Sipro Lab

Download as pdf or txt
Download as pdf or txt
You are on page 1of 16

SIPRO LAB

B.1.1 :
So the given wave is
t cos o t

a) Fourier Transform of x and plot the spectrum :-

t cos o t

o t o t
t t o

o t o t
t o t

o t o
t t o

Using complex sinusoid for imaginary axis. Thus for complex sinusoid,
o

By using this,

t t t
We get, o
t
o
t

So, magnitude of fourier transform is :-


t
t t
o
And the plot of magnitude spectrum is:-

Fig 1:- Plot of Magnitude Spectrum of |Fcc x(f)|


b)This sine wave is sampled with sampling frequency fs > 2f0. The expression of
sampled signal xs :-
t cos o t

 t cos o t

Now let us consider that,


So that, t cos o t

o t o t
=t t o

o t o
=t t o

Complex sinusoid in case of unit circle. So for complex sinusoid,


o
晦晦晦
By using this,
DFT of this function is:-
t t
晦晦晦 t 晦晦晦 t
o

So, magnitude of fourier transform is :-


t
晦晦晦 t 晦晦晦 t
o

The DFT is periodic with t and t

And the plot of magnitude spectrum of


t
o
t
晦晦晦 t 晦晦晦 t is:-

Fig 2:- plot of magnitude spectrum of |(1/fs)FdcXs(f/fs)|

c)Fourier transform of the windowed sampled signal y xsrectNt by means of


Dirichlet kernel DNt is as follows:-

Fdc as,
For discrete time signals, the convergence domain is restricted to the unit circle: the
unit constant, the complex sinusoid, the kronecker comb. So for complex sinusoid,
o
晦晦晦
By using this,
DFT of this function is:-
t t
晦晦晦 t 晦晦晦 t
o
t t
晦晦晦 t 晦晦晦 t
o
Kronecker comb is defined for integral values. So the circular convolution of this with
will consider the values dependent of integral values and disregard the cardinal
sine value as sinc is dependent on the values where the values are non-integral
t t
t t
o

Magnitude of :- t
= o t t
So the plot of magnitude spectrum of
t
is:-

Fig 3:- Plot of magnitude spectrum of |(1/fs)FdcYs(f/fs)|


d) The equation of sampled waveform YNf by means of discrete comb 1Nf,0 where Nf=
Nt and there exists an integer k0 such that λ =k/Nf & λ0=k0/Nf is:-
t t
t t
o

t t
o
t t ]

B.1.2 :
We saw above that the numerical computation to perform, for a Nt-points data
sequence (y[n])0≤n≤Nt−1 and for Nf frequencies, is, for all k ∈{0,...,Nf −1}:
t
o
y

a)If Nf=Nt, we recognize the Discrete Fourier Transform. So this transform has been
implemented by transffourier() function in the matlab code.

b)If Nf=Nt, numerical implementation uses Fast Fourier Transform algorithms.

In case of Nf>Nt i.e., zero-padding technique, we have considered nЄ[0,Nt-1] &


kЄ[0,Nf1] . In the question, exp (-j2пkn/Nf) denotes the kernel of the DFT.

Mathematically, we can find that,


o t
y e t t t t t
o t t
o t
t y e t t t t t
.
.
.
o t t
o t
t y e t t t t t
t t t
tt t t t

.
.
.
t t t t
For different values of k from 0 to Nt-1, we can observe that the summation can be
carried out. But after Nt-1, n will have no values. Now from Nt to Nf-1, we have to use
Zero padding to avoid the data loss.
For k≤Nt-1, the algorithm will work and the algorithm will also work for k>Nt-1 but
the values of will be zero.
Thus the new sequence has the same number of data points in the time domain as well
as in the frequency domain. Now after zero padding, Nt(New)=Nf.

c)For Nf<Nt, the FFT algorithm can also be used. But in this case,
o t
y e t t t t t

o t t
o t
t y e t t t t t
.
.
.
o t t
o t
t y e t t t t t
Now the rest of the values of i.e., after t cannot be taken into
account for DFT. So the signal is truncated. Thus the Time Domain Aliasing occurs
which can be understood as the effect of insufficient zero padding, which can be
thought of as under-sampling in the frequency domain. So important information may
be lost form the signal when it’s reproduced in the frequency domain.

Thus the necessary condition in case of DFT is Nf≥Nt.

