Program 1 (impulse & step response of a given LTI system)

clc;
close all;
N=input('enter the length of the response:');
b=input('enter the numerator co-efficients:');
a=input('enter the denominator co-efficients:');
x1=[1,zeros(1,N-1)];
x2=[ones(1,N)]
h1=filter(b,a,x1);
h2=filter(b,a,x2);
n=0:N-1;
n1=-2:(N-2)-1;
subplot(2,2,1)
stem(n,x1)
title('impulse input')
xlabel('n')
ylabel('delta n')
subplot(2,2,2)
stem(n,x2)
title('step input')
xlabel('n')
ylabel('u(n)')
subplot(2,2,3)
stem(n1,h1)
title('impulse response')
xlabel('n')
ylabel('h(n)')
subplot(2,2,4)
stem(n1,h2)
title('step response')
xlabel('n')
ylabel('h(n)')
 Inputs given

enter the length of the response: 5

Output
 Program 2 (exponential response of a given LTI system)

a)using filter

clc;
close all;
N=input('enter the length of the response:');
n=0:N-1;
x=[2*exp(j*.26*n)];
b=input('enter the numerator co-efficients:');
a=input('enter the denominator co-efficients:');
h=filter(b,a,x);
subplot(2,2,1)
stem(n,x)
title('exponential input')
xlabel('n')
ylabel('exp ip')
subplot(2,2,2)
stem(n,h)
title('exponential response')
xlabel('n')
ylabel('h(n)')

b)using convolution

clc;
close all;
N=input('enter the length of the response:');
n=0:N-1;
x=[2*exp(j*.26*n)];
h=[1 1 1];
y=conv(x,h);
stem(y)
 Inputs given
enter the length of the response: 128

enter the numerator co-efficients: [1 1 1]

enter the denominator co-efficients: 1

Output
 Program 3 (FIR filter design using pole-zero concept)

clc;
close all;
N=input('Enter the length of the response:');
b=input('enter the numerator co-efficients:');
a=1;
fs=300;
[mag,freq]=freqz(b,a,N);
freq1=freq*fs/(2*pi);
mag1=20*log(abs(mag));
plot(freq1,mag1);
Inputs given:

( for k = 1, i.e, without gain )

Enter the length of the response: 256

enter the numerator co-efficients: [1 -1 1]

Output:

40

20

-20
mag in dB

-40

-60

-80

-100
0 50 100 150
freq in Hz

Inputs given:

( for k = 2, i.e, with gain )

Enter the length of the response: 256

enter the numerator co-efficients: 2 * [1 -1 1]

Output:

40

20

-20
mag in dB

-40

-60

-80

-100
0 50 100 150
freq in Hz
 Program 4 (IIR filter design using pole-zero concept)

clc;
close all;
N=input('Enter the length of the response:');
b=input('enter the numerator co-efficients:');
a=input('enter the denominator co-efficients:');
fs=300;
[mag,freq]=freqz(b,a,N);
freq1=freq*fs/(2*pi);
mag1=20*log(abs(mag));
plot(freq1,mag1);
Inputs given:

( for k = 1, i.e, without gain )

Enter the length of the response: 256

enter the numerator co-efficients: [1 -1 1]

enter the denominator co-efficients: [.8015 -0.895 1]

Output:

10

-10

-20
mag in dB

-30

-40

-50

-60

-70
0 50 100 150
freq in Hz

Inputs given:

( for k = 2, i.e, with gain )

Enter the length of the response: 256

enter the numerator co-efficients: 2 * [1 -1 1]

enter the denominator co-efficients: [.8015 -0.895 1]

Output:

20

10

-10
mag in dB

-20

-30

-40

-50
0 50 100 150
freq in Hz
 Program 5 (To study the effect of time domain windowing)

clc;
close all;
N=input('enter the length of the data points:');
A=input('enter the value of freq shift A in kHz:');
wp=2*pi*((1+A)/8);
ws=2*pi*(8/2)/8;
p=ceil(4*pi/(ws-wp));
q=ceil(8*pi/(ws-wp))
f=[0:255]*8000/256;
f1=[0:(255+p-1)]*8000/(255+p-1);
f2=[0:(255+q-1)]*8000/(255+q-1);
n=0:N-1;
x=3*cos(2*pi*1.2*n/8)+2*cos(2*pi*(1+A)*n/8);
h1=fir1(p,(1+A)/4);
h2=fir1(q,(1+A)/4,hamming(q+1));
y=conv(x,h1);
z=conv(x,h2);
m=fft(x,256);
n=fft(y,256+p-1);
l=fft(z,256+q-1);
subplot(3,2,1)
plot(x)
title('x')
xlabel('t')
ylabel('amp')
subplot(3,2,2)
stem(f,abs(m))
title('xk')
xlabel('n')
ylabel('amp')
subplot(3,2,3)
plot(y)
title('yk')
xlabel('t')
ylabel('amp')
subplot(3,2,4)
stem(f1,abs(n))
title('yk')
xlabel('n')
ylabel('amp')
subplot(3,2,5)
plot(z)
title('z')
xlabel('t')
ylabel('amp')
subplot(3,2,6)
stem(f2,abs(l))
title('zk')
xlabel('n')
ylabel('amp')
Inputs given:
enter the length of the data points: 256
enter the value of freq shift A in kHz: 1
Output:
q =16

