%% Complex Numbers: Phasors
% Table of Contents
% Euler's Formula: Phasor Form
%     Exercise: Show that
%     Phasor Notation
%     MATLAB Complex Number Phasor Declaration
%     Phasors are stunning!
%         Applications of Phasors
%     Trigonometric Example of Phasors
%     de Moivre's Theorem
%     Polar Form Multiplication
%         The trigonometric way:
%         Phasor Multiplication
%     Dividing with Phasors
%         Exercises
% A Geometric or Phase Representation with Phasors
%         Example:
%     Important Phase Example
% Visualising Waveforms: Understanding a basic Waveform
%         The Sine Wave and Sound
%     Adding in Phase to a Sinusoid
%     MATLAB Sine Wave of a Given Ampltude and Frequency
%     MATLAB Cos v Sin Wave
% Amplitude of a Sine Wave
% Frequencies of a Sine Wave
% Phase of a Sine Wave
% Phasors: Sum of Two Sinusoids of Same Frequency
%         Example:
%         Another example
% Supporting functions
%
% We now hopefully a good basic understanding of the foundations of Complex Numbers.
% This section is a critical section for the understanding of the whole of part 1 of this module.
% As we will we can use COMPLEX numbers in a very useful, mathematically compact way to model the basic phenomenom of waveforms. We will study this with respect tp 1D Audio and then extend thois to 2D imagery. However, note that similar apporaches are used accross many subject areas with the displines of Mathematics, Engineering and Physics.
%
%
%% Euler's Formula: Phasor Form
%
% Euler's Formula states that we can express the trigonometric form as:
%
% e^{i\phi} = \cos \phi + i \sin \phi}, \:\:\phi \in \Re
% We won't prove this here. Proof here if interested.
%
%% Exercise: Show that \mathbf{e^{-i\phi} = \cos\phi - i \sin \phi}
% (Hint: Recall from CM2104 that \cos\phi is an even function and that \sin \phiis odd)
% This is also known as phasor form or a Phasor, for short.
% Note: Phasors and the related trigonometric form are very important to Fourier Theory which we study later and is integral to the first part of this module..
%
%% Phasor Notation
% The General Phasor Form is: r e^{i\phi}
% So:
% r e^{i\phi} = r(\cos\phi + i \sin\phi)
%
%% MATLAB Complex Number Phasor Declaration
% MATLAB speaks phasors natively, for example:

z = exp( 1i*(pi/4) )

[abs(z), angle(z)]  % note pi/4 =  0.7854 radians

