I have a few sinusoidal signals at approximately the same frequency. I would like to know the phase difference between the three signals. What is the best way to do this? It would also be helpful to attach an uncertainty to these phase measurements. I've attached a picture of the signals
The obvious answer would be to take the FFT of each signal. The phase angle of the highest (magnitude) bin gives you the phase of the fundamental sinewave buried in each signal. Compare and contrast as needed.
Use an appropriate window function to minimize edge effects.
I recommend transfer-function estimation techniques such as the Welch Periodogram, used in conjunction with coherence estimates, which can be used to put confidence intervals around your estimate. I recommend Bendat and Piersol, Random Data: Analysis and Measurement Procedures for a very thorough treatment.
The Welch periodogram involves taking the spectrum of overlapping segments of data to yield the periodogram.
Coherence pretty much indicates the noise and nonlinearities inherent in the relationship between input and output. A coherence of 1.0 reflects a perfectly linear noise-free system. Both noise and nonlinearity reduce coherence.
You can also plot one signal against the other, and calculate phase from the resulting Lissajou figure
Another option would be to examine the peaks in the cross correlation between the two signals.
That appears to be a matplotlib graph correct? This would imply you are using python to post-process. Equally you want to post-process already collected data.
What I would suggest is apply a moving average filter to each signal. While such a filter will "phase shift" an individual waveform with respect to the original, relative to each other (post-filter) there will be no difference as the same "shift" will have been applied
1/4x + 1/4(x-1) + 1/4(x-2) + 1/4(x-3) should do or a simple numpy.mean with a fix slice window and varying index.
Either way once you have lightly filtered the data to minimise zero crossing jitter you can utilise a couple of numpy functions:
zero_crossings = numpy.where(numpy.diff(numpy.sign(DATA)))
This will return an array of all the zero crossings. You may need to increase the filtering if there are multiple very close together.
If you have matlab a similar scheme can be achieved using the filter command & there is a crossing.m file publically available from mathworks exchange.
NOTE this won't give a distinction for positive going or negative going crossings.
Take this piece of python code:
#!/usr/bin/env python3 import numpy as np from pylab import * import matplotlib matplotlib.pyplot.switch_backend('QT4Agg') def filt(data, l=4): tmp = np.zeros(len(data)) for i in range(l,len(data)): tmp[i] = np.mean(data[i-l:i]) return tmp t = np.arange(0,10*np.pi,0.1e-3) s = np.sin(t) + np.random.normal(0,0.01,len(t)) crossings = np.where(np.diff(np.sign(s))) print(len(crossings)) s2 = filt(s,l=4) crossings = np.where(np.diff(np.sign(s2))) print(len(crossings)) s3 = filt(s,l=10) crossings = np.where(np.diff(np.sign(s3))) print(len(crossings)) s4 = filt(s,l=150) crossings = np.where(np.diff(np.sign(s4))) print(len(crossings)) filt(s) plot(t,s) grid(True) hold(True) plot(t,s2) plot(t,s3) plot(t,s4,linewidth='2') show()
This generates a 300,000 datapoint sin wave and couples onto it a normally distributed noise source with a sd of 0.01
As you can see a lot of zero crossing noise, an annoying amount.
The rest of the script does a really crude moving average (a better implementation would be a movingAverageFilter & this is easy in a matlab script via the "filter" command. (http://www.mathworks.co.uk/help/matlab/data_analysis/filtering-data.html)
Each iteration uses a different window size & equally prints out the total number of zero crossings detected.
[jrb@FluidMotion:~/tmp] 27s 1 $ ./test.py 1126 272 120 12
The number of detected zero crossings has fallen from 1126 to 12
As mentioned there is an absolute phaseshift (would be significantly reduced if a better filter was used, but this was a 5min example), but for comparing relative differences it is valid as all datasets are shifted by the same amount.
Once you know the indexes of the zero crossings and you know the timestep of one sample, you are able to calculate the time difference between zero crossings. Couple that with the number of samples between three zero crossings of the same data & you know the number of samples per 360degrees
This took 150samples out of 300,000 (to produce a reasonable zero-crossing detection but the source had alot of random zero crossings. The OP data appears to only have a couple of zero-crossing events.
You could aggressively band-pass filter the signals, since it sounds like you know the approximate frequency you're interested in. Low pass to get rid of the noise, and high pass because different peaks of a single signal are at wildly differing amplitudes. Then compare the times of the zero crossings.