Add --multiple-threshold for returning all peaks above a score

README.md

@@ -51,10 +51,12 @@ The following command-line options can be provided to alter the behaviour of the
 | --within audio file  |  ...within this file |
 | --sr sample rate |  Target sample rate in Hz during downsampling (default: 8000) |
 | --trim seconds  |  Only use the first n seconds of each audio file |
+| --start seconds  |  Skip the first n seconds of each audio file before processing.  Combined with `--trim`, the audio considered is the window `[start, start + trim]`. Default: 0. |
 | --resolution samples  |  Resolution (maximum accuracy) of search in samples (default: 128) |
 | --show-plot  |  Display a plot of the cross-correlation results |
 | --save-plot filename |  Save a plot of the cross-correlation results to a file (in a format that matches the extension you provide - png, ps, pdf, svg) |
 | --json  |  Output in JSON for further processing |
+| --multiple-threshold score  |  Return all correlation peaks with a standard score above this threshold (instead of just the single best peak).  JSON output becomes an array of `{time_offset, standard_score}` objects, sorted by score descending. |
 
 You can fine-tune the results for your application by tweaking the sample rate, trim and resolution parameters:
 * The _sample rate_ option refers to a resampling operation that is carried out before the audio offset search is carried out.  It does not refer to the sample rate(s) of the audio files being compared.  Resampling at a higher sample rate retains higher audio frequencies, but increases the time required to search for an offset.  The default sample rate is 8000Hz, which is a good compromise for most audio.

audio_offset_finder/audio_offset_finder.py

@@ -16,6 +16,7 @@
 
 from subprocess import Popen, PIPE
 from scipy.io import wavfile
+from scipy.signal import find_peaks
 import librosa
 import os
 import tempfile
@@ -38,7 +39,9 @@ def mfcc(audio, win_length=256, nfft=512, fs=16000, hop_length=128, numcep=13):
     ]
 
 
-def find_offset_between_files(file1, file2, fs=8000, trim=None, hop_length=128, win_length=256, nfft=512, max_frames=2000):
+def find_offset_between_files(
+    file1, file2, fs=8000, trim=None, start=0, hop_length=128, win_length=256, nfft=512, max_frames=2000
+):
     """Find the offset time offset between two audio files.
 
     This function takes in two file paths, and (assuming they are media files with a valid audio track)
@@ -55,6 +58,9 @@ def find_offset_between_files(file1, file2, fs=8000, trim=None, hop_length=128,
         The sampling rate that the audio should be resampled to prior to MFCC calculation, in Hz
     trim: int
         The length to which input files will be truncated before processing, in seconds.  A value of "None" indicates no trimming.
+    start: float
+        The number of seconds to skip at the beginning of each input file before processing.  Defaults to 0.
+        When combined with "trim", the audio considered is the window [start, start + trim].
     hop_length: int
         The number of samples (at the resampled rate "fs") to skip between each calculated MFCC frame
     win_length: int
@@ -78,8 +84,8 @@ def find_offset_between_files(file1, file2, fs=8000, trim=None, hop_length=128,
     ------
     InsufficientAudioException if the audio supplied is too short to analyse.
     """
-    tmp1 = convert_and_trim(file1, fs, trim)
-    tmp2 = convert_and_trim(file2, fs, trim)
+    tmp1 = convert_and_trim(file1, fs, trim, start=start)
+    tmp2 = convert_and_trim(file2, fs, trim, start=start)
     a1 = wavfile.read(tmp1, mmap=True)[1].astype(float)
     a2 = wavfile.read(tmp2, mmap=True)[1].astype(float)
     offset_dict = find_offset_between_buffers(a1, a2, fs, hop_length, win_length, nfft)
@@ -164,6 +170,50 @@ def find_offset_between_buffers(buffer1, buffer2, fs, hop_length=128, win_length
     }
 
 