%
%% Phasors are stunning!
% Phasers on stun!
%
% Phasors are stunning!
%% Applications of Phasors
% Can simplify Trigonometric proofs, Trig. expression manipulation etc
%            Can do Trigonometry without Trigonometry (well almost!)
%            We'll see some examples below
% Electrical Signals: Can apply simplify AC circuits to DC circuit theory (e.g. Ohm's Law)!
% Power engineering: Three phase AC power systems analysis
%              See final example for an illustration of this
% Signal Processing: Fourier Theory, Filters
%              The main focus of Phaors for this module: plenty of examples in coming sections.
%
%% Trigonometric Example of Phasors
% From Euler's Formula we can write:
% {\cos \phi = \frac{e^{i\phi} + e^{-i\phi}}{2}}
% {\sin \phi= \frac{e^{i\phi} - e^{-i\phi}}{2i}}
% Prove the above
% Proof of \cos \theta:
%
%     \mbox{We know:}\\ ~~~~e^{i\phi} = \cos\phi + i \sin \phi~~~~~~~(1)\\
% \mbox{and}\\
% ~~~~e^{-i\phi} = \cos\phi - i \sin \phi~~~~~~~(2)\\~\\
% \mbox{Adding} (1) \mbox{and} (2):\\
% ~~~~~~~e^{i\phi} + e^{-i\phi} = 2\cos \phi\\
% \mbox{Rearranging, we get:}\\
% ~~~~~~~\cos \phi = \frac{e^{i\phi} + e^{-i\phi}}{2}}
% .
% Proof of \sin \theta: As above but subtract (1) and (2).
%
%% de Moivre's Theorem
%
% Phasors provide an elegant and powerful way to achieve some proofs that would otherwise be quite tedious to evaluate!
%         This is an important concept to grasp and a prime reason why we use phasors – They also elegantly model phase, which is the key concept but that is story we will tell very soon.
%
% de Moivre's Theorem states that: If n is an integer then show that:
% {(\cos \theta + i \sin \theta)^{n} = \cos n\theta + i \sin n\theta}.
%
% Proof:
%     1. Long tedious way: Look up or recall the 	Binomial Expansion  (See CM2208 Formulae Sheet), which states:
% Expand using the terms \cos \theta and  i \sin \theta, gather terms and after a lot of simplification you shoud get the answer.
% This is not that super hard but it is long winded and challenging!
% We will no go this way!
%     2. Proof by induction: It's true for a trivial case, n =1, by Euler's formula. Assume it's true for a general case k and expand to k+1by muliplying by \cos \theta + i \sin \theta.
%         The proof is here: but there's a small amount  gathering of trigonometric terms involved.
%     3. The Phasor Way
% \mbox{From Euler's formula we know:}\\
% ~~~~~~~~~~e^{i\phi} = \cos\phi + i \sin \phi~~~~~~~~~(1)\\~\\
%
% \mbox{Now}:\\
% ~~~~~~~~~{(\cos\theta + i \sin \theta)}^{n} = {(e^{i\theta}})^{n}, ~~\mbox{and}\\
%  ~~~~~~~~~{(e^{i\theta})}^{n}=e^{in\theta}\\~\\
% \mbox{Substituting for,}~\phi = n \theta ~\mbox{in}~ (1) ~\mbox{we get}\\
%  ~~~~~~~~~e^{in\theta} = \cos n\theta + i \sin n\theta\\~~\\
% \mbox{So we get}, \\
%  ~~~~~~~~~(\cos \theta + i \sin \theta)^{n} = \cos n\theta + i \sin n\theta\\
% \mbox{Q.E.D}
% Quite elegant and efficeint!  Much quicker than 1 and less Trig knowledge (or look up) than 2.
%
%% Polar Form Multiplication
% This is also much easier and quicker with Phasors:
%% The trigonometric way:
% Let z_1 = [r_1, \phi_1] and z_2 = [r_2, \phi_2] then
%         z_1 = r_1(\cos \phi_1 + i \sin \phi_1) and
%  z_2 = r_2(\cos \phi_2 + i \sin \phi_2)
%
% Therefore:
% z_1z_2 = [r_1(\cos \phi_1 + i \sin \phi_1)]\times[r_2(\cos \phi_2 + i \sin \phi_2)]\\ ~~~= r_1r_2[(\cos \phi_1 + i \sin \phi_1)\times (\cos \phi_2 + i \sin \phi_2)]\\ ~~~= r_1r_2[ \cos \phi_1 \cos \phi_2 - \sin \phi_1\sin \phi_2 \\ ~~~~~ + i ( \cos \phi_1\sin \phi_2 + \sin \phi_1 \cos \phi_2)]
% From trigonometry we have the following relations:
% \sin (A \pm B) = \sin A \cos B \pm \cos A \sin B,
% \cos (A \pm B) = \cos A \cos B \mp \sin A \sin B,
% So finally we have:
% z_1z_2 = r_1r_2[\cos(\phi_1 + \phi_2) + i \sin(\phi_1 + \phi_2)]\\ {\bf z_1z_2 } = {\bf [r_1r_2,\phi_1+\phi_2]}
%% Phasor Multiplication
% Alternatively, we can multiply complex numbers via Phasors:
% z_1 = r_1e^{i\phi_1} and z_2 = r_2e^{i\phi_2}.
% Therefore:
% z_1z_2 = r_1e^{i\phi_1} \times r_2e^{i\phi_2}\\ = r_1r_2e^{i\phi_1}e^{i\phi_2}
% Now in general,  e^x e^y = e^{(x+y)}
% So we get:  z_1z_2 = r_1r_2e^{i(\phi_1 + \phi_2)} which (as we should expect) gives:
% \mathbf{z_1z_2  = \mathbf{[r_1r_2,\phi_1+\phi_2]}
% This is a much easier way to prove this fact — Agree?
%
%% Dividing with Phasors
% Sticking with the Phasor formulation, we can divide two complex numbers:
% z_1 = r_1e^{i\phi_1} and z_2 = r_2e^{i\phi_2}.
% Therefore:
% \frac{z_1}{z_2} = \frac{r_1e^{i\phi_1}}{r_2e^{i\phi_2}}\\ ~~= \frac{r_1}{r_2} \frac{e^{i\phi_1}}{e^{i\phi_2}}\\ ~~= \frac{r_1}{r_2} e^{i\phi_1}e^{-i\phi_2},\:\:\: \mbox{by a same argument as in multiplication}\\ ~~= \frac{r_1}{r_2} e^{i(\phi_1-\phi_2)}\\{\bf \frac{z_1}{z_2} }= {\bf [\frac{r_1}{r_2},\phi_1-\phi_2]}
% Exercise: Prove this formula via the trigonometric polar form — yawn!.
%
%% Exercises
% If z_1 = 3\sqrt{2} + 3\sqrt{2}i and z_2 = \frac{3\sqrt{3}}{2} + \frac{3}{2}i, find z_1z_2 and \frac{z_1}{z_2}, leave your answer in polar form.
% Evaluate (\frac{1}{2} + \frac{\sqrt{3}}{2}i)^3, give your answer in Cartesian form.
%
%% A Geometric or Phase Representation with Phasors
%
% Multiplying a complex number z = x + iy by i rotates the vector representing z through 90^\circ anticlockwise
%% Example:
% Let z_1 = 1.
% Then z_2 = iz_1 = i.
% Polar form of z_1 = [1,0^\circ].
% Polar form of z_2 = [1,90^\circ], Q.E.D.
%
%% Important Phase Example
%
% W.r.t  this module's topics that following (All of Digital Signal Processing and Fourier Theory): This is the most important concept we have met so far
% A phasor is a complex number used to represent a sinusoid.
% In particular, not the following duality:
%                 Sinusoid: x(t) = M\cos(\omega t + \phi), \:\: -\infty < t < \infty — a function of time
%                 Phasor: X = Me^{i\phi} = M\cos(\phi) + i M\sin(\phi) — a complex number
% So phasors and Sinusoids are related:
% \Re[Xe^{i\omega t}]  =  \Re[Me^{i\phi}e^{i\omega t}]\\ ~~~ =  \Re[Me^{i(\omega t + \phi})]\\ ~~~ =  \Re[M (\cos(\omega t + \phi) + i \sin(\omega t + \phi))\\ ~~~ =  M \cos(\omega t + \phi)\\ ~~~=  {\bf x(t)}
% The above example is importnat to grasp. In the above not that we convert the sinusoid into a Phasor, a special complex number, manipulate the Phasor – this can anything from a change or phase to more complex processing – and then take the real part of the Phasor to obtain the new, clearly real, resultant sinusoid.
%%
%
%% Visualising Waveforms: Understanding a basic Waveform
% We'll revisit this topic and study it in more depth in the next section.
% Simple Waveforms
%
% Frequency is the number of cycles per second and is measured in Hertz (Hz)
% Wavelength is  inversely proportional to frequency, i.e. Wavelength varies as \frac{1}{frequency}
%% The Sine Wave and Sound
%
% The general form of the sine wave we shall use (quite a lot of) is as follows: y = A. sin( 2\pi.n.F_w/F_s)
% where:
% A is the amplitude of the wave,
% F_w is the frequency of the wave,
% F_s is the sample frequency,
% n is the sample index.
% MATLAB function: sin() used — works in radians
% For more information look at:

% MATLAB Doc and Help

showDoc = 0;
if showDoc
    doc sin;
end

showHelp = 1;

if showHelp
    help sin;
end

%
% Example – Simple Basic sin()function and plot, run the following code:
%

figure

amp = 0.5;
freq = 2;
Fs = 256;
n = 0:floor(Fs/freq);

y = amp.*sin(2*pi*n*freq/Fs);

plot(y); % Figure above used stem() instead of plot();
title ('1 Cycle Duration Sine Wave');

% Note some convenience functions are supplied (see code at bottom, also in zip supplied) to plot sin() and other waveforms
%% Adding in Phase to a Sinusoid
%
% Now we need to add in Phase, which relates to the main adoption of the phasor representation (a clue is the name)
% Even though the Sinusoid is the simplest of waveforms, it it has many fundamental uses as we shall see more in following sections. It is critical the concepts of amplitude, frequency and, now, phase are well understood. Amplitude and frequency should be relatiuvely intiuitive but phase may need a little explaining.
%
%
% MATLAB is a great help in visualising and undertsnding these concepts
%% MATLAB Sine Wave of a Given Ampltude and Frequency
%

figure

% Natural frequency is 2*pi radians
% If sample rate is F_s HZ then 1 HZ is 2*pi/F_s
% If wave frequency is F_w then frequency is 
%       F_w* (2*pi/F_s)
% set n samples steps up to sum duration nsec*F_s where 
% nsec is the duration in seconds
% So we get y = amp*sin(2*pi*n*F_w/F_s);
amp = 0.5;

F_s = 11025;
F_w = 440;
nsec = 2;
dur= nsec*F_s;
n = 0:dur;



y = amp*sin(2*pi*n*F_w/F_s);
figure(1)
plot(y(1:500));
title('N second Duration Sine Wave');

%% MATLAB Cos v Sin Wave

figure
% Cosine is same as Sine (except 90 degrees out of phase)

yc = amp*cos(2*pi*n*F_w/F_s);

figure;
hold on;
plot(yc(1:50),'b');
plot(y(1:50),'r');
title('Cos (Blue)/Sin (Red) Plot (Note Phase Difference)');
hold off;

%
% Note some convenience functions are supplied (see code at bottom, also in zip supplied) to plot sin() and other waveforms
%
% As we will often call sin() with a need for specific parameters (amplitude, frequency, phase, duration and sample frequency)  we will also often use a convenient fuction, mysin():

figure
% mysin() example call

amp = 0.5;
freq = 2;
phase = 0; % See below for example of phase
samp_freq = 256;
dur = samp_freq/freq; % Number of Samples so on cycle again here.

y = mysin(amp,freq,phase,dur, samp_freq);

stem(y);

%
%
%% Amplitude of a Sine Wave
%
% Code for sinampdemo.m:

% Simple Sin Amplitude Demo
figure

samp_freq = 400;
dur = 800; % 2 seconds
amp = 1;
phase = 0;
freq = 1;
s1 = mysin(amp,freq,phase,dur,samp_freq);

axisx = (1:dur)*360/samp_freq; % x axis in degrees

plot(axisx,s1);
set(gca,'XTick',[0:90:axisx(end)]);

fprintf('Initial Wave: \t Amplitude = %d\n \t\t Frequency = %d\n \t\t Phase = %d\n \t\t Duration (Samples) = %d\n \t\t Sample Frequency = %d\n\n', amp, freq, phase, dur, samp_freq);

% change amplitude

amp = 2;

s2 = mysin(amp,freq,phase,dur,samp_freq);
hold on
plot(axisx, s2,'r');
set(gca,'XTick',[0:90:axisx(end)]);
hold off;

%
%% Frequencies of a Sine Wave
%
% Code for sinfreqdemo.m
%

% Simple Sin Frequency Demo

figure

samp_freq = 400;
dur = 800; % 2 seconds
amp = 1;
phase = 0;
freq = 1;
s1 = mysin(amp,freq,phase,dur,samp_freq);

axisx = (1:dur)*360/samp_freq; % x axis in degrees

plot(axisx,s1);
set(gca,'XTick',[0:90:axisx(end)]);

fprintf('Initial Wave: \t Amplitude = %d\n \t\t Frequency = %d\n \t\t Phase = %d\n \t\t Duration (Samples) = %d\n \t\t Sample Frequency = %d\n\n', amp, freq, phase, dur, samp_freq);

% change frequency

freq = 10;

s2 = mysin(amp,freq,phase,dur,samp_freq);
hold on
plot(axisx,s2,'r');

set(gca,'XTick',[0:90:axisx(end)]);
hold off;

%%
%
%% Phase of a Sine Wave
%
% The code is sinphasedemo.m:

figure

% Simple Sin Phase Demo
samp_freq = 400;
dur = 800; % 2 seconds
amp = 1;
phase = 0;
freq = 1;
s1 = mysin(amp,freq,phase,dur,samp_freq);

axisx = (0:dur-1)*360/samp_freq; % x axis in degrees

plot(axisx,s1,'b');
set(gca,'XTick',[0:90:axisx(end)]);


% change phase

phase = 135;

s2 = mysin(amp,freq,phase,dur,samp_freq);
hold on
plot(axisx,s2,'r');
set(gca,'XTick',[0:90:axisx(end)]);
hold off;

%
%% Phasors: Sum of Two Sinusoids of Same Frequency
% Back to some Phasor maths, for the last time for now.
% Hopefully we now have a good understanding and can visualise Sinusoids of different phase, amplitude and frequency.
%
% Back to Phasors:
% X = Me^{i\phi} = M\cos(\phi) + i M\sin(\phi)
%
% Consider two sinusoids: Same frequency, \omega but different phase, \thetaand \phi and amplitude, A and B
% {\bf  A\cos(\omega t + \theta)}and {\bf  B \cos(\omega t + \phi)
%  Let's add them together
%  A\cos(\omega t + \theta) + B \cos(\omega t + \phi) = \Re[Ae^{i(\omega t + \theta)} + Be^{i(\omega t + \phi)}]\\ ~~~~~~~~~~~~~~~~~~~~~~~~~~~
% = \Re[e^{i\omega t}(Ae^{i\theta} + Be^{i\phi})]
% Now let Ae^{i\theta} + B e^{i\phi} = C e^{i\gamma} for some C and \gamma, then
% \Re[e^{i\omega t}(Ae^{i\theta} + Be^{i\phi})] = \Re[e^{i\omega t}(C e^{i\gamma})]\\ ~~~~~~~~~~~~~~~~= C \cos(\omega t + \gamma)
% So we have derived the Trigonometric Equation
% \mathbf{A\cos(\omega t + \theta) + B \cos(\omega t + \phi) = C \cos(\omega t + \gamma)}
% The equivalent phasor equation is given by:
%  \mathbf{Ae^{i\theta} + B e^{i\phi} = C e^{i\gamma}
% Which version is neater? Let's see
%% Example:
%
% Simplify
% {\bf 5 \cos(\omega t + 53^\circ) + \sqrt{2} \cos(\omega t + 45^\circ)}
% Hard way via trigonometry
% Use the cosine addition formula three times
%    See maths formula sheet handout again for useful formulae
% For our third time to simplify the result.
% Not difficult but tedious once more!
%
% Easy Way Phasors
%  \Re[5e^{i53^\circ} + \sqrt{2}e^{i45^\circ}] = (3 + 4i) + (1 + i)\\ ~~~~~~~~~~~~~~~= (4 + 5i)\\
% ~~~~~~~~~~~~~~~= 6.4e^{i51^{\circ}}
% So:
% 5 \cos(\omega t + 53^\circ) + \sqrt{2} \cos(\omega t + 45^\circ) = \Re[6.4e^{i(wt + 51^{\circ})}] \\ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% = 6.4 \cos(\omega t + 51^{\circ})
% This is a very important example - make sure you understand it.
%% Another example
% Simplify {\bf \cos(\omega t + 30^\circ) + \cos(\omega t + 150^\circ) + \sin(\omega t)}
% First trick to note: \sin(\omega t) = \cos(\omega t - 90^\circ) So now simplify: {\bf \cos(\omega t + 30^\circ) + \cos(\omega t + 150^\circ) + \cos(\omega t - 90^\circ)}
% Hard way via trigonometry
% Use the cosine addition formula three times
% see maths formula sheet handout for formula
% Not difficult but tedious!
% Easy Way Phasors
% e^{i30^\circ} + e^{i150^\circ} + e^{-i90^\circ} =e^{i90} (= i)
%
% So we get:
% \Re[e^{i90}] = \cos(90^\circ)\\ ~~~~~~= 0
% So:
% {\bf \cos(\omega t + 30^\circ) + \cos(\omega t + 150^\circ) + \cos(\omega t - 90^\circ)} =0
% or
% {\bf \cos(\omega t + 30^\circ) + \cos(\omega t + 150^\circ) + \sin(\omega t)} =0
% This fact is used in three-phase AC (alternating current) to conserve current flow and balance loads.
%
%% Supporting functions
% mysin():

function s = mysin(amp,freq,phase,dur, samp_freq)
% example function to so show how amplitude,frequency and phase
% are changed in a sin function
% Inputs: amp - amplitude of the wave
%         freq - frequency of the wave
%         phase - phase of the wave in degree
%         dur - duration in number of samples
%         samp_freq - sample frequncy


x = 0:dur-1;
phase = phase*pi/180;

s = amp*sin(2*pi*x*freq/samp_freq + phase);

end