SimuLib: def integer_beacon_sync

lib/timing.py

@@ -95,12 +95,21 @@ if not hasattr(signal, 'correlation_lags'):
 ##### end of monkey patch correlation_lags
 
 def beacon_time_delay(samplerate, ref_beacon, beacon):
+    """
+    Determine the time delay between two beacons using correlation.
+
+    """
     grid = correlation_grid(in1_len=len(ref_beacon), in2_len=len(beacon), mode='full')
     time_lag, errs = lag_gridsearch(grid, samplerate, ref_beacon, beacon)
 
     return time_lag, errs
 
 def beacon_phase_delay(samplerate, f_beacon, ref_beacon, beacon):
+    """
+    Determine total phase delay between two beacons using correlation.
+
+    Internally uses beacon_time_delay.
+    """
     time_delay, errs = beacon_time_delay(samplerate, ref_beacon, beacon)
 
     phase = 2*np.pi*f_beacon*time_delay
@@ -108,6 +117,32 @@ def beacon_phase_delay(samplerate, f_beacon, ref_beacon, beacon):
     
     return phase, phase_err
 
+def beacon_integer_period(samplerate, f_beacon, ref_impulse, impulse, k_step=1):
+    return _beacon_integer_period_sum(samplerate, f_beacon, ref_impulse, impulse, k_step=k_step)
+
+def _beacon_integer_period_sum(samplerate, f_beacon, ref_impulse, impulse, k_step=1):
+    """
+    Use the maximum of a coherent sum to determine
+    the best number of periods of f_beacon.
+    """
+    max_k = int( len(ref_impulse)*f_beacon/samplerate )
+    ks = np.arange(0, max_k, step=k_step)
+
+    maxima = np.empty(len(ks))
+
+    best_i = 0
+
+    for i,k in enumerate(ks, 0):
+        augmented_impulse = util.time_roll(impulse, samplerate, k/f_beacon)
+
+        maxima[i] = max(ref_impulse + augmented_impulse)
+
+        if maxima[i] > maxima[best_i]:
+            best_i = i
+
+    return ks[best_i], (ks, maxima)
+
+
 def lag_gridsearch(grid, sample_rate, reference, signal_data):
     """
     Return the best time shift found when doing a grid search.