Find when two lists are parallel with autocorrelation

Question

I am working on a signal processing project with two oscillators. The first oscillator has a frequency $f_1$ that is independent, and the second oscillator has a frequency $f_2$ that increases, and will become 'frequency entrained' (i.e., 'synchronized', or 'phase locked') when its near to $f_1$.

The best way to determine if the two oscillators are phase locked is to compare the phases of both signals. The real data is messy and a bit cumbersome, but a minimum working example might be two data sets like this:

phases1 = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20}
phases2 = {0., 0.5, 1., 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.5, 7.5, 8.5, 9.5, 10.5, 12.5, 14.5, 16.5, 18.5}

There is a spot in the middle where the two lines are parallel. One way to find where they are parallel is to take the derivative of both lists and subtract:

Differences[phases1] - Differences[phases2]

And then use a list manipulation function to locate the zeros. However, for the system I am working with, this isn't good enough because of noise and other issues. Is there a better way to determine when two lists are parallel in mathematica?

Should I use autocorrelation or some other function?

Answer 1

Here is the basic analog engineer approach. Multiply the outputs of your oscillators and apply a lowpass filter:

    phases1 = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 
  17, 18, 19, 20}
    phases2 = {0, 0.5, 1., 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.5, 
  7.5, 8.5, 9.5, 10.5, 12.5, 14.5, 16.5, 18.5}  

k = 10; (* k = it scales the frequencies; your frequencies are certainly  
different; can be a problem *)

(* ouputs of the oscillators :*)
osc1output[t_] = 
  Sin[k Interpolation[phases1, InterpolationOrder -> 1][t]];
osc2output[t_] = 
  Sin[k Interpolation[phases2, InterpolationOrder -> 1][t]];
tfin = 18.4;

(* graphics rendering : *)
options = {Frame -> True, ImageSize -> 300, Axes -> False, 
   ImagePadding -> 40};
Column[{Plot[osc1output[t], {t, 1, tfin}, Evaluate @ options, 
   FrameLabel -> {"time", "", Style["osc1output", FontSize -> 20]}],
  Plot[osc2output[t], {t, 1, tfin}, Evaluate @ options, 
   FrameLabel -> {"time", "", Style["osc2output", FontSize -> 20]}],
  Plot[osc1output[t]* osc2output[t], {t, 1, tfin}, Evaluate @ options,
    FrameLabel -> {"time", "", 
     Style["osc1output * osc2output", FontSize -> 20]}]}, 
 Dividers -> None]

Low pass filter, adjusted a posteriori:

mySampleRate = 100; (* should be larger than your frequencies, say at least  
                       ten times *)
l1 = Table[osc1output[t] osc2output[t], {t, 1, tfin, 1./mySampleRate}];
Manipulate[
 ListPlot[LowpassFilter[l1, ωc, n, SampleRate -> mySampleRate],
   PlotRange -> {Automatic, {-1, 1}}, 
  DataRange -> {1, tfin}], {{ωc, 8.2}, 0, 100}, {{n, 139}, 0, 
  400, 1}]

When the frequencies are the same, the filtered signal is constant, and not necessarily zero (here, it's a coincidence).

Answer 2

This is very similar to eldo's answer, but IMO a little more robust.

data1 = Range[0, 20];
data2 = Join[Range[0., 5.5, .5], Range[6.5, 10.5], Range[12.5, 18.5, 2.]];

getIndices[data1_, data2_, threshold_: 1.*^-3] :=
  Module[{indices},
    indices =
      MapIndexed[
        If[Abs[#1] < threshold, #2, Nothing] &, 
        Differences[Thread[Subtract[data1, data2]]]];
    Append[indices, Last[indices] + 1]]

With[{indices = getIndices[data1, data2]},
  ListLinePlot[
    {Transpose[{data1, data1}],
     Transpose[{data1, data2}],
     Transpose[{Flatten[indices] - 1, Extract[data2, indices]}]},
    AspectRatio -> Automatic,
    PlotStyle -> {Automatic, Automatic, AbsoluteThickness[5]}]]

Answer 3

p1 =
  {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20};

Some noise added

p2 = 
  {0., 0.5, 1., 1.5, 2., 2.5, 3., 3.5, 4., 4.5, 5., 5.5, 6.6, 7.4, 
   8.5, 9.51, 10.45, 12.5, 14.5, 16.5, 18.5};

dif = Abs[Differences@p1 - Differences@p2];

Position of differences <= n to allow for noise

pos = Position[dif, a_Real /; a <= 0.25];

"Parallel" section

p3 = Transpose[{Flatten@pos, Extract[p2, pos]}]

 ListLinePlot[{p1, p2, p3},
   GridLines -> Automatic,
   PlotStyle -> {Automatic, Automatic, [email protected]}]