167 views (last 30 days)

Hi,

I need to construct my 24 point discrete time series by using its first 5 harmonics. I know I can use x=ifft(y) to obtain the original signal in which x and y should have the same dimensions (24x1 and 24x1). My question is that is it possible to obtain the same size of x (24x1) by using let's say first 5 element of y?

Thanks,

Star Strider
on 8 Mar 2017

You could use only the first five frequencies, if you wanted an approximation. Reconstructing from less than all the available frequency components loses detail in the reconstruction. The reconstructed signal will be the size of the frequencies you used to reconstruct it. You can zero-pad the reconstructed signal out to the length of the original vector (I did not do that here), but you cannot recover the information you lost.

Example:

t = linspace(0, 1, 24); % Time Vector

x = (t >= 1/3) & (t <= 2/3); % Signal Vector

y = fft(x); % Fourier Transform

ys1 = fftshift(y); % Shift To Centre

ys2 = fliplr(fftshift(ys1(8:18))); % Take Centre 5 Frequencies, Shift Back, And Flip

xinv = ifft(ys2);

figure(1)

plot(x,'b')

hold on

plot(xinv,'r')

plot(ifft(y), '--g')

hold off

legend('Original Waveform', 'Reconstructed From First 5 Freqencies', 'Reconstructed From All Frequencies')

Star Strider
on 10 Mar 2017

That’s an interesting concept.

If you want to do the frequency domain analysis, with only 24 points, you can easily use a loop to calculate Eqns (4) and (6). I doubt that there would be anything to be gained by using the fft function, considering the data you want to recover and the length of the records only being 24 samples.

That said, this seems to be taking the long way round to avoid using a simple filter. The only advantage I can see to this approach is if the filters to isolate the various harmonics are longer than the data record. That would be the situation with a FIR filter, but might not be with an IIR filter. The paper didn’t discuss the advantages of their approach over using a discrete filter, and since I have no experience in doing load forecasting, I’ve not dealt with the practical aspects of designing filters for such data.

See if this filtering approach does what you want:

The Code —

t = linspace(0, 1, 24); % Time Vector

x = (t >= 1/3) & (t <= 2/3); % Signal Vector

Fs = 24; % Sampling Frequency (Samples/Day)

Fn = Fs/2; % Nyquist Frequency

Wp = [0.6 1.4]/Fn; % First Harmonic Passband

Ws = [0.2 1.8]/Fn; % First Harmonic Stopband

Rp = 10;

Rs = 30;

[n,Ws] = cheb2ord(Wp, Ws, Rp, Rs);

[z,p,k] = cheby2(n,Rs,Ws);

[sos1,g1] = zp2sos(z,p,k); % Filter For First Harmonic

figure(1)

freqz(sos1, 2^16, Fs)

Wp = (1+[0.6 1.4])/Fn; % Second Harmonic Passband

Ws = (1+[0.2 1.8])/Fn; % Second Harmonic Stopband

[n,Ws] = cheb2ord(Wp, Ws, Rp, Rs);

[z,p,k] = cheby2(n,Rs,Ws);

[sos2,g2] = zp2sos(z,p,k); % Filter For Second Harmonic

figure(2)

freqz(sos2, 2^16, Fs)

H1 = filtfilt(sos1,g1,+x);

H2 = filtfilt(sos2,g2,+x);

figure(3)

plot(t,x, 'LineWidth',1)

hold on

plot(t, H1)

plot(t, H2)

hold off

grid

set(gca, 'YLim',[-0.5 1.5])

legend('Original Signal', '1^{st} Harmonic', '2^{nd} Harmonic')

The Bode plots are stated in units of Hz. In reality, they are in units of ‘cycles/day’.

The ‘+’ in front of the ‘x’ in ‘+x’ converts the logical vector to a numeric vector. It would not make any difference with real data.

The Plot —

Star Strider
on 10 Mar 2017

My pleasure.

Note that all I did to design the second harmonic filter was to add 1 to the bandpass and bandstop vectors. Do the same for as many harmonics (up to 12) as you want. The Chebyshev Type II design will provide good separation with filters short enough to work with a 24-element data vector.

Opportunities for recent engineering grads.

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

Start Hunting!
## 1 Comment

## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/328801-does-anyone-know-how-to-construct-original-discrete-time-signal-by-using-ifft-x-function#comment_435286

⋮## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/328801-does-anyone-know-how-to-construct-original-discrete-time-signal-by-using-ifft-x-function#comment_435286

Sign in to comment.