[signal, freq]=wavread('orig1.wav'); %read orig1.wav file and return 2 parameters
n=length(signal);

figure(1);
plot((1:n)/freq,signal);
title('original signal in the time domain');
xlabel('second');
grid;

figure(2);
signal10=signal-mean(signal);
fsignal=fft(signal10);
plot((1:n/2)/n*freq,abs(fsignal(1:n/2)));

title('original signal in frequency domain');
xlabel('Hz');
ylabel('magnitude');
grid;

figure(3);
plot((1:n/2)/n*freq,unwrap(angle(fsignal(1:n/2))));
title('original signal in frequency domain');
xlabel('Hz');
ylabel('phase');
grid;

sound(signal,freq); %playing sound

Plotting the noisy signal for time , frequency (magnitude and phase) domain.

[signal, freq]=wavread('noisy1.wav'); %read noisy.wav file and return 2 parameters
n=length(signal);

figure(4);
plot((1:n)/freq,signal);
title('noisy signal in the time domain');
xlabel('second');
grid;

figure(5);
signal10=signal-mean(signal);
fsignal=fft(signal10);
plot((1:n/2)/n*freq,abs(fsignal(1:n/2)));

title('noisy signal in frequency domain');
xlabel('Hz');
ylabel('magnitude');
grid;

figure(6);
plot((1:n/2)/n*freq,unwrap(angle(fsignal(1:n/2))));
title('noisy signal in frequency domain');
xlabel('Hz');
ylabel('phase');
grid;

sound(signal,freq); %playing sound

Filterring process for experiment Part 1 using the following function ( butter , buttord , freqz , filter ) , the calculation of MSE and SNR values and the rest of the procedure in part 2

[signal2, freq]=wavread('orig1.wav');
[signal1, freq]=wavread('noisy1.wav');
a=signal1;

fs=3417;
fp=243;
ft=8000;
ws=2*fs/ft; %calculate the normalized stopband frequency
wp=2*fp/ft;  %calculate the normalized passband frequency

[N,Wn]=buttord(wp,ws,0.03,40);
[B,A]=butter(N,Wn);

figure(7);
stem(b); %stem the impluse response
title('impluse response');
xlabel('samples');

figure(8);
freqz(B,A,N);
y=filter(B,A,a);
M=length(y);

sound(y,freq); %play the filtered signal sound
wavwrite(y,'parttwo'); %save the filtered signal as parttwo.wav file

MSE = sum((y(1:M)-signal2(1:M)).^2)/M %calculation of MSE(mean squared error)
NSR = 10*log((sum(signal2(1:M).^2)/M)/MSE) %calculation of SNR(signal to noise ratio)
%calculation for SNR for the noisy signal before filtering
MSE = sum((signal1(1:M)-signal2(1:M)).^2)/M ;
NSR = 10*log((sum(signal2(1:M).^2)/M)/MSE)

Plotting the filtered signal for time , frequency (magnitude and phase) domain.

[signal, freq]=wavread('parttwo.wav'); %read the parttwo.wav file which is the filtered signal
n=length(signal);

figure(9);
plot((1:n)/freq,signal);
title('filtered signal in the time domain');
xlabel('second');
grid;

figure(10);
signal10=signal-mean(signal);
fsignal=fft(signal10);
plot((1:n/2)/n*freq,abs(fsignal(1:n/2)));

title('filtered signal in frequency domain');
xlabel('Hz');
ylabel('magnitude');
grid;

figure(11);
plot((1:n/2)/n*freq,unwrap(angle(fsignal(1:n/2))));
title('filtered signal in frequency domain');
xlabel('Hz');
ylabel('phase');
grid;