Monday, 17 February 2014

Plot the spectrum of the real sinuoids and complex sinuoids.

Magnitude spectrum:
Magnitude spectrum of the sinusoids

Phase spectrum:

Phase spectrum (-fs/2 to +fs/2)

3D-Plots of Sinusoids: [note that Sin(-wt) should be read as Sin(wt)]
Sinusoids on Real and Imaginary axes

Complex Sinusoids 


%Code used to generate above plots

%% Clear the workspace
clear all;
close all;
clc;
home;

%% Parameters
fSamp = 1e3; % Hz
ts = 1/fSamp;
T = 1; % seconds - Data duration
t = 0:ts:T-ts;

fSig = 10; % Hz - Signal frequency
%% Real - sinusoids
% 1. Sine
s = sin(2*pi*fSig*t);

% 2. Co-sine
c = cos(2*pi*fSig*t);

% 3. Complex exponential
compExpoPositive = exp(1i*2*pi*fSig*t);
compExpoNegative = exp(-1i*2*pi*fSig*t);

%% Their Magnitude spectrum
df = fSamp / length(s);
freqIndex = -fSamp/2:df:fSamp/2 - df;

sFft = fftshift(fft(s)) ./ length(s);
% phaseSine = unwrap(angle((sFft))); %output is in radians but continuos (unwrap)
phaseSine = phase(sFft); %equal to unwrap and angle

cFft = fftshift(fft(c))./ length(s);
phaseCosine = unwrap(angle((cFft)));

fftCompPositive = fftshift(fft(compExpoPositive))./ length(s);
phaseCompPositive = phase(fftCompPositive);
fftCompNegative = fftshift(fft(compExpoNegative))./ length(s);
phaseCompNegative = phase(fftCompNegative);
% Sine
figure;
subplot(3,1,1);
stem(freqIndex,real(sFft),'b*'); hold on;
stem(freqIndex,imag(sFft),'ro');


xlim([-15 15]);
title('Sinusoid''s magnitude spectrum','color','b');
xlabel('Frequency (Hz)');
ylabel('Mag.');
legend('Real','Imag','location','best');

% Cosine
subplot(3,1,2);
stem(freqIndex,real(cFft),'g-'); hold on;
stem(freqIndex,imag(cFft),'m^');
xlim([-15 15]);
title('Co-sinusoid''s magnitude spectrum','color','b');
xlabel('Frequency (Hz)');
ylabel('Mag.');
legend('Real','Imag','location','best');

% Complex sinusoids
subplot(3,1,3);
stem(freqIndex,real(fftCompPositive),'b*'); hold on;
stem(freqIndex,imag(fftCompPositive),'ro');

stem(freqIndex,real(fftCompNegative),'g*');%,'markersize',5);
stem(freqIndex,imag(fftCompNegative),'mo');%,'markersize',5);
xlim([-15 15]);
title('Complex exponentials / Complex sinuoids','color','b');
xlabel('Frequency (Hz)');
ylabel('Mag. (arbitrary unit)');
legend('Real(exp(jwt))','Imag(exp(jwt))','Real(exp(-jwt))','Imag(exp(-jwt))','location','best');

% Setting plot properties
set(findall(gcf,'Type','line'),'linewidth',2);
set(findall(gcf,'type','axes'),'fontsize',12);
set(findall(gcf,'Type','text'),'FontSize',12);
saveas(gcf,'MagSpectrumFreqDomain.jpg');




%% Their Phase spectrum
figure;
plot(freqIndex,phaseSine*180/pi,'b');
hold on;
plot(freqIndex,phaseCosine*180/pi,'r');

plot(freqIndex,phaseCompPositive*180/pi,'g');
plot(freqIndex,phaseCompNegative*180/pi,'m');


xlim([-15 15]);
% Setting plot properties
title('Phase plot');
xlabel('Frequency (Hz)');
ylabel('Phase (degrees)');



legend('Sine','Co-sine','exp(+jwt)','exp(-jwt)','location','best');
set(findall(gcf,'Type','line'),'linewidth',2);
set(findall(gcf,'type','axes'),'fontsize',12);
set(findall(gcf,'Type','text'),'FontSize',12);

saveas(gcf,'PhaseSpecFreqDomain.jpg');




***************************************************************** 







%% 3D plot of the real and complex sinusoids

% Real sinusoids

% Complex sinusoids

% Code used to generate 3-D plots

%% Clear the workspace

clc;
clear all;
close all;

figure;

%% Plotting components of Co-Sinusoids

x=ones(1,6) .* 0.5;

y=0:0.1:0.5;

z=zeros(1,length(x));

plot3(x,y,z,'g');
hold on;

x=ones(1,6) .* -0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));

plot3(x,y,z,'g+');

%% Plotting components of Sinusoids

x=ones(1,6) .* 0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));

plot3(x,z,-y,'m');

x=ones(1,6) .* -0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));
plot3(x,z,y,'m+');

%% Plotting axis lines

xL = [-1 1];%xlim;
yL = xL;
zL = yL;
line([0 0], yL, [0 0],'linewidth',3,'color','y'); %x-axis
line(xL, [0 0], [0 0],'linewidth',3,'color','k'); %y-axis
line([0 0], [0 0],zL,'linewidth',3,'color','b'); %y-axis
hold on;
grid on;

