How might I go about detecting phase changes in a sine wave of a certain frequency?

Question

I want to detect occasional instantaneous changes in the phase of a sine wave.

I have a signal that is a sine wave whose frequency is in the kHz range and experiences random changes in its phase. I want to detect these phase changes and count the occurrences. Ideally, I would also be able to extract some approximate phase change in degrees for each occurrence. How might I go about that?

I implemented Jason R's suggestion in Python:

# generating a 30KHz sine wave sampled at 500KHz with a pi/2 phase shift halfway through the signal
F = 30e3
Fs = 500e3
t = np.arange(0, 1e-3, 1/Fs)
S = np.append(np.sin(F*2*np.pi*t[:int(len(t)/2)] + 0), np.sin(F*2*np.pi*t[int(len(t)/2):] - np.pi/2))

S_analytic =  hilbert(S) # calculating the analytic signal

def GetRealImagArray(Array):
    ImagArray = np.array([num.imag for num in Array])
    RealArray = np.array([num.real for num in Array])
    return RealArray, ImagArray

R, I = GetRealImagArray(S_analytic)
phi = np.arctan2(I, R)
f = phi[1:] - phi[:-1]

I then get the following by plotting the signal and FM discriminator:

The large positive spike indicates the change in phase as required in my question. What is the cause of the downward spikes at the minima of the sine wave? Is this expected behaviour or an error in my Python implementation?

Answer 1

The easiest way to do this would be to apply a digital FM discriminator to the signal. At a high level, you would do the following:

• Form the analytic signal of your input. This also can be termed as ensuring it is at complex baseband.

• Calculate the phase versus time of your analytic signal $x[n]$: $$ \phi[n] = \text{atan2}\left(\frac{\text{Im}\{x[n]\}}{\text{Re}\{x[n]\}}\right) $$

• The instantaneous frequency of the signal versus time can be estimated by differentiating $\phi[n]$ using any discrete differentiation method (first-order difference is usually used, i.e. $f[n] = \phi[n]-\phi[n-1]$).

Why does this work? If your signal is typically a constant frequency, your FM discriminator's output $f[n]$ should be approximately constant. If there is a phase jump at a particular sample index, this would appear to the discriminator to be a sudden change to a much different frequency. After the jump, the discriminator output $f[n]$ will again track the sinusoid's frequency.

Therefore, you can run a discriminator continuously on your input signal and apply a threshold, looking for large changes in frequency from sample to sample. That would correspond to a phase discontinuity.

Edit: The above is subject to spurious outputs when the phase value output by the arctangent function wraps around $\pm \pi$ (see the comment below). A more robust way of implementing the above would instead be:

$$ \hat{x}[n] = x[n] x^*[n-1] \\ f[n] = \text{atan2}\left(\frac{\text{Im}\{ \hat{x}[n]\}}{\text{Re}\{ \hat{x}[n]\}}\right) $$

That is, form the product of $x[n]$ with the conjugate of $x[n-1]$. At each sample instant $n$, the result will have a phase that is equal to the difference in phase between $x[n]$ and $x[n-1]$ (which is what you're trying to measure). The subtraction alleviates any wrapping issues when the absolute phase of $x[n]$ is near $\pm \pi$. Just find the phase of the resulting product to get $f[n]$.