# Confirmation on how to calculate phase differences at every frequency point between two wideband signals?

I have a problem that I thought was going to be simple, but it has become surprisingly stubborn and I question my method...

I have used the method (described below), but I wanted to confirm that this was the right approach, and see if there was another method I could use.

What I basically have is two wideband signals, (on the order of 100 Mhz BW), in quadrature. So I have the I and Q of those two signals. What I am trying to find is the phase difference between those two signals as a function of frequency.

The method I am currently trying is the following:

1. I take my quadrature complex signals, and window each of them with a hamming window 3 times over.
2. I take the FFT of the windowed complex signals, at some arbitrary length, say 2048. Call the results $S_1(w)$ and $S_2(w)$.
3. I multiply $S_1(w)$ by the conjugate of $S_2(w)$.
4. I take the argument out of the complex exponential vector from the above step (inverse tangent of imag over real). This gives me the phase differences between the two signals as a function of frequency.

My questions on this method:

• Is this the right approach for this simple type of problem?
• When I tried this method for a simulation with 1, 2 or a couple tones, I can get good results. I find the index of where my tones exist, (look at the abs of the FFT), I take that index, and look at the $\Delta$phase value from step 4 above. This gives me the right phase.
• However, when I look the $\Delta$phase values for frequency indicies where said those frequencies dont exist, I get wild answers. Why is this? I realize those frequencies dont exist but shouldn't the phases just be 0 - 0 = 0 in that case?
• Based on this (with a couple tones), can I readily extend this to a wideband signal and also expect good results?
• My last question regards the nature of the phase delay here - physically what is happening is that one of the signals (optical) is going through a device that delays different frequencies differently. Of course the delays are manifested as phase offsets between all frequencies... but at the same time there is a time delay associated with this as well. Even if we were to measure the $\Delta$phase between the two signals properly, wouldnt it still be impossible to measure the time delay between them at those frequencies? I can have a time delay of 10 wavelenghts, but a phase offset of 0.

Thoughts?

Thanks in advance!

EDIT FROM COMMENT FEEDBACKS:

I think I am clearer now as to what my question really is: For a case where you have well separated tones and high SNR, my method can work.

However, what to do in the wideband case? In this sense, suppose $s_1[n]$ and $s_2[n]$ both contain two tones that are very close to each other in the frequency space. But your FFT can only give you so much frequency resolution. So instead of applying my method to two peaks, you only have one peak to go on. Now what?

• What relation would the $\Delta$phase of that peak have to the actual $\Delta$phase of the 'true' peaks? Is it their average?
• One method I am thinking of is correlating the two signals in the time domain, and working backwards to get the phase offsets per frequency taking advantage of the fourier translation property, although part of me thinks it might not work since there is no new information, just looking at the problem from a different domain.
• Another thought I had for the wideband case would be to put the whole thing through many narrow band filters, and then apply my original method. Is there a name for something like this?

Thanks again

• @DilipSarwate "I am not sure the phases if you are extracting are representative of the phase difference between the tones corresponding to those bins." I think that the phases difference between them is in fact this. At the end of the day the way I see it, I am comparing the phase of one phasor rotating at some frequency, against a phasor of the same frequency but with different starting phase. I dont see why this wouldn't be true. – Spacey Jan 12 '12 at 4:07
• @DilipSarwate "And what happens if the frequencies of the tone(s) you have chosen fall between bins?" Here I agree with you - I too 'worry' about this scenario, that I have some frequency out there that cannot be represented by the FFT and the number of samples I have. This I think is the biggest problem. That being said, I had heard anecdotaly from somewhere that you can orthogonalise different frequency bins by gram-schmit? Or something along those lines?... – Spacey Jan 12 '12 at 4:10
• If you rereference phase to the center of the FFT aperture (via fftshift; or flipping every other result bin) you can easily interpolate the phase of frequencies that fall between between bins. – hotpaw2 Jan 12 '12 at 7:43
• I am trying to do the same thing. can Hilbert transform help us out in this regard? – user5686 Oct 16 '13 at 19:09
• @Spacey Did you ever successfully finish this project? – JZYL Aug 30 '19 at 21:43

## 2 Answers

What you seem to want is the phase response of the device that your optical system is passing through. Based on that, your high-level approach makes intuitive sense overall; I'm not sure why you chose a Hamming window and applied it three times. In any event, though, if you have a linear discrete-time system with (complex-valued) frequency response $H(e^{j\omega})$, then:

$$Y(e^{j\omega}) = H(e^{j\omega})X(e^{j\omega})$$