%%
axis equal;
set(findall(gcf,'Type','line'),'linewidth',2);
legend('Cos(wt) @ w', 'Cos(wt) @ -w', 'Sin(wt) @ w','Sin(-wt) @ -w');
ylabel('Real');
xlabel('Frequency (rad)');
zlabel('Imaginary')
view(36,25);
saveas(gcf,'RealSinusoids.jpg');

%% Plotting complex exponentials

figure;

x=ones(1,6) .* +0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));
plot3(x,y,z,'r');
hold on;
x=ones(1,6) .* -0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));
plot3(x,y,z,'b');

%
x=ones(1,6) .* 0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));
plot3(x,z,-y,'r+');

x=ones(1,6) .* -0.5;
y=0:0.1:0.5;
z=zeros(1,length(x));
plot3(x,z,y,'b+');

%% Plotting lines
xL = [-1 1];%xlim;
yL = xL;
zL = yL;
line([0 0], yL, [0 0],'linewidth',3,'color','y'); %x-axis
line(xL, [0 0], [0 0],'linewidth',3,'color','k'); %y-axis
line([0 0], [0 0],zL,'linewidth',3,'color','b'); %y-axis
hold on;
grid on;
axis equal;
set(findall(gcf,'Type','line'),'linewidth',2);

legend('exp(+jwt) ','exp(-jwt)','-j*exp(+jwt)','j*exp(-jwt)');
ylabel('Real');
xlabel('Frequency (rad)');
zlabel('Imaginary')
whitebg('white');
% whitebg('blue')
view(36,25);
saveas(gcf,'ComplexExpo.jpg');    






Thursday, 6 February 2014

Euler's formula and How to plot the XY axis at the center of the plot?


%% Clearing workspace
clc;
clear all;
close all;
home;

%% Just a circle -- Euler's circle
f = 1; tm = 1/f;
numberOfSamplesPerCycle = 6; % define the sides of the polygon, if large the polygon appears as circle. If = 6, we get hexagon.

dt = tm / numberOfSamplesPerCycle;

t = 0:dt:(1/f); % the number of sides in the circle depend on the 'dt' step - this define the sample time interval

s = exp(1i*2*pi*f*t);


%% Another way to plot a circle -- Trigonometric relation with complex numbers

% We can get the same polygon using this below code:
r = abs(s);
theta = phase(s);
x = r .* cos(theta);
y = r .* sin(theta);
figure; plot(s,'g'); hold on; plot(x,y,'r*');% hold on;
legend('Euler''s style', 'polar style','Location','best');
xlabel('Real part');
ylabel('Imaginary part');
set(findall(gcf, '-property', 'FontSize'), 'FontSize', 12, 'fontWeight', 'normal')
set(findall(gcf,'Type','line'),'linewidth',2);
saveas(gcf,'fig1.jpg');
Circle plot using Euler's formula and polar form formula

%% To draw the XY axes at the centre of the plot
xL = xlim;
yL = ylim;
line([0 0], yL); %x-axis
line(xL, [0 0]); %y-axis


xlabel('Real part');
ylabel('Imaginary part');
set(findall(gcf, '-property', 'FontSize'), 'FontSize', 12, 'fontWeight', 'normal')
set(findall(gcf,'Type','line'),'linewidth',2);
saveas(gcf,'fig2.jpg');

Plotting Center XY coordinate lines

%% To draw a circle
f = 1; tm = 1/f;
numberOfSamplesPerCycle = 60; % define the sides of the polygon, if large the polygon appears as circle. If = 6, we get hexagon.

dt = tm / numberOfSamplesPerCycle;

t = 0:dt:(1/f); % the number of sides in the circle depend on the 'dt' step - this define the sample time interval
s = exp(1i*2*pi*f*t);

figure; plot(s,'r');

xL = xlim;
yL = ylim;
line([0 0], yL); %x-axis
line(xL, [0 0]); %y-axis

xlabel('Real part');
ylabel('Imaginary part');
set(findall(gcf, '-property', 'FontSize'), 'FontSize', 12, 'fontWeight', 'normal')
set(findall(gcf,'Type','line'),'linewidth',2);

saveas(gcf,'fig3.jpg');
Smooth Cirle



%% Concentric circles
r1 = 0:0.1:50;
theta1 = r1;
x1 = r1 .* cos(theta1);
y1 = r1 .* sin(theta1);

figure;plot(x1,y1,'g')
xL = xlim;
yL = ylim;
line([0 0], yL); %x-axis
line(xL, [0 0]); %y-axis

xlabel('Real part');
ylabel('Imaginary part');
set(findall(gcf, '-property', 'FontSize'), 'FontSize', 12, 'fontWeight', 'normal')
set(findall(gcf,'Type','line'),'linewidth',2);
saveas(gcf,'fig4.jpg');
Spiral plot using Euler's formula


%% Complex plane: s-plane
sigma = 0:0.1:1;
jOmega = sigma;
s = sigma + 1i * jOmega;
figure;plot(s)
xlabel('Real part');
ylabel('Imaginary part');
set(findall(gcf, '-property', 'FontSize'), 'FontSize', 12, 'fontWeight', 'normal')
set(findall(gcf,'Type','line'),'linewidth',2);

saveas(gcf,'fig5.jpg');
S-plane plot of a complex number