+def find_peaks_in_correlation(results, threshold):
+    """Find all local maxima in a cross-correlation curve with a standard score above the given threshold.
+
+    Parameters
+    ----------
+    results: dict
+        A results dictionary as returned by find_offset_between_files() or find_offset_between_buffers().
+    threshold: float
+        The minimum standard score that a peak must have to be returned.
+
+    Returns
+    -------
+    A list of dicts, sorted by standard_score in descending order.  Each dict contains:
+    time_offset (float), frame_offset (int), standard_score (float)
+    """
+    c = results["correlation"]
+    time_scale = results["time_scale"]
+    latest_frame_offset = results["latest_frame_offset"]
+
+    mean = np.mean(c)
+    std = np.std(c)
+    if std < 1e-10:
+        return []
+
+    height = mean + threshold * std
+    peak_indices, _ = find_peaks(c, height=height)
+
+    peaks = []
+    for idx in peak_indices:
+        frame_offset = int(idx)
+        if frame_offset > latest_frame_offset:
+            frame_offset -= len(c)
+        score = float((c[idx] - mean) / std)
+        peaks.append(
+            {
+                "time_offset": frame_offset * time_scale,
+                "frame_offset": frame_offset,
+                "standard_score": score,
+            }
+        )
+    peaks.sort(key=lambda p: p["standard_score"], reverse=True)
+    return peaks
+
+
 # returns an array in which the first half represents an offset of mfcc2 within mfcc1,
 # and the second half (accessed by negative indices) vice-versa.
 def cross_correlation(mfcc1, mfcc2, nframes):
@@ -209,7 +259,7 @@ def std_mfcc(array):
     return (array - np.mean(array, axis=0)) / np.std(array, axis=0)
 
 