B.1.3 :
Let’s have a continuous-time signal x, recorded between time 0 and T, sampled with
frequency fs. We obtain the discrete-time signal y such that:
, t
, ݁
where Nt is the rounded value of T fs towards −∞. Fourier transform of x is
approximated by:
t

With
t
o

Write a Matlab function whose input-output syntax is:


[f, tfx] = transffourier(y, Nf, fs)
Inputs are y is a vector containing the Nt values of the signal to process, Nf is the
number of frequencies, fs the sampling frequency.
Outputs are f is a vector containing the frequencies (in Hz) for which the transform
was calculated, i.e.,

tfx is the approximation of the Fourier transform of x, i.e.,


t
t

With this function, we can do plot(f, abs(tfx)) to plot the magnitude spectrum.
Matlab script for transfourier() function :-
function [f,tfx] = transfourier(y,Nf,fs)
f = (0:Nf-1)*fs/Nf
tfx=fft(y,Nf)/fs;
end

Line 1:- The function is defined here. This function is used to compute the DFT of the
given waveform. Here, inputs are y,Nf,fs and the outputs are f, tfx.
Line 2:- Here limit of f is from 0 to Nf -1 and multiplied with fs/ Nf.
Line 3:- tfx is the variable which stores the value of the approximation of the Fourier
transform of x.
B.1.4 :
Let x be a sine wave: t cos o t .
This cosine wave is sampled with frequency fs and windowed. We obtain the signal y:
, t
, ݁
Initialization
Nt = 128;
Nf = 128;
fs = 1000;
f0 = 5/128*fs;
a = 1;
phi = pi/4;
n = 0:Nt-1;
t = n/fs;
y = a*cos(2*pi*f0*t + phi);
[f, tfx] = transfourier(y, Nf, fs)

Line 1:- Nt is the variable which stores the number of sampling points i.e., 128.
Line 2:- Nf is the number of frequencies which is equal to Nt here. The value is 128.
Line 3:- fs is the sampling frequency which is defined as 1000 Hz.
Line 4:- f0 is the natural frequency of the waveform which is dependent on fs.
Line 5:- a is the amplitude of the cosine wave.
Line 6:- phi is the phase of the sine wave.
Line 7:- n is the row vector which stores the number of sampling points where
nЄ{0,…, Nt-1}.
Line 8:- t is the time-period for which the sampling is being computed.
Line 9:- The cosine wave is defined. This wave will be sampled.
Line 10:- This function is called here and it works out the FFT.

a)Matlab Script :-

figure(1)
plot (t, y)
title('Plot of sinusoidal waveform with respect to time');
xlabel('Time(Seconds)');
ylabel('Amplitude');

Line 1:- It gives the plot a certain number i.e., Figure 1, Figure 2 etc.
Line 2:- This function will generate a plot of the sampled signal w.r.t. the sampling
time period.
Line 3:- It gives the plot a title.
Line 4:- It gives the abscissa a certain name.
Line 5:- It gives the ordinate a certain name.
Fig 4 : The plot of the magnitude spectrum

b)In this case, we are taking Nf=128 which is equal to Nt.

Matlab Script:-

figure(2)
plot (f, abs(tfx))
title('Resultant signal computed by FFT');
xlabel('Frequency(Hz)');
ylabel('Magnitude spectrum');

Line 2:- Plot of the Magnitude spectrum vs Frequency of FFT


(abs because tfx may be negative or possitve but we need only absolute)
Fig 5 : Plot of the Magnitude spectrum vs Frequency of FFT

It can be seen from the figure above that the magnitude spectrum of this waveform is
like a comb. It is also noted that the peaks of the spectrum occur at X = 39.06, 960.9
and Y = 0.064 which means that the Discrete-time Fourier transform is periodic with
the period equal to sampling frequency fs. We would have obtained a perfect
kronecker comb if the signal had been given for infinite duration but such signals do
not exist in nature. As the condition is Nf= Nt, the Time Domain Aliasing occurs
which can be understood as the effect of under-sampling in the frequency domain.As
there is a sudden raise and fall in the signalwhich implies that the signal is
unstable .So important information may be lost form the signal when it’s reproduced
in the frequency domain.

c)Plot the magnitude spectrum where Nf=4096.

figure(3)
Nf = 4096; % plot of Magnitude spectrum vs Frequency of FFT
[f, tfx] = transfourier(y, Nf, fs)
plot (f, abs(tfx),'x')
title('Resultant signal computed by FFT');
xlabel('Frequency(Hz)');
ylabel('Magnitude spectrum');

