FFT and IFFT: Change the values of complex number

Hi everyone,

I am in the midst of analyzing a time series data of tidal current speed.

I have performed FFT analysis on the data and would like to change few of the complex number values at certain time steps to zero. However, upon changing the complex number to zero, when I perform an IFFT, it return a complex double values instead of double values of speed.

Any clue what causes this? What is the correct way to change complex number to 0+0i?

What I am trying to do here is to remove the tidal constituent from the current speed data and keep the non - harmonic residual speed time series.

Thanks.

Jack

Answer 1

[Star Strider]

I am not quite certain what you are doing. The fft of a time series presents data as amplitudes at each frequency. If you want to have the amplitudes at specific frequencies (not time steps) to be zero (which is what it seems that you are doing), it is necessary to filter your time series data with a ‘notch’ or ‘bandstop’ digital filter. This is not difficult to design (with the Signal Processing Toolbox). You have to have regularly-sampled data, and know what the sampling interval (‘Ts’) or sampling frequency (‘Fs’) are.

[Jack Vun Zac Lee]

My apology if my question is misleading.

First, I transformed (FFT) a time series of tidal speed. Within the output complex number (a + bi), I would like to set the amplitude (a) of certain time steps to zero. Then, I am trying to inverse transform (IFFT) the corrected time series of complex number back into the speed. However, the inverted time series contain complex double instead of the desired double values of speed.

Hence I am wondering where did I gone wrong.

Thank you.

[Star Strider]

If I understand correctly, what you are doing wrong is manipulating the fft of your data. The fft does not present data in amplitude as a function of time steps but as amplitude as a function of frequency.

If you want to filter out the twice-daily tidal signal, you have to use a digital filter for that, and filter your time series data. (It sounds counter-intuitive, but is the correct procedure.) If you attach a sample of your data (time and amplitude, preferably as a .mat file) for one or more days, I can help you design your filter. (It is not difficult, but it can be confusing if you have no experience in signal processing.)

(It is GMT-6 here — 23.00 — and the end of my day, so I will continue in the morning.)

[Jack Vun Zac Lee]

Thank you for your kind effort.

As what you said are correct, I am trying to remove the sinusoidal tidal speed from my raw time series data. In other words, I am trying to get the residual speed (non harmonic components).

I hereby attached the raw time series data (xlsx).

[Star Strider]

My pleasure!

Good morning.

This is my initial fft analysis. What peaks do you want to remove? (We can also remove the baseline easily be setting ‘CrntFT(1)=0’. It is always a real value, so it will affect the ifft only by setting the mean of the time series data to zero.)

The code:

[d,s,r] = xlsread('Jack Vun Zac Lee data.xlsx');

t_dn = datenum(s(2:end,1));
dc = find(diff(t_dn) > 1/143);  % Check For Missing Intervals

Ts = 1/144;                                 % Sampling Ts = 10 min
Fs = 1/Ts;                                  % Sampling Frequency
Fn = Fs/2;                                  % Nyquist Frequency
Crnt = d;                                   % Tidal Current
CrntFT = fft(Crnt)/size(Crnt,1);            % Fourier Transform
Fv = linspace(0, 1, size(Crnt,1)/2+1)*Fn;   % Frequency Vector
Ix = 1:length(Fv);                          % Index Vector

[pkv,loc] = findpeaks(abs(CrntFT(Ix)), 'MinPeakHeight',0.01);

figure(1)
plot(Fv, abs(CrntFT(Ix)))
grid
xlabel('Frequency (Day^{-1})')
ylabel('Amplitude')
axis([0  12    ylim])
pklbls = strsplit(sprintf('%d\n', 1:length(loc)));
text(Fv(loc), pkv, pklbls(1:end-1), 'HorizontalAlignment','center', 'VerticalAlignment','bottom')

The figure:

[Jack Vun Zac Lee]

Thank you once more.

I managed to obtain the same results as yours. I would like to remove peaks 1 to 9 then perform IFFT to produce a current speed which does not contain these peak events.

[Star Strider]

My pleasure.

You don’t need to go through the fft - ifft step with a digital filter. It does everything in the time domain. I designed a bandpass filter with a lower frequency cutoff of 5 cycles/day and an upper frequency cutoff of 24 cycles/day, so it also works to eliminate high-frequency noise. (If you change the high-frequency cutoff, keep the stopband frequency at about 1.25 times the passband frequency.)

You need the Signal Processing Toolbox for this code. It uses the ‘t_dn’, ‘Ts’, ‘Fs’, ‘Fn’, and ‘Crnt’ variables I created in my previous code, so be sure to include those:

[d,s,r] = xlsread('Jack Vun Zac Lee data.xlsx');
t_dn = datenum(s(2:end,1));
Ts = 1/144;                                 % Sampling Ts = 10 min
Fs = 1/Ts;                                  % Sampling Frequency
Fn = Fs/2;                                  % Nyquist Frequency
Crnt = d;                                   % Tidal Current

Fco = [6.0 22]/Fn;                          % Passband
Fsb = [5.0 24]/Fn;                          % Stopband
Rp  =  1;                                   % Passband Ripple
Rs  = 10;                                   % Stopband Ripple
[n,Wn]  = buttord(Fco, Fsb, Rp, Rs);        % Order of Butterworth Filter
[b,a]   = butter(n,Wn);                     % Butterworth Transfer Function
[sos,g] = tf2sos(b,a);                      % ‘SOS’ For Stability

figure(3)                                   % Check Filter Performance (Optional)
freqz(sos, 512, Fs)

CrntFilt = filtfilt(sos,g,Crnt);            % Filter ‘Current’ (Time Domain)

figure(2)
plot(t_dn, Crnt, '-b')
hold on
plot(t_dn, CrntFilt, '-r')
hold off
grid
datetick('x', 'dd:mm HH.MM')
xlabel('Time (min)', 'FontSize',10, 'FontWeight','bold')
ylabel('Current (units)', 'FontSize',10, 'FontWeight','bold')
legend('Raw Data', 'Filtered', 'Location','NE')
set(gca, 'FontSize',8)
axis([min(t_dn)  max(t_dn)    ylim])

The figure(3) freqz call lets you see how the filter works. I used it to be sure the filter was stable and does what I want it to , but it is not necessary for the code.

[Jack Vun Zac Lee]

Thank you for the hard works. I will give it a try!

[Star Strider]

My pleasure!

Answer 2

[Honglei Chen]

Most likely you didn't maintain the symmetry when you setting this to zero. It would be helpful if you show the snippet where you set the spectrum to 0.

[Jack Vun Zac Lee]

Thank you for the feedback. What do you mean by setting the symmetry to zero?