-def convert_and_trim(afile, fs, trim=None):
+def convert_and_trim(afile, fs, trim=None, start=0):
     """Converts the input media to a temporary 16-bit WAV file and trims it to length.
 
     Parameters
@@ -221,6 +271,9 @@ def convert_and_trim(afile, fs, trim=None):
     trim: float
         The length to which the output audio should be trimmed, in seconds.  (Audio beyond this point will be discarded.)
         A value of "None" implies no trimming.
+    start: float
+        The number of seconds to skip at the beginning of the audio.  Audio before this point will be discarded.
+        Defaults to 0.
 
     Returns
     -------
@@ -235,7 +288,7 @@ def convert_and_trim(afile, fs, trim=None):
     ffmpeg_command += ["-i", afile]
     ffmpeg_command += ["-ac", "1"]
     ffmpeg_command += ["-ar", str(fs)]
-    ffmpeg_command += ["-ss", "0"]
+    ffmpeg_command += ["-ss", str(start)]
     if trim:
         ffmpeg_command += ["-t", str(trim)]
     ffmpeg_command += ["-acodec", "pcm_s16le"]

audio_offset_finder/cli.py

@@ -14,7 +14,7 @@
 # See the License for the specific language governing permissions and
 # limitations under the License.
 
-from .audio_offset_finder import find_offset_between_files
+from .audio_offset_finder import find_offset_between_files, find_peaks_in_correlation
 import argparse
 import sys
 
@@ -34,6 +34,13 @@ def main(argv):
     parser.add_argument("--within", metavar="audio file", type=str, help="...within this file.")
     parser.add_argument("--sr", metavar="sample rate", type=int, default=8000, help="Resample to this rate before searching")
     parser.add_argument("--trim", metavar="seconds", type=int, help="Only consider the first n seconds of the audio files")
+    parser.add_argument(
+        "--start",
+        metavar="seconds",
+        type=float,
+        default=0,
+        help="Skip the first n seconds of each audio file before processing. Combined with --trim, considers the window [start, start+trim].",
+    )
     parser.add_argument(
         "--resolution", metavar="samples", type=int, default=128, help="Resolution (maximum accuracy) of search in samples"
     )
@@ -46,6 +53,16 @@ def main(argv):
         help=("Save a plot of cross-correlation results to a file " "(format matches extension - png, ps, pdf, svg)"),
     )
     parser.add_argument("--json", action="store_true", dest="output_json", help="Output in JSON for further processing")
+    parser.add_argument(
+        "--multiple-threshold",
+        metavar="score",
+        type=float,
+        dest="multiple_threshold",
+        help=(
+            "Instead of returning the single best peak, return all peaks with a standard score above this threshold. "
+            "JSON output becomes an array of objects; plain output lists each peak."
+        ),
+    )
     args = parser.parse_args(argv)
     if not (args.find_offset_of and args.within):
         parser.error("Please provide input audio files")
@@ -56,23 +73,40 @@ def main(argv):
             trim = int(args.trim)
 
         results = find_offset_between_files(
-            args.within, args.find_offset_of, fs=int(args.sr), trim=trim, hop_length=int(args.resolution)
+            args.within,
+            args.find_offset_of,
+            fs=int(args.sr),
+            trim=trim,
+            start=args.start,
+            hop_length=int(args.resolution),
         )
     except Exception as e:
         print(e, file=sys.stderr)
         return 1
 
+    peaks = None
+    if args.multiple_threshold is not None:
+        peaks = find_peaks_in_correlation(results, args.multiple_threshold)
+
     if args.output_json:
         import json
 
-        json_results = {"time_offset": results["time_offset"], "standard_score": results["standard_score"]}
+        if peaks is not None:
+            json_results = [{"time_offset": p["time_offset"], "standard_score": p["standard_score"]} for p in peaks]
+        else:
+            json_results = {"time_offset": results["time_offset"], "standard_score": results["standard_score"]}
         print(json.dumps(json_results))
     else:
-        print("Offset: %s (seconds)" % str(results["time_offset"]))
-        print("Standard score: %s" % str(results["standard_score"]))
+        if peaks is not None:
+            print("Found %d peak(s) with standard score >= %s:" % (len(peaks), str(args.multiple_threshold)))
+            for p in peaks:
+                print("  Offset: %s (seconds), Standard score: %s" % (str(p["time_offset"]), str(p["standard_score"])))
+        else:
+            print("Offset: %s (seconds)" % str(results["time_offset"]))
+            print("Standard score: %s" % str(results["standard_score"]))
 
     if args.show_plot or args.plot_file is not None:
-        plot_results(args, results)
+        plot_results(args, results, peaks=peaks)
 
 
 # Re-order the cross-correlation array so that the index of the earliest frame offset is at one end of the range
@@ -83,7 +117,7 @@ def reorder_correlations(cc, earliest_frame_offset):
     return concatenate((cc[earliest_frame_offset:], cc[:earliest_frame_offset]))
 
 
-def plot_results(args, results):
+def plot_results(args, results, peaks=None):
     import matplotlib.pyplot as pyplot
     import matplotlib.ticker as ticker
 
@@ -109,8 +143,12 @@ def plot_results(args, results):
     plot_title = "Offset of %s in %s" % (args.find_offset_of, args.within)
     pyplot.title(plot_title, fontsize=14)
 
-    peak_xvalue = results["frame_offset"]
-    pyplot.axvline(x=peak_xvalue, color="red", linestyle="dotted")
+    if peaks is not None:
+        for p in peaks:
+            pyplot.axvline(x=p["frame_offset"], color="red", linestyle="dotted")
+    else:
+        peak_xvalue = results["frame_offset"]
+        pyplot.axvline(x=peak_xvalue, color="red", linestyle="dotted")
 
     if args.plot_file is not None:
         pyplot.savefig(args.plot_file)

tests/audio_offset_finder_test.py

@@ -16,7 +16,7 @@
 
 import pytest
 from audio_offset_finder.audio_offset_finder import find_offset_between_files, std_mfcc, cross_correlation
-from audio_offset_finder.audio_offset_finder import InsufficientAudioException
+from audio_offset_finder.audio_offset_finder import InsufficientAudioException, find_peaks_in_correlation
 import numpy as np
 import os
 
@@ -86,6 +86,44 @@ def test_std_mfcc():
     np.testing.assert_array_equal(std_mfcc(m), np.array([[-1.0 / s1, -1.0 / s2, -0.5 / s3], [1.0 / s1, 1.0 / s2, 0.5 / s3]]))
 
 
+def test_find_offset_with_start():
+    # Skipping the first N seconds identically in both files preserves the relative offset.
+    results = find_offset_between_files(path("timbl_1.mp3"), path("timbl_2.mp3"), hop_length=160, trim=25, start=5)
+    assert results["time_offset"] == pytest.approx(12.26)
+    assert results["standard_score"] > 10
+
+    # Auto-correlation with a non-zero start still finds offset 0.
+    results = find_offset_between_files(path("timbl_1.mp3"), path("timbl_1.mp3"), hop_length=160, trim=20, start=10)
+    assert results["time_offset"] == pytest.approx(0.0)
+    assert results["standard_score"] > 10
+
+
+def test_find_peaks_in_correlation():
+    results = find_offset_between_files(path("timbl_1.mp3"), path("timbl_2.mp3"), hop_length=160, trim=35)
+
+    # A very high threshold returns just the dominant peak
+    peaks = find_peaks_in_correlation(results, threshold=20.0)
+    assert len(peaks) >= 1
+    assert peaks[0]["time_offset"] == pytest.approx(12.26)
+    assert peaks[0]["standard_score"] == pytest.approx(28.99, rel=1e-2)
+    # Sorted by standard_score descending
+    for i in range(len(peaks) - 1):
+        assert peaks[i]["standard_score"] >= peaks[i + 1]["standard_score"]
+    # Above the threshold
+    for p in peaks:
+        assert p["standard_score"] >= 20.0
+
+    # A lower threshold returns more peaks
+    more_peaks = find_peaks_in_correlation(results, threshold=5.0)
+    assert len(more_peaks) >= len(peaks)
+    for p in more_peaks:
+        assert p["standard_score"] >= 5.0
+
+    # An impossibly-high threshold returns nothing
+    none_peaks = find_peaks_in_correlation(results, threshold=1000.0)
+    assert none_peaks == []
+
+
 def test_cross_correlation():
     m1 = np.array([[-0.5, -0.4, -0.4], [0.5, 0.5, 0.4], [0.1, -0.1, 0.1]])
     m2 = np.array([[0.5, 0.5, 0.4], [0.1, -0.1, 0.1], [-0.6, 0.0, -0.3]])

tests/tool_test.py

@@ -71,3 +71,53 @@ def test_json():
         assert len(json_array) == 2
         assert pytest.approx(json_array["time_offset"]) == 12.26
         assert pytest.approx(json_array["standard_score"], rel=1e-2) == 28.99
+
+
+def test_multiple_threshold_json():
+    import json
+
+    args = (
+        "--find-offset-of tests/audio/timbl_2.mp3 --within tests/audio/timbl_1.mp3 --resolution 160 "
+        "--trim 35 --json --multiple-threshold 20.0"
+    )
+    with patch("sys.stdout", new=StringIO()) as fakeStdout:
+        main(args.split())
+        output = fakeStdout.getvalue().strip()
+        peaks = json.loads(output)
+        assert isinstance(peaks, list)
+        assert len(peaks) >= 1
+        # First peak is the dominant one
+        assert pytest.approx(peaks[0]["time_offset"]) == 12.26
+        assert pytest.approx(peaks[0]["standard_score"], rel=1e-2) == 28.99
+        for p in peaks:
+            assert "time_offset" in p
+            assert "standard_score" in p
+            assert p["standard_score"] >= 20.0
+
+
+def test_start():
+    import json
+
+    args = (
+        "--find-offset-of tests/audio/timbl_2.mp3 --within tests/audio/timbl_1.mp3 --resolution 160 "
+        "--trim 25 --start 5 --json"
+    )
+    with patch("sys.stdout", new=StringIO()) as fakeStdout:
+        main(args.split())
+        output = fakeStdout.getvalue().strip()
+        result = json.loads(output)
+        assert pytest.approx(result["time_offset"]) == 12.26
+
+
+def test_multiple_threshold_plain():
+    args = (
+        "--find-offset-of tests/audio/timbl_2.mp3 --within tests/audio/timbl_1.mp3 --resolution 160 "
+        "--trim 35 --multiple-threshold 20.0"
+    )
+    with patch("sys.stdout", new=StringIO()) as fakeStdout:
+        main(args.split())
+        output = fakeStdout.getvalue().strip()
+        assert "peak" in output.lower()
+        assert "12.26" in output
+        assert "Offset:" in output
+        assert "Standard score:" in output