How to find the correct magnitude/frequency value for FFT when there has noise?

49 views (last 30 days)
Hello. I am running FFT for some research on small unknown magnitude current generated in my experiment. Those currents are non periodic and their magnitude vary everytime. I used the Matlab FFT example code to analyze the current. For the sampling frequency Fs, I used the 2* Nyquist frequency which is same as the 1/(interval between two time points). I did not try to make 2^n for the first trial run. I got some result but not as ideal as what I expected. I am also little bit confused about the the way I was doing. Can you guys help me check see what I did is correct?
Avg_freq = 1/((Time(end)-Time(1))/length(Timex));
Fs = Avg_freq;
L = length(Time);
Y = abs(fft(Current)/L);
Y1 = Y(1:L/2+1);
Y1(2:end-1) = 2*Y1(2:end-1);
f = Fs*(0:(L/2))/L;
figure(4)
plot(f,Y1,'-')
grid on
xlim([0 1*10^8])
  1 Comment
Hantao
Hantao on 13 Apr 2023
There is also one thing that really makes me confused is that since the value of f is depending on the value of Fs, if I intentionally set the Fs very high, my Y1 result will go to wrong place. What way should I defind the Fs most reasonablly?

Sign in to comment.

Accepted Answer

Paul
Paul on 14 Apr 2023
Hi Hantao,
Load the data
x = readtable('Data.xlsx');
Time = x.Time;
Current = x.Current;
Plot it
figure
plot(Time,Current),grid
It looks like the response has some pulses. Is that linear interval between -10 and -9.8 expected. Maybe somemthing like spectrogram would give more insight?
Check the time between samples. Ideally this would be constant
figure
plot(Time(1:end-1),1./diff(Time))
axis padded
Looks like the data was collected at 10e9 Hz, assuming Time in seconds, but there are a few data points that may need editing. Maybe just the time stamp is corrupted?
Avg_freq = 1/((Time(end)-Time(1))/length(Time));
Set the Avg_freq based on the apparent sampling frequency. This will only affect the scaling of the x-axis when plotting the spectrum.
%Fs = Avg_freq;
Fs = 10e9;
L = length(Time)
L = 50001
Y = abs(fft(Current)/L);
Y1 = Y(1:L/2+1);
Warning: Integer operands are required for colon operator when used as index.
With L odd, we should have
Y1 = Y(1:((L-1)/2+1));
%Y1(2:end-1) = 2*Y1(2:end-1);
Y1(2:end) = 2*Y1(2:end);
f = Fs*(0:(L-1)/2)/L;
figure
plot(f,abs(Y1),'-')
grid on
xlim([0 1*10^8])
Expand the frequency axis
figure
plot(f,abs(Y1),'-')
grid on
ylim([0 0.02])
What's the difference between these plots and what you were expecting?
  2 Comments
Hantao
Hantao on 14 Apr 2023
Hi Paul, thank you for answering.
For the data has a weird linear part, I checked and I think there was something I did incorrectly that causing that werid straight line now I have fixed it. The graph does meet what I was expecting as I derive the FFT by myself. However, there are two parts still confused me.
First, I do knew the device measuring frequency is 1 GHz but I do want to verify whether the data is using that frequency, that is why I kind of using the avg frequency interval to determine. Fun thing is that if you run my data with the code 'Avg_freq = 1/((Time(end)-Time(1))/length(Time))', it would give me 1e10 which is 10 times bigger than the device parameter but as I put Avg_freq into the fft function it actually still working same as I directly put 10^9 into the function. Do you know what kind of math mistakes did I made?
Second, when I run my old code and I did get approximately same result as yours in both magnitude and frequency. Nevertheless, when I include more data into the fft function (I added the current data for longer time), the general fft curve result still looked similar but the frequency peak started shifting a little bit (which is understandable if its small amount) and so does the magnitude, but I don't know does that level of changing is consider normal or not... (for my case I increased the size of data points from 50001 to 149999, see the atttached file) (and the magnitude peak around 11900000 HZ from 0.5 to 0.32, kinda big change)
Again, I really appreciate your time and help!
Paul
Paul on 14 Apr 2023
As to the first part, 1e10 == 10e9, so your equation for Avg_freq should give essentially the same results as above (however, the computation for Avg_freq doesn't yield exactly 10e9, if I recall correctly). Of course I have no idea why the sampling frequency in the data file is 10 times higher than what the expected sampling frequency of the device.
As to the second part, I can't offer any insight as to why adding more data changes the spectrum. However, if this new data has the same characteristic as the first with a series of pulses, you might want to isolate each pulse and compute the spectrum of each to see if the frequency content is changing with time, or try something like spectrogram to get a time-frequency plot to see if the spectral content changes with time.

Sign in to comment.

More Answers (0)

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!