How do I find the phase for a periodic cos wave using fft?
6 次查看(过去 30 天)
显示 更早的评论
I have a simple cos wave that has a period, in this case 6. I am offsetting the wave by an angular amount with phase. How can I use the fft function to recover back my initial phase offset? I tried converting the fourier transformed value with the highest amplitude using abs (which does correclty identify the period) back to a phase using angle, but this never seems to work. If someone could enlighten me as to how to do this, I would be so grateful! (Edit: I would like to see how to correctly recover the original period as well)
n = 1e3;
a = linspace(0, 2*pi, n);
period = 6;
phase = rand * 2*pi;
f = rescale(cos((a-phase)*pr));
3 个评论
Mathieu NOE
2024-6-28
sure
NB that in the other post, the author also wanted to use fft, but did not really succeed
采纳的回答
David Goodmanson
2024-6-28
编辑:David Goodmanson
2024-6-28
Hi Noam,
here is a version using time and frequency in the usual manner. A couple of ponts here. If you construct a grid as you do with linspace, the cosine will end up having value 0 at the first and last array point. But that does not give a satisfactory result. The array needs to fall one point short at the top end, so that if you think of a set of windows side by side, cos=0 falls into the first point of the next window up. Then the wave is periodic in the desired manner.
In the argument of cosine, the phase is a separate entity and stands apart from multiplication by any factor.
The phase here is set into the equivalent range -pi<= phi<= pi for direct comparison to the angle function.
After using fftshift (which shifts the point corresponding to f=0 to the center of the array; you need to construct a frequency array that matches), positive freq +6 corresponds to array index 507 and negative freq -6 corresponds to array index 495. Since
cos(2*pi*f0*t+phi) = (1/2)*e^(i*(2*pi*f0*t+phi)) + (1/2)*e^(-i*(2*pi*f0*t+phi))
the phase of the positive frequency peak is the same as the orignal phase and the the phase of the negative frequency peak comes in with a minus sign.
n = 1e3;
t = ((0:n-1)/n);
f0 = 6;
phi = -.3
y = cos(2*pi*f0*t + phi);
figure(1)
plot(t,y)
grid on
f = (-n/2:n/2-1);
z = fftshift(fft(y)/n);
figure(2)
stem(f,abs(z))
grid on
xlim([-20 20])
ind2 = max(find(abs(z)>.4)) % find the array index of the positive freq peak
phi2 = angle(z(ind2))
t0 = -phi2/(2*pi*f0) % since cos peak occurs when (2*pi*f0*t0 + phi) = 0
3 个评论
David Goodmanson
2024-6-28
编辑:David Goodmanson
2024-6-28
Hi Noam,
I modified the answer to show the time t0 at which the peak occurs. The example has phi = -.3.
The peak in the cos function occurs at the time when its argument = 0. So
2*pi*f0*t0 + phi = 0
and you can then determine t0. The result is t0 = .008 and if you zoom in on the plot in figure 1, that's where the peak is.
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Fourier Analysis and Filtering 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!