Continue to Site

Welcome to EDAboard.com

Welcome to our site! EDAboard.com is an international Electronics Discussion Forum focused on EDA software, circuits, schematics, books, theory, papers, asic, pld, 8051, DSP, Network, RF, Analog Design, PCB, Service Manuals... and a whole lot more! To participate you need to register. Registration is free. Click here to register now.

Help me solving this long-time-haunted problem...

Status
Not open for further replies.

powersys

Advanced Member level 1
Joined
Nov 29, 2005
Messages
439
Helped
3
Reputation
6
Reaction score
2
Trophy points
1,298
Activity points
4,981
This problem has been haunting me for several months. Kindly advise.

Let's say a user is given a dataset of signal as shown in Figure 1.
Assume that he doesn't know how the signal is created (i.e. he doesn't know the codes in Part-A).
He performs FFT on the signal (using Part-B code) and the results he obtains are given in Figure 2 and Table 1 (only harmonics 2-10 are shown).
Clearly, he finds out there are 4 major harmonics in the signal, i.e. 1st, 3rd, 7th, and 9th.
Although there are some little values at 4th, 8th, and 10th, he decides to assume that the magnitudes of these harmonics are zero.
He was taught that a signal can be constructed by a series of sinusoids.
Since he knows the magnitudes of major harmonics, he thinks he can reconstruct the signal using code in Part-C (assume there is no Matlab ifft function).
However, he soon finds out the reconstructed signal (as shown in Figure 3) is not similar to that in Figure 1.
Why?

Questions:
[1] Is the code performed in Part-C to reconstruct the signal correct?
[2] Are those information (especially Y(2:10) and YPHASE(2:10)) given in Table 1 useful in signal reconstruction? If YES, how can we use them correctly?

Thank you very much


Code:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-A: To create signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x = (0:1:359);
y1 = 2*sin(x*pi/180);
y3 = 1*sin(3*x*pi/180 + 30/180*pi);
y7 = 0.5*sin(7*x*pi/180);
y9 = 0.1*sin(9*x*pi/180);
y=y1+y3+y7+y9;

figure(1);
plot(x,y); xlim([0 360]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-B: To compute the FFT of the signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N=length(y);
Y=fft(y)/N;

YMAG=2*abs(Y);
YPHASE=angle(Y);
realY=real(Y);
imagY=imag(Y);

figure(2);
bar(YMAG(2:20));

[Y(2:10)' YMAG(2:10)' YPHASE(2:10)']
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-C: To reconstruct the signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rcy1 = YMAG(1+1)*sin(1*x*pi/180);
rcy3 = YMAG(1+3)*sin(3*x*pi/180);
rcy7 = YMAG(1+7)*sin(7*x*pi/180);
rcy9 = YMAG(1+9)*sin(9*x*pi/180);
rcy  = rcy1+rcy3+rcy7+rcy9;

figure(3);
plot(x,rcy);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


 

Hi
There is some mistake in the code. You did not take in account of phase term. Following code will take care of that. You can extend this concept for more harmonics.
Hope it will be helpful:).
Code:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-A: To create signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x = (0:1:359);
y1 = 2*sin(x*pi/180);
y3 = 1*sin(3*x*pi/180 + 30/180*pi);
y7 = 0.5*sin(7*x*pi/180);
y9 = 0.1*sin(9*x*pi/180);
y=y1+y3+y7+y9;

figure(1);
plot(x,y); xlim([0 360]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-B: To compute the FFT of the signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N=length(y);
Y=fft(y)/N;

Ycos=(Y(2:N/2)+Y(N:-1:N/2+2));
Ysin=j*(Y(2:N/2)-Y(N:-1:N/2+2));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PART-C: To reconstruct the signal
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rcy1 = Ycos(1)*cos(1*x*pi/180)+Ysin(1)*sin(1*x*pi/180);
rcy3 = Ycos(3)*cos(3*x*pi/180)+Ysin(3)*sin(3*x*pi/180);
rcy7 = Ycos(7)*cos(7*x*pi/180)+Ysin(7)*sin(7*x*pi/180);
rcy9 = Ycos(9)*cos(9*x*pi/180)+Ysin(9)*sin(9*x*pi/180);
rcy  = rcy1+rcy3+rcy7+rcy9;

figure(2);
plot(x,rcy);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 

    powersys

    Points: 2
    Helpful Answer Positive Rating
Hi MHanif, thanks for your reply. I tested your code and it works great! However, I don't understand some parts in the code. For examples:

[1] In part B, may I know what's the purpose of Ycos and Ysin?

[2] In part C, may I know why cosine function is used to contruct the signal instead of sine function?

[3] Is your computation based on inverse DFT formulation?

Thank you very much
 

Hi Power sys (Most probably its not your name :D ),
This formulation is based on relationship between complex exponential and cosine and sine form of Fourier series/Transform. If you know the relationships which is present in every standard text and much of it is also available in internet, then it is easy to follow.
Basically the Euler Equation is used for expanding exp(j*w*t) which is cos(w*t)+j*sin(w*t). Rearranging the coefficients will yield the coefficients of cosine and sine terms. Hope you comprehend above statements :).
Lastly you asked about the sine and cosine terms, these come always in the representation of any periodic signal. The phase term in the sine term introduces two terms i-e Asin(w*t)+Bcos(w*t).
I do not know how much is this reply helpful for you :).
 

    powersys

    Points: 2
    Helpful Answer Positive Rating
I used 'powersys' as username because at the time I registered as a member of Edaboard, I was attending a course on power system. Ha Ha...

All your replies do give me lot of inputs. However, I still don't really understand the following statement:

"The phase term in the sine term introduces two terms, i.e. Asin(w*t) + Bcost(w*t)"

I'm not able to visualize how the phase in the sine term related to those two terms. My apology for my poor background on DFT.

Thanks.
 

Hi
Let us see this simple statement.
Code:
y3 = 1*sin(3*x*pi/180 + 30/180*pi);
Here you used the term sin(ωt+θ)=sin(ωt)cos(θ)+cos(ωt)sin(θ). Letting A=cos(θ) and B=sin(θ) will make sin(ωt+θ)=Asin(ωt)+Bcos(ωt).
Hope this explanation will help you understand why to use sin(ωt) and cos(ωt) terms in general.
(By the way we are dealing with discrete signals so actually continuous time variable 't' should not be used).
 

    powersys

    Points: 2
    Helpful Answer Positive Rating
Status
Not open for further replies.

Similar threads

Part and Inventory Search

Welcome to EDABoard.com

Sponsor

Back
Top