Inputs given:
enter the length of the data points: 256
enter the value of freq shift A in kHz: 2
Output:
q =32
 Program 6_1 (Butterworth Filter)

clc;
close all;
clear all;
Fp=input ('enter the passband edge frequency:');
Fs=input ('enter the stopband edge frequency:');
fsamp=input ('enter the sampling frequency:');
Ap=input ('enter passband attenuation:');
As=input ('enter stopband attenuation:');
ft=1;
wp=2*pi*Fp/fsamp;
ws=2*pi*Fs/fsamp;
o1=2*ft*tan(wp/2);
o2=2*ft*tan(ws/2);
[N,wc]=buttord(o1,o2,Ap,As,'s');
disp('order:')
disp(N)
wc1=wc*fsamp/(2*pi);
disp('cutoff freq:')
disp(wc1)
[b,a]=butter(N,wc,'s');
[num,den]=bilinear(b,a,ft);
freqz(num,den,256);
Inputs given:

enter the passband edge frequency:1500

Output:

order: 11
cutoff freq: 1.8604e+003

200
Magnitude (dB)

-200

-400
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)

0
Phase (degrees)

-500

-1000
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)
 Program 6_2 Kaiser window

close all;

clear all;

As =input('Enter the stopband attenuation in db:');

f1 =input('Enter the passband edge frequency in hz:');

f2 =input('Enter the stopband edge frequency in hz:');

fs =input('Enter the sampling frequency in hz:');

wp =2*pi*(f1/fs);

ws =2*pi*(f2/fs);

wc=(wp+((ws-wp)/2))/(2*pi);

N=ceil((8*pi)/(ws-wp));

k=kaiser(N+1,0.5);

figure(1);

stem(k);

b=fir1(N+1,wc);

disp('Order')

disp(N)

disp('Cutoff freq normalized:')

disp(wc)

figure(2);

freqz(b,1,256);
Input

Enter the stopband attenuation in db:30

Enter the passband edge frequency in hz:260

Enter the stopband edge frequency in hz:340

Enter the sampling frequency in hz:2000

Output

Order : 100

Cutoff freq normalized : 0.1500

0.9

0.8

0.7

0.6

0.5

0.4

0.3

0.2

0.1

0
0 20 40 60 80 100 120
 Program 6_3 Hamming

clear all;
close all;
As= input('Enter the stop band attenuation in db:');
wp= input('Enter the pass band frequency in radian:');
ws= input('Enter the stop band frequency in radian:');
wc= ((((ws-wp)/2)+wp)/(2*pi))
N= ceil((8*pi)/(ws-wp))
b= fir1(N,wc,hamming(N+1))
disp('order')
disp('N')
wc1=wc*2000
disp('Cutoff freq')
disp('wc1')
[mph,freq1]=freqz(b,1,128)
freq=freq1*(2000/(2*pi))
m=20*log(abs(mph))
subplot(2,1,1)
plot(freq,m)
xlabel('freq rad/s')
ylabel('gain in db')
m1= angle(180*mph*(2*pi))
subplot(2,1,2)
stem(freq,m1)
xlabel('freq rad/s')
ylabel('phase in degree')

Input hamming
Enter the stop band attenuation in db:50
Enter the pass band frequency in radian:0.8168
Enter the stop band frequency in radian:1.06814

wc = 0.1500
N = 100
Output hamming
 Program 7 Sampling of an image
clc;
close all;
clear all;
a=imread('C:\Program Files\MATLAB\R2007b\work\penguin.jpg');
[m n]=size(a);
p=input('enter the size you want:');
for i=1:m
for j=1:n
for k=1:p
b(i,(j-1)*p+k)=a(i,j);
end
end
end
C=b;
[m n]=size(C);
for i=1:n
for j=1:m
for k=1:p
b((j-1)*p+k,i)=C(j,i);
end
end
end
imshow(a),title('original image'),impixelinfo;
figure,imshow(b),title('zoomed image'),impixelinfo;
xlabel(sprintf('zooming factor is %g',p))
 Input given
enter the size you want: 5

Output

original image

zoomed image

zooming factor is 5
 Program 8 Comparison of DFT and DCT
