0
$\begingroup$

I have a continuous-time signal that is of the form $$x(t)=A \sin(2 \pi f t+\phi_1)+B \sin(2 \pi (5f) t+\phi_2)+C \sin(2 \pi (7f) t+\phi_3)+n(t)$$ sampled at 50us where $f$ is the frequency I want to find out A,B and C are constants where $$0.15 \cdot A>B+C+|n(t)|$$. I don't need to know the value of any of the phases or amplitudes. I also know that the value of $f$ is 50Hz +/-15 Hz though i want to know it to at least +/- 1 Hz accuracy.

What ive tried so far:

  • A PLL-because $f$ is quite a low frequency and my frequency range is so wide and harmonics were present, the bandwidth of my loop filter was quite limited, I couldn't beat a settling time of 100ms(within 1 HZ).
  • DFT-Sampled over 100ms, the spectra was only really clear when $f$=50Hz.enter image description here
  • Measuring the time between zero crossings(Though this falls apart when the amplitudes of the harmonics and noise gets quite high)

Code used for DFT

T = 0.1; % sample time 100ms
Ts = 50e-6; % Sampling period
fs = 1/Ts; % Sampling rate
t = 0:Ts:T-Ts; % Time vector of length 100ms
f1=50;
f2=35;
f3=65;
x1 = sin(f1*2*pi*t)+0.1*sin(3*f1*2*pi*t+pi/3)+0.1*sin(5*f1*2*pi*t+pi/3)+0.01*sin(50*f1*2*pi*t+pi/3); 
x2 = sin(f2*2*pi*t)+0.1*sin(3*f2*2*pi*t+pi/3)+0.1*sin(5*f2*2*pi*t+pi/3)+0.01*sin(50*f2*2*pi*t+pi/3);
x3 = sin(f3*2*pi*t)+0.1*sin(3*f3*2*pi*t+pi/3)+0.1*sin(5*f3*2*pi*t+pi/3)+0.01*sin(50*f3*2*pi*t+pi/3);
N = length(t); % Number of elements in the time vector
F = fs/(N); % Frequency step for FFT display
f = (-fs/2):F:(fs/2)-F; % Frequency vector(same length as the time vector)

X=fftshift(fft(x1))/N; % FFT of the signal x1  (blue) 
X2=fftshift(fft(x2))/N; % FFT of the signal x2 from task 2 (red)
X3=fftshift(fft(x3))/N; % FFT of the signal x3 from task 3 (black)

figure(2);
subplot(2,1,1), plot(t, x1, 'b', t, x2, 'r', t, x3, 'k'), axis([0 .1 -1.5 1.5]);
legend('x=50Hz','x=35Hz','x=65Hz')
subplot(2,1,2), plot(f, abs(X), 'b', f, abs(X2), 'r', f, abs(X3), 'k'), axis([-100 100 0 0.5]);
|improve this question
$\endgroup$
  • $\begingroup$ Why does your DFT have so few points (from what I can see in the plot)? With a sampling interval of $T=50\mu s$ you get $2000$ samples in your time interval of $100$ms. How did you compute the DFT from these $2000$ points? $\endgroup$ – Matt L. Oct 12 '14 at 18:37
  • $\begingroup$ @MattL. I've attached the code i used to the bottom of my question. I'm really hoping you can tell me i did something terribly wrong, because i was a bit dissapointed with the result. Thanks $\endgroup$ – Darren Lanigan Oct 12 '14 at 18:46
  • $\begingroup$ I see, you zoomed in on a small range, haven't noticed that. That's why you have so few points because your sampling frequency is so high compared to the frequency you're interested in. I guess the sampling frequency is given and can't be changed? $\endgroup$ – Matt L. Oct 12 '14 at 19:06
  • $\begingroup$ @MattL. Oh yeah i zoomed right in, i shoulda mentioned that. I did that because the frequency I'm looking for is only in the range of 35-65Hz. The sampling frequency can be decreased from that value just not increased, I'm pretty new to DFTs so i just set it as high as possible as it seemed intuitive that a higher sampling rate would offer a better resolution result. $\endgroup$ – Darren Lanigan Oct 12 '14 at 19:12
  • $\begingroup$ You should probably filter your signal with a lowpass filter that leaves the fundamental frequency unchanged and attenuates the rest. In this way you will also reduce the noise. Then you can downsample your signal, zero-pad it and compute an FFT. The zero-padding will give you more points in your FFT output in the frequency range that you're interested in. $\endgroup$ – Matt L. Oct 12 '14 at 19:34