Line 4 : ‘x’ is given in order to get plot in the form of x marks


Fig 6 : Plot of the Resultant signal computed by FFT

Here Nf>Nt, there are a lot of non zero peaks in the above plot so by increasing the Nf
accrate magnitude spectrum can be achived but there is a cost of computation time So
for different values of k from 0 to Nt-1, we can observe that the summation can be
carried out. But after Nt-1, n will have no values. So from Nt to Nf-1, we have to use
zero padding to avoid the data loss.
So for k≤Nt-1, the algorithm will work.
And for k>Nt-1, the algorithm will also work but the values of will be zero.
Thus the new sequence has the same number of data points in the time domain as well
as in the frequency domain. So after zero padding, Nt(New)=Nf.
In addition to that, the signal is more stable as there is not so sudden jump and fall
like in the previous case. This smooth transition from high to low value can be easily
understood by plotting it with ‘x’ marks.

d) Plot for f0=(5.7/128)fs :-

i)Where Nf=128 :-
figure(4)
Nf = 128;
f0 = 5.7/128*fs;
[f, tfx] = transfourier(y, Nf, fs)
plot (f, abs(tfx))
title('Resultant signal computed by FFT');
xlabel('Frequency(Hz)');
ylabel('Magnitude spectrum');
Fig 7 : Plot of the Resultant signal computed by FFT

ii) Where Nf=4096 :-

figure(5)
Nf = 4096;
f0 = 5.7/128*fs;
[f, tfx] = transfourier(y, Nf, fs)
plot (f, abs(tfx))
title('Resultant signal computed by FFT');
xlabel('Frequency(Hz)');
ylabel('Magnitude spectrum');

Fig 8 : Plot of the Resultant signal computed by FFT

For both the plots of Fig 7 and Fig 8 we get peaks at X = 39.06 and Y = 0.064
B.2 :
A Basic Codec:-
B.2.1 :
Preliminary Spectral Analysis:-
Matlab Script :-
Nf=2^18;
[y,srv]=audioread('E:\Sipro\Lab\musiqcue\musique.wav');
info=audioinfo('E:\Sipro\Lab\musiqcue\musique.wav');
bits=info.BitsPerSample;
sound(y,srv)
[f,tfx]=transfourier(y,Nf,srv);
figure(1)
plot(f,abs(tfx))
title('Wave form of encoder signal');
xlabel('Frequency points');
ylabel('Amplitude');
Line 1:- It is the number of frequency points.Here it’s taken as in a form of Nf=2k
where k>1 is an integer. This is used in Radix-2 FFT algorithm. It’s used here to
make the DFT very fast and computational complexity is very less.
Line 2:- It reads data from the given file and returns sampled data y and a sample rate
for that data srv.
Line 3:- It returns the information about the contents of the audio file.
Line 4:- It displays the bits per sample from the set of samples collected in the
variable info.
Line 5:- It sends the audio signal y to the speaker at sampling rate srv.
Line 6:- It calls the function transffourier() to do the FFT of the given sound wave.
Line 7:- It generates the name of the figure.
Line 8:- It plots the Magnitude Spectrum wrt frequency.
Line 9:-It gives the title of the plot.
Line 10:- It gives the name for the abscissa.
Line 11:- It gives the name for the ordinate.

Fig 9 : Plot of the Wave form of encoder signal


From the plot it’s also observed that there are two peaks: one is at 0 Hz and another
one is around 4.5x104 Hz.
B.2.2 :

Coder-Decoder

Coder Function:-

This function is defined as [npt,echelle]=codeur(y,fs,bits,fmin, fmax, fichier). Here


y,fs,bits,fmin, fmax, fichier are inputs and npt,echelle are outputs.

Codeur() defined on Matlab :-

function [npt, echelle] = codeur(y, fe, bits, fmin, fmax, fichier)

tfy = fft(y);
tfy = tfy(:);%stores as a column vector

npt = length(y);%npt stores the length of the matrix Y


kmin = round(npt*fmin/fe) + 1;%kmin and kmax are converting masque's
frequency value into matlab index
kmax = round(npt*fmax/fe) + 1;
tfymasq = tfy(kmin:kmax);%creatiing masque matrix

tfymasq = [real(tfymasq) imag(tfymasq)];% dividing the values of masque


into real and imaginary