clc;
close all;
clear all;
N=128;
n=0:N-1;
x=n-64;
subplot(3,1,1)
stem(n,x)
title('discrete time signal x(n)')
X=fft(x);
L=input('enter the value L=');
k=64-L:64+L;
X(k)=0;
Xmdft=ifft(X);
subplot(3,1,2)
stem(n,Xmdft)
title('Inverse DFT xmdft(n)')
for m=0:N-1;
z(m+1)=x(m+1)-Xmdft(m+1);
v=abs(z);
end
edft=sum(v)/N;
edft
W=dct(x);
L1=input('enter the value L1=');
K=127-L1:127;
W(K)=0;
xmdct=idct(W);
subplot(3,1,3)
stem(n,xmdct)
title('Inverse DCT xmdct(n)')
for P=0:N-1;
Zm(P+1)=x(P+1)-xmdct(P+1);
Vm=abs(Zm);
end
edct=sum(Vm)/N;
edct
 Inputs given

enter the value L= 2

enter the value L1= 2

Output
edft = 0.8219

edct = 1.8332e-004

discrete time signal x(n)

50

-50

-100
0 20 40 60 80 100 120 140

Inverse DFT xmdft(n)

50

-50

-100
0 20 40 60 80 100 120 140

Inverse DCT xmdct(n)

50

-50

-100
0 20 40 60 80 100 120 140
 Program 9 Sampling rate up conversion

clc;

clear all;

B=fir1(26,0.2708);

fs=8000;

N=2048;

L=3;

n=0:1:N-1;

x=5*sin(n*pi/4)+cos(5*n*pi/8);

w=zeros(1,L*N);

n=0:1:N-1;

w(L*n+1)=x(n+1);

NL=length(w);

W=2*abs(fft(w,NL))/NL;

W(1)=W(1)/2;

f=[0:1:NL/2-1]*fs*L/NL;

y=filter(B,1,w);

Y=2*abs(fft(y,NL))/NL;

Y(1)=Y(1)/2;

fSL=[0:1:NL/2-1]*fs*L/NL;

subplot(2,1,2);

plot(f,W(1:1:NL/2));

grid;

xlabel('frequency(hz)');

subplot(2,1,1);

plot(fSL,Y(1:1:NL/2));

grid;

xlabel('frequency(hz)');
Output
 Program 10 Frequency transformation
clc;

clear all;

fs=input('enter the sampling frequency:');

wd=input('enter passband edge frequency in rad/sec:');

t=1/fs;

wa=(2/t)*tan(wd*t/2);

[B,A]=lp2hp([1],[1 1 1],wa);

[b,a]=bilinear(B,A,fs)

freqz(b,a,512,fs)
Input

enter the sampling frequency:300

enter passband edge frequency in rad/sec:10

Output

b = .9833 -1.9667 0.9833

a = 1.0000 -1.9661 0.9672

100
Magnitude (dB)

-100

-200

-300
0 50 100 150
Frequency (Hz)

200
Phase (degrees)

150

100

50

0
0 50 100 150
Frequency (Hz)
 Program 11 Fft & Ifft of image
clc;

clear all;

im1=imread('C:\Program Files\MATLAB\R2007b\work\winter.jpg');

figure(1),imshow(im1),title('original image'),impixelinfo;

im2=im2bw(im1);

im3=fft2(im2);

figure(3),imshow(im3),title('fft of black &white image'),impixelinfo;

im4=ifft2(im3);

figure(4),imshow(im4),title('ifft of black &white image'),impixelinfo;

Output fft & idft

original image
fft of black &white image

ifft of black &white image

clear all;

close all;

fs= input('Enter the sampling frequency:');

fp= input('Enter the PassBand frequency:');

fst= input('Enter the StopBand frequency:');

Ap= input('Enter the Passband attenuation:');

Ast= input('Enter the StopBand attenuation:');

d=fdesign.lowpass('fp,fst,Ap,Ast',fp,fst,Ap,Ast,fs);

hd=design(d,'equiripple');

fvtool(hd);
Input

Enter the sampling frequency:10000

Enter the PassBand frequency:500

Enter the StopBand frequency:700

Enter the Passband attenuation:1

Enter the StopBand attenuation:50

Output
 Program 13 N-point DFT
#include<stdio.h>

#include<math.h>

float x[10],re[10],im[10],w;

main()

int N,k,i,n,l;

printf("enter no. of points in dft\n");

scanf("%d",&N);

printf("enter sequence length\n");

scanf("%d",&l);

printf("enter sequence\n");

for(i=0;i<N;i++)

scanf("%f",&x[i]);

if(l<N)

for(i=l;i<N;i++)

x[i]=0;

if(l>N)

printf("invalid sequence restart program\n");

goto finish;

for(k=0;k<N;k++)

re[k]=0;

im[k]=0;

for(n=0;n<N;n++)