2
$\begingroup$

if you know that the 5th and 7th harmonic are precise harmonics of the fundamental at frequency $f$, you can use them along with the fundamental to help against the noise. (the signal is not "buried" in its own harmonics, the harmonics are part of the signal, unlike the noise.)

this is the "pitch detection" problem (even if you're only interested in frequencies around 50 Hz). perhaps considering computing the Average Squared Difference Function (ASDF)

$$ Q(L) = \frac{1}{M} \sum\limits_{i=0}^{M-1} \left(x[n-i] - x[n-L-i] \right)^2 $$

try computing that for lags $L$ around the period length ($\frac{1}{50}$ second). pick the $L$ such that $Q(L)$ is minimum. the noise will keep it from getting to zero, but it will get close to zero.

as to the question in the title of this thread, the shortest time to estimate the frequency is about 1 period ($\frac{1}{50}$ second) plus a little more.

|improve this answer
$\endgroup$
  • $\begingroup$ Ok so going on from your answer I would want to compute that 30 times (35-65Hz with 1Hz accuracy) and the value of Q(L) with the lowest value would be the frequency of my signal. So if i put my signal through a lowpass filter to eliminate some of the noise I should get a value of the frequency in a 1/35 seconds+filter time constant+a little more? I understand that the greater the value of "M" the greater the accuracy of my results, but how much should the "little more" be to get an accurate result. 10 samples? a quarter of the period? $\endgroup$ – Darren Lanigan Oct 13 '14 at 15:01
  • $\begingroup$ i wouldn't necessarily LPF the signal. if you know the noise is at the high frequencies, fine, LPF it. it's just that if the noise is uncorrelated to the signal (as it should be), then delaying the noise by one period should leave you with totally different noise. but delaying the signal by one period should look like the same signal. $\endgroup$ – robert bristow-johnson Oct 13 '14 at 15:15
  • $\begingroup$ so, if your sampling period is $$T=\frac{1}{f_s}= 0.000050 \text{ sec}$$ then, just to make sure we're on the same page $$ x[n] \triangleq x(nT) $$. so start with lags $L$ of about $$ L \approx f_s \cdot \frac{1}{50 \pm 15} $$ or $$ \frac{f_s}{35} < L < \frac{f_s}{65} $$. i can't tell you how many integer values of $L$ to use, perhaps all of them. evaluate $Q(L)$ for each $L$ and see which $Q(L)$ is minimum. try $M$ about the same value as the largest $L$, but it could be smaller. i dunno. $\endgroup$ – robert bristow-johnson Oct 13 '14 at 15:25
  • $\begingroup$ Sorry this should be another question, but given that there is no DC offset in my signal, what would be the disadvantage to rearanging the ASDF like so $$ Q(L) = \frac{1}{M} \sum\limits_{i=0}^{M-1} \left(x[n+i] + x[n+L-i] \right)^2 $$ $\endgroup$ – Darren Lanigan Oct 13 '14 at 16:22
  • $\begingroup$ is not your rearrangement equivalent? if this is running real time and "$n$" is the present time, if $n+L-i>n$, then you'll be looking at "future" samples. if it's not real-time or if it's looking at a buffer of samples already delayed, then it's not a problem. $\endgroup$ – robert bristow-johnson Oct 13 '14 at 18:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.