where $X(e^{j\omega})$ and $Y(e^{j\omega})$ are the discrete-time Fourier transforms (DTFTs) of the input and output signals, respectively. If you express the above in polar format (magnitude and phase), then the phase response of the system is just:

$$\angle H(e^{j\omega}) = \angle Y(e^{j\omega}) - \angle X(e^{j\omega})$$

that is, the phase response of the system just expresses the amount of phase shift that the system applies to the input signal as a function of (normalized) frequency $\omega$. Your approach takes advantage of the fact that for two complex numbers $X(e^{j\omega})$ and $Y(e^{j\omega})$:

$$\angle Y(e^{j\omega}) - \angle X(e^{j\omega}) = \angle \left(X^*(e^{j\omega})Y(e^{j\omega})\right)$$

So, given the DTFTs of the system's input and output signals, you can straightforwardly calculate the phase response that you seek. The devil is in the details of how one would obtain the DTFT representations of these signals. The nice part about the DTFT is that it is continuous in frequency; the bad part is that the DTFT sum is infinite in length. So, if you really want to evaluate the phase response at arbitrarily many frequencies, you better be patient, because you'll need to collect an infinite amount of data first.

The discrete Fourier transform (DFT) is often used as an approximation to the DTFT, which is valid under certain conditions. If the signals $x[n]$ and $y[n]$ are finite in length (which for a practical scenario like yours, they are), then the DFT can be thought of as evaluation of the DTFT on the finitely-spaced frequency grid $\frac{2\pi k}{N}, k = 0, 1, \ldots , N-1$ (where $N$ is the signal length). So, the larger you make $N$, the smaller your sample spacing across the frequency axis.

Most physical devices will have a frequency response that is smooth to some degree (taken with the caveat that I'm not knowledgeable on optical device physics), so you should be able to find a DFT length that gives a decent approximation to the DTFTs for the input and output signals of your system. You can then apply your approach to estimate the phase response of the device.

You noted in your question that you see incorrect results for the estimated phase response for DFT bins in the output signal that didn't contain any content in the input signal; this is to be expected. Your assumption that you would just get zero phase for those bins is not likely to hold in a practical scenario; noise present in the sampled signals, spectral leakage, or other interference can cause you to observe some small amount of energy at frequencies that you wouldn't expect it. These extraneous sources can introduce distortion between the input and output that is not due to the device's response, and for complex numbers with small magnitude, small amounts of additive noise can drastically affect the number's resulting phase.

Therefore, you can only consider the phase measurement valid at frequencies where there is significant spectral content in the input. Stated differently, if you want to make a precise phase measurement, you want high signal-to-noise ratio in the complex value that you're using to estimate the phase.

• "So, the larger you make N, the smaller your sample spacing across the frequency axis." Yes, I concur. In this case, the 'N' you are referring to here corresponds to the number of actual data points in my original time vector, and NOT the N being the size of my zero-padded time vector before FFT correct? (The former will have more data, whereas the latter just interpolates my frequency bins with no additional information). Have I understood you correctly? contd... – Spacey Jan 12 '12 at 20:59
• " I'm not sure why you chose a Hamming window and applied it three times." Well, I really wanted to minimize any effects from artifacts due to the FFT grid spacing being the way it was, so I windowed the signals once, and then when windowing again I got better results. So I windowed a third time for giggles although I dont think it made a difference compared to doing it twice. – Spacey Jan 12 '12 at 21:02
• (See edit to question) – Spacey Jan 12 '12 at 21:04
• Your ability to resolve peaks close to one another in frequency is limited by the time duration of your observation (i.e. FFT zero-padding will not help; you need more actual samples). The DFT is nothing but a bank of uniformly-spaced filters with rectangular impulse response, so given two tones at specific frequency and phase offsets from a bin center, you can calculate what the effective phase of the DFT output would be. You're not going to be able to make any meaningful conclusions as to each tone's phase from the composite measurement, however. – Jason R Jan 12 '12 at 21:36
• Regarding your last question in your edit, yes, there is a name for such a method, where you pass the signal through a bunch of narrow filters: a longer FFT. The DFT can be thought of exactly as a critically-sampled filter bank; as you increase the transform length, you lengthen the filters' impulse responses (effectively narrowing each filter's bandwidth) and space them more finely on the frequency axis. – Jason R Jan 12 '12 at 21:38

Noise pretty much has a random phase. So it's usually reasonable to just clamp the phase of any bin with a magnitude below some noise threshold to zero. Make sure this threshold is above the (unavoidable) numerical noise level.