{
w=-(2*3.1416*k*n/N);

re[k]=re[k]+x[n]*cos(w);

im[k]=im[k]+x[n]*sin(w);

printf("%f j %f\n",re[k],im[k]);

finish:;

Input

enter no. of points in dft : 4

enter sequence length :4

enter sequence :1234

Result

: re[k]= 10 -2 -2 -2

: im[k]= 0 2 0 -2
real part

Imaginary part
 Program 14 Linear convolution

#include<stdio.h>

#include<math.h>

int m,n,i,j,x[30],h[30],y[30];

main()

printf("enter the length of I seq\n");

scanf("%d",&m);

printf("enter the length of II seq\n");

scanf("%d",&n);

printf("enter first sequence of length m+n-1 by padding n-1 zeroes\n");

for(i=0;i<m+n-1;i++)

scanf("%d",&x[i]);

printf("enter second sequence of length m+n-1 by padding m-1 zeroes\n");

for(i=0;i<m+n-1;i++)

scanf("%d",&h[i]);

for(i=0;i<m+n-1;i++)

y[i]=0;

for(j=0;j<=i;j++)

y[i]+=x[j]*h[i-j];

for(i=0;i<m+n-1;i++)

printf("%d\n",y[i]);

}
Input

Enter the length of 1st seq : 4

Enter the length of 2nd seq : 4

Enter the 1st seq :1 2 3 4 0 0 0

Enter the 2nd seq :4 5 6 7 0 0 0

Result: 4 13 28 50 52 45 28
 Program 15 Impulse response

#include<stdio.h>

#define order 2

#define Len 10

float y[Len]={0,0,0},sum;

main()

int j,k;

float a[order+1]={0.1311,0.2622,0.1311};

float b[order+1]={1,-0.7478,0.2722};

for(j=0;j<Len;j++)

sum=0;

for(k=1;k<=order;k++)

if((j-k)>=0)

sum=sum+(b[k]*y[j-k]);

if(j<=order)

y[j]=a[j]-sum;

else

{
 y[j]=-sum;

printf("Response[%d] = %f\n",j,y[j]);

Output :
Response[0] = 0.131100

Response[1] = 0.360237

Response[2] = 0.364799

Response[3] = 0.174741

Response[4] = 0.031373

Response[5] = -0.024104

Response[6] = -0.026565

Response[7] = -0.013304

Response[8] = -0.002718

Response[9] = 0.001589
 Program 16 DTMF
#include<stdio.h>
#include<math.h>
#define order 2
#define len 10
float
y1[len]={0,0,0},
y2[len]={0,0,0},
y3[len]={0,0,0},
sum1,sum2,sum3;
main()
{
int j,k;
float a1[order+1]={0,0.6202,0},a2[order+1]={0,0.8131,0};
float b1[order+1]={1,-1.5687,1},b2[order+1]={1,-1.1641,1};
for(j=0;j<len;j++)
{
sum1=0;
for(k=1;k<=order;k++)
{
if((j-k)>=0)
sum1+=(b1[k]*y1[j-k]);
}
if(j<=order)
y1[j]=a1[j]-sum1;
else
y1[j]=-sum1;
printf("Tone1[%d]=%f\n",j,y1[j]);
}
for(j=0;j<len;j++)
{
sum2=0;
for(k=1;k<=order;k++)
{
if((j-k)>=0)
sum2+=(b2[k]*y2[j-k]);
}
if(j<=order)
y2[j]=a2[j]-sum2;
else
y2[j]=-sum2;
printf("Tone2[%d]=%f\n",j,y2[j]);
}
for(k=1;k<=len;k++)
{
y3[k]=y1[k]+y2[k];
}
for(j=1;j<=len;j++)
{
printf("DTMF[%d]=%f\n",j,y3[j]);
}
}
 Output
Tone1[0]=0.000000
Tone1[1]=0.620200
Tone1[2]=0.972908
Tone1[3]=0.906000
Tone1[4]=0.448335
Tone1[5]=-0.202697
Tone1[6]=-0.766306
Tone1[7]=-0.999407
Tone1[8]=-0.801464
Tone1[9]=-0.257849
Tone2[0]=0.000000
Tone2[1]=0.813100
Tone2[2]=0.946530
Tone2[3]=0.288755
Tone2[4]=-0.610390
Tone2[5]=-0.999310
Tone2[6]=-0.552907
Tone2[7]=0.355671
Tone2[8]=0.966944
Tone2[9]=0.769948
DTMF[1]=1.433300
DTMF[2]=1.919437
DTMF[3]=1.194756
DTMF[4]=-0.162055
DTMF[5]=-1.202007
DTMF[6]=-1.319213
DTMF[7]=-0.643736
DTMF[8]=0.165480
DTMF[9]=0.512099
DTMF[10]=0.000000