echelle = max(max(abs(tfymasq)))*1.01;%echelle is scale and 1.01 is used


fo keeping the error between -1 & +1
tfymasq = tfymasq/echelle;

audiowrite(fichier,tfymasq, fe, 'BitsPerSample',bits)%writing a new file


by encoding the original file
end

Decoder Fuction:-

This function is defined as: - [y, fs,bits]=decodeur(fichier, fmin, fmax, npt, echelle).

Decodeur() defined in Matlab:-

function [y,srv,bits] = decodeur(fichier,fmin,fmax,npt,echelle)%


Function definition
[tfymasq,srv]=audioread(fichier);% reads the encoded file
info = audioinfo(fichier);% gets the information of file
bits = info.BitsPerSample;% gets the information
tfymasq=tfymasq(:,1)+i*tfymasq(:,2);%getting a complex waveform
tfymasq = tfymasq*echelle;
kmin = round(npt*fmin/srv)+1;
kmax = round(npt*fmax/srv)+1;
tfy = zeros(npt,1);%creating a matrix having value zero
tfy(kmin:kmax)=tfymasq;% zero padding the tfymasq
tfy=tfy(:);%all the elements of tfy in a single column
y=ifft(tfy,'symmetric');%symmetic is done to get a two sided signal
end
Line 1 :Function definition
Line 2 :reads the encoded file
Line 3 :gets the information of file
Line 4 :gets the information
Line 5 :getting a complex waveform
Line 9 :creating a matrix having value zero
Line 10 :zero padding the tfymasq
Line 11:all the elements of tfy in a single column
Line 12 : symmetic is done to get a two sided signal

This function works as reverse process of encoder function.


It goes through the sampled data from the encoded wave file.The sampled data is
stored in two different channels in this file. The first channel stores the real
component, whereas the second one stores the imaginary component of the data by
analysing the file.
The function scales back the volume using the value of echelle & combines the real
and imaginary parts of the data together from two different channels.
We can see that there are more data points in npt thant the original file so zero
padding should be used to get the real wave form therefore all the data points before
and after kmin and kmax are zero padded.
At the end, the data is transformed using IFFT. It’s noted that while doing IFFT, the
option ‘symmetric’ must be chosen to get the waveform a two sided signal.

Matlab Script where codeur() and decodeur() is implemented:-

%encorder
[y,srv]=audioread('E:\Sipro\Lab\musiqcue\musique.wav');
info = audioinfo('E:\Sipro\Lab\musiqcue\musique.wav');
bits = info.BitsPerSample;
sound(y,srv)
Nf = 2^18;
[f,tfx] = transfourier(y,Nf,fs);
figure(1);
plot(f, abs(tfx))
title('Waveform of encoded signal');
xlabel('Frequency points');
ylabel('Magnitude spectrum');
fmin=0;
fmax=8000;
[npt, echelle] = codeur(y, srv, bits, fmin,
fmax,'E:\Sipro\Lab\compressed.wav')
%decorder
[y,srv,bits]=decodeur('E:\Sipro\Lab\compressed.wav',fmin,fmax,npt,ech
elle);
[f,tfx] = transfourier(y,Nf,fs);
figure(2)
plot(f, abs(tfx))
title('Waveform of decoded signal');
xlabel('Frequency points');
ylabel('Magnitude spectrum');
sound(y,srv)

Fig 10 : Plot of the Wave form of decoded signal

Fig 11 : Plot of the Wave form of encoded signal


The above plots show that it can be observed that the values are zero in the
center part of the graph due to zero-padding at 8000 Hz. Similarly , by using
0-6000Hz frequency interval instead of 0-8000Hz, we can get a similar frequency
response in the center where there is zero padding.

We can observe that there is not much change in the audio quality when we
change the frequency from 0-6000Hz to 0-8000Hz. Here fmin and fmax should be taken
sucha way that spectrum in frequency should not go out during FFTor IFFT .
The size of the encoded file will be lesser than the original file and the size of
the file will decrease if the difference between fmax and fmin is decreased. If fmin is
taken more than 0 Hz, there will loss of amplitude near the the lower frequencies. So
we have to choose fmin=0 and human ear’s cannot hear the frequency above 8000Hz
so it is better to choose a lesser value than that.
The conclusion is that the original file is encoded into a new one which has a
lower storage space.the only thing we noticed is that there was no noticeable
difference between the original and created file.

You might also like