Monday, 13 January 2014

dsp.Write a c program for FIR filter design using windows

/* Write a c program for FIR filter design using windows
 
   Input :1.Number of coefficients of the filter (M),
          2.Cutoff frequency of digital filter(wc),
          3.Choice : high pass or low pass.
          4.Choice of the window.
   
  Output :Coefficients of the corresponding filter

  Roll No :-20
  Batch :-A          */

 
 
#include<conio.h>
#include<stdio.h>
#include<math.h>

void main()
 {
   float wc,tou,M,hd[50],h[50],wn,pi,n;
   int ch,p;
   char HPF;
   clrscr();
   printf("\t\tFIR filter design using windows\n\n");
   printf("n\nEnter the length(M) of the filter(coefficient):");
   scanf("%f",&M);
   p=(int)M;
   printf("\n\nEnter the cutoff Frequency(Discrete Frequency) Wc:");
   scanf("%f"&wc);
   printf("\n\nChoice for the filter:\n");
   printf("\n\nEnter 'Y' for High Pass Filter=");
   HPF=toupper(getch());
   printf("%c"HPF);
   tou=(M-1)/2;
   pi=22.0/7.0;
   for(n=0;n<=M-1;n++)
     {
       hd[n]=(sin(wc*(n-tou)))/(pi*(n-tou));
       if(n==tou)&&((p/2)*2!=p))hd[n]=wc/pi;
         {
           if(HPF=='Y')
            {
              for(n=0;n<=M;n++)
               {
                 hd[n]=(sin(pi*(n-tou))-sin(wc*(n-tou)))/(pi*(n-tou));
                 if(n==tou)&&((p/2)*2!=p))hd[n]=1(wc/pi);
               }
            }
          printf("\n\nEnter the windows you want to use");
          printf("\n1.Rectangular window \n2.Triangular(bartlett) window\n3.Hamming Window\n4.Hanning window");
          printf("\nEnter your choice");
          scanf("%d",&ch);
          switch(ch)
           {
             case 1:for(n=0;n<=M;n++)
                     {
                       wn=1;
                       h[n]=hd[n]*wn;
                     }
                     break;
            case 2:for(n=0;n<=M-1;n++)
                     {
                       wn=1-(2*abs(n-tou)/(M-1));
                       h[n]=hd[n]*wn;
                     }
                     break;
            case 3:for(n=0;n<=M-1;n++)
                    {
                      wn=0.54=0.46*cos((2*pi*n)/(M-1));
                      h[n]=hd[n]*wn;
                    }
                     break;
           case 4:for(n=0;n<=M-1;n++)
                   {
                     wn=(1-cos((2*pi*n)/(M-1)));
                     h[n]=hd[n]*wn;
                   }
                   break;
         }
        if(HPF=='Y')
          printf("\n\nCoefficient of High Pass FIR Filter are as follows");
        else
          printf("\n\nCoefficient of Low Pass FIR Filter are as follows");
          for(n=0;n<=M-1;n++)
           {
             printf("\n\nh[%1.0f]=%f",n,h[n]);
           }
          getch();
 }

***********************************************************
                         FIR filter design using windows

 Enter the length(M) of the filter (coefficient):4

 Enter the cutoff frequency(descrete frequency) wc =50

 Choice for the filter
 Enter 'y' for High Pass Filter
 Press key 'l' for Low Pass Filter = Y

 Enter the Window you want to use
 rectangle window (Enter 1)
 triangular(bartlett) window (Enter 2)
 hamming window (Enter 3)
 hanning window (Enter 4)
 choice =1

 coefficient of lowpass FIR filter are as follows...
 h0=-0.082257
 h1=-0.084224
 h2=-0.084224
 h3=-0.082257 

dsp.C program to design Butterworth filter design using Bilinear Transformation

/* C program to design Butterworth filter design using Bilinear
   Transformation.

   Inputs: 1.Order of the butterworth filter (N)
  2.Cutoff frequency of digital filter (Wc).
 Outputs: Coefficients of digital filter(cascaded second order).
   The cascaded second order sections
 H(z)=H1(z)*H2(z)*H3(z)...
               b00+b01z^-1+b02z^-2     b10+b11z^-1+b12z^-2
   H(z)=B0-------------------- *B1---------------------*...
  a00+a01z^-1+a02z^-2     a10+a11z^-1+a12z^-2
=====================================================================*/
#include<stdio.h>
#include<conio.h>
#include<math.h>
void main()
{
 int i;
 float bi[20],b[20][3],a[20][3],p_real,p_imag;
 float wc,N,theta,OmegaC,Den;
 clrscr();
 printf("\t\t Butterworth filter design using Billinear Transformation\n");
 printf("\n\n Enter the order of the filter N= ");
 scanf("%f",&N);
 printf("\n Enter the cutoff frequency of digital filter Wc:");
 scanf("%f",&wc);
 OmegaC=tan(wc/2);
 for(i=0;i<N/2;i++)
 {
  theta=((N+2*i+1)*M_PI)/(2*N);
  p_real=-1*OmegaC*cos(theta);
  p_imag=OmegaC*sin(theta);
  Den=1+2*p_real+p_real*p_real+p_imag*p_imag;
  bi[i]=(OmegaC*OmegaC)/Den;
  b[i][0]=1;
  b[i][1]=2;
  b[i][2]=1;
  a[i][0]=1;
  a[i][1]=(2*(p_real*p_real+p_imag*p_imag)-2)/Den;
  a[i][2]=(1-2*p_real+p_real*p_real+p_imag*p_imag)/Den;
 }
 if((N/2)!=i)
 {
  i--;
  Den=OmegaC+1;
  bi[i]=OmegaC/Den;

  b[i][0]=1;
  b[i][1]=1;
  b[i][2]=0;
  a[i][0]=1;
  a[i][1]=(OmegaC-1)/Den;
  a[i][2]=0;
 }
 printf("\n The Coefficients of cascaded second order");
 printf(" \n\n sections of digital filter as follows...\n");
 for(i=0;i<N/2;i++)
 {
  printf("\n B[%d]= %f\tb[%d][0]= %f\ta[%d][0]= %f",
   i,bi[i],i,b[i][0],i,a[i][0]);
  printf("\n\t\tb[%d][1]= %f\ta[%d][1]= %f",i,b[i][1],i,a[i][1]);
  printf("\n\t\tb[%d][2]= %f\ta[%d][2]= %f",i,b[i][2],i,a[i][2]);
 }
 getch();
}

****o*p***
                Butterworth filter design using Billinear Transformation


 Enter the order of the filter N= 2

 Enter the cutoff frequency of digital filter Wc:90

 The Coefficients of cascaded second order

 sections of digital filter as follows...

 B[0]= 0.443609 b[0][0]= 1.000000       a[0][0]= 1.000000
                b[0][1]= 2.000000       a[0][1]= 0.549059
                b[0][2]= 1.000000       a[0][2]= 0.225377

dsp.Fourier Transform of the sequence and computation of transfer function. This program computes the fourier transform of the sequence x(n) and plots its magnitude and phase transfer function characteristics on the screen.

/* Fourier Transform of the sequence and computation of transfer function.
   This program computes the fourier transform of the sequence x(n) and
   plots its magnitude and phase transfer function characteristics on the
   screen.

   Inputs: 1.Number of samples of x(n)
  2.Values of samples of x(n)
  3.Frequency w at which fourier transform is to be evaluated.

   Outputs: X(w) at given value is displayed.
   Magnitude and phase plots for 0 to pi values of w are
   displayed on screen

   Assumptions: This program is for written real values of sequence x(n)
   Magnitude and phase transfer function plots are computed
   for x(n) for o to pi values of w.
==========================================================================*/

#include<stdio.h>
#include<conio.h>
#include<math.h>
#include<graphics.h>

void main()
{
 float x[100],w,xwreal,xwimag,pi,wstep;
 float static mag[640],phase[640],ymag,yphase;
 int k,N,n,gd,gm;

 clrscr();
 printf("\n Fourier Transform and Computation of transfer function.\n");
 printf("\n Enter the number samples in the Sequence x(n): ");
 scanf("%d",&N);
 printf("\n\n Enter the samples of Sequence x(n):\n");

 for(n=0;n<N;n++)
 {
  printf("\n x(%d)=",n);
  scanf("%f",&x[n]);
 }

 printf("\n Enter the Frequency W at which fourier transform is");
 printf("\n\n to be evaluated w(between 0 to pi)= ");
 scanf("%f",&w);
 xwreal=xwimag=0.0;

 for(n=0;n<N;n++)
 {
  xwreal=xwreal+x[n]*cos(w*n);
  xwimag=xwimag+x[n]*sin(w*n);
 }
 xwimag=xwimag*(-1.0);

 printf("\n The value of fourier transform is:");
 printf("\n\n x(w=%f)= %f+j(%f)",w,xwreal,xwimag);
 printf("\n\n Press any key to see magnitude and phase.....");
 getch();

 w=0.0;
 pi=22.0/7.0;
 wstep=pi/640.0;

 for(k=0;k<640;k++)
 {
  xwreal=xwimag=0.0;
  w=w+wstep;
  for(n=0;n<N;n++)
  {
   xwreal=xwreal+x[n]*cos(w*n);
   xwimag=xwimag+x[n]*sin(w*n);
  }
  xwimag=xwimag*(-1.0);
  mag[k]=sqrt(xwreal*xwreal+xwimag*xwimag);
  phase[k]=atan2(xwimag,xwreal);
 }

 detectgraph(&gd,&gm);
 initgraph(&gd,&gm,"c:\\tc\\bgi");
 setlinestyle(DOTTED_LINE,1,1);
 line(0,250,640,250);
 line(0,350,640,350);

 for(k=0;k<640;k++)
 {
  ymag=250-mag[k]*200;
  putpixel(k,ymag,YELLOW);
  yphase=350-phase[k]*50;
  putpixel(k,yphase,GREEN);
 }

 outtextxy(500,260,"Magnitude plot");
 outtextxy(500,450,"Phase plot");
 getch();
 closegraph();
}

***op***

                           Fourier Transform and Computation of transfer function.


Enter the number samples in the Sequence x (n): 1


 Enter the samples of Sequence x (n):

 x (0)=1

 Enter the Frequency w at which Fourier transform is

 to be evaluated w(between 0 to pi)= pi

 The value of Fourier transform is:

 x (w=-0.000000)= 1.000000+j(-0.000000)

 Press any key to see magnitude and phase.....









                                                                                                 

dsp.Write a C program to find DFT of a given sequences using Goertzel algorithm

/* Write a C program to find DFT of a given sequences using Goertzel algorithm. */

#include<stdio.h>
#include<conio.h>
#include<math.h>

void main()
{
int k,n,N;
float static X[100],X_Real[100],X_Imag[100];

clrscr();
printf("\tDiscrete Fourier Transform(DFT)\n");

printf("\n Enter the number samples in the sequence X(n)=");
scanf("%d",&N);
printf("Enter the number samples of sequence X(n)\n");

for(n=0;n<N;n++)
{
printf("X(%d)=",n);
scanf("%f",&X[n]);
}

for(k=0;k<N;k++)
 {
  X_Real[k] = X_Imag[k]=0.0;
    for(n=0;n<N;n++)
     {
       X_Real[k]=X_Real[k]+X[n]*cos((2*M_PI*k*(n-N))/N);
       X_Imag[k]=X_Imag[k]+X[n]*sin((2*M_PI*k*(n-N))/N);
     }
     X_Imag[k]=X_Imag[k]*(-1.0);
 }

printf("\nThe %d point DFT of given sequence is:\n",N);
printf("\n\n\tReal X(k)\t\tImaginary X(k)\n");
   for(k=0;k<N;k++)
   printf("\nX(%d)= %f\t\t\t%f\t\t",k,X_Real[k],X_Imag[k]);

getch();
}

***op**

        Discrete Fourier Transform(DFT)

 Enter the number samples in the sequence X(n)=8
Enter the number samples of sequence X(n)
X(0)=-1
X(1)=0
X(2)=2
X(3)=0
X(4)=-4
X(5)=0
X(6)=2
X(7)=0

The 8 point DFT of given sequence is:


        Real X(k)               Imaginary X(k)

X(0)= -1.000000                 -0.000000               
X(1)= 3.000000                  -0.000000               
X(2)= -9.000000                  0.000000                
X(3)= 3.000000                  -0.000000               
X(4)= -1.000000                  0.000000                
X(5)= 3.000000                  -0.000000               
X(6)= -9.000000                  0.000000                

dsp..C program to compute N-point Radix-2 DIT FFT algorithm.

/*  C program to compute N-point Radix-2 DIT FFT algorithm.

 Inputs: 1.Number of samples of x(n)
2.Values of samples of x(n)

 Output: Real and imaginary part of DFT X(k).
==================================================================*/

#include<stdio.h>
#include<conio.h>
#include<math.h>
#define SWAP(a,b)var=(a);(a)=(b);(b)=var;

void main()
{
 int N,n,m,j,k,i,p;
 float data[200],real1,imag1,real2,imag2,var;
 float costheta,sintheta,t,Theta;

 clrscr();
 printf("\n\t\t Readix-2 DIT FFT algorithm\n\n");

 printf("\n\n Enter the number of samples in the sequence x(n),N=");
 scanf("%d",&N);

 printf("\n\n Enter the Samples of the Sequence x(n):\n");
 printf("\n Real part  Imaginary part");

 for(n=1;n<=N;n++)
 {
   printf("\n x(%d)=",n-1);
   scanf("%f%f",&data[2*n-1],&data[2*n]);
 }
 n=N<<1;
 j=1;

 for(i=1;i<n;i=i+2)
 {
  if(j>i)
  {
   SWAP(data[j],data[i]);
   SWAP(data[j+1],data[i+1]);
  }
  m=n>>1;
  while(m>=2 && j>m)
  {
   j-=m;
   m>>=1;
  }
  j+=m;
 }

 k=1;m=1;t=0.0;
 while((N/(2*k))>=1)
 {
  p=pow(2,m);
  n=1;
  Theta=((2*M_PI)/(float)p)*t;
  costheta=cos(Theta);
  sintheta=sin(Theta);

  for(i=1;i<=2*N;)
  {
   real1=data[i]+costheta*data[i+p]+sintheta*data[i+1+p];
   imag1=data[i+1]+costheta*data[i+1+p]-sintheta*data[i+p];
   real2=data[i]-costheta*data[i+p]-sintheta*data[i+1+p];
   imag2=data[i+1]-costheta*data[i+1+p]+sintheta*data[i+p];

   data[i]=real1;
   data[i+1]=imag1;
   data[i+p]=real2;
   data[i+p+1]=imag2;

   if(n<k)
   {
    t=t+1;
    Theta=((2*M_PI)/(float)p)*t;
    costheta=cos(Theta);
    sintheta=sin(Theta);
   }
   else
   {
    i=i+p+2;
    n=1;
    t=0;
    Theta=((2*M_PI)/(float)p)*t;
    costheta=cos(Theta);
    sintheta=sin(Theta);
   }
  }
  k=k<<1;
  m++;
 }

 printf("\n\n Output of DIT FFt is as follows:\n");
 printf("\n\n Real part of X[k]       Imaginary part of X[k]");
 for(n=1;n<=N;n++)
 {
  printf("\n%f\t\t %f ",data[2*n-1],data[2*n]);
 }
 getch();
}

op*********


                 Readix-2 DIT FFT algorithm



 Enter the number of samples in the sequence x(n),N=8

 Enter the samples of the sequence x(n):

       Real Part   Imaginary Part

 x(0)=   0.5           0

 x(1)=   0.5           0

 x(2)=   0.5           0

 x(3)=   0.5           0

 x(4)=    0            0   

 x(5)=    0            0   

 x(6)=    0            0   

 x(7)=    0            0   


 Output of DIT FFT is as follows

       Real part of X(k)   Imaginary part of X(k)

         2.000000             0.000000

         0.500000             -1.207107  

         0.000000             0.000000

         0.500000             -0.207107

         0.000000             0.000000

         0.500000             0.207107

         0.000000             0.000000
 
         0.500000             1.207107 

dsp.CIERCULAR CONVOLUTION USING DFT AND IDFT

/*-------------------CIERCULAR CONVOLUTION USING DFT AND IDFT---------------
 This program computes the circular convolution of two causal sequences x(n)
 and h(n) using DFT and IDFT.

 Inputs: 1.Length of the sequence i.e.N
2.Samples of the two sequences to be convolved

 Outputs:  Circular convolution sequence of x(n) and h(n)

 Assumptions: The given sequences are real and causal.
============================================================================*/
#include<stdio.h>
#include<conio.h>
#include<math.h>
void main()
{
 float h[20],x[20],y[20];
 float h_real[20],h_imag[20],x_real[20],x_imag[20];
 float N,M,n,k,m,y_real[20],y_imag[20];

 clrscr();
 printf("\t\tCircular Convolution using DFT and IDFT\n\n");

 printf("\n Enter the value of N=");
 scanf("%f",&N);

 printf("\n Enter the sequence h(n):\n");
 for(n=0;n<N;n++)
 {
  printf("\n h[%1.0f]= ",n);
  scanf("%f",&h[n]);
 }

 printf("\n\n Enter the Sequence x(n):\n");
 for(n=0;n<N;n++)
 {
  printf("\n x[%1.0f]= ",n);
  scanf("%f",&x[n]);
 }

 for(k=0;k<N;k++)
 {
  h_real[k]=h_imag[k]=0.0;
  for(n=0;n<N;n++)
  {
   h_real[k]=h_real[k]+h[n]*cos((2*M_PI*k*n)/N);
   h_imag[k]=h_imag[k]+h[n]*sin((2*M_PI*k*n)/N);
  }
  h_imag[k]=h_imag[k]*(-1.0);
 }

 for(k=0;k<N;k++)
 {
  x_real[k]=x_imag[k]=0.0;
  for(n=0;n<N;n++)
  {
   x_real[k]=x_real[k]+x[n]*cos((2*M_PI*k*n)/N);
   x_imag[k]=x_imag[k]+x[n]*sin((2*M_PI*k*n)/N);
  }
  x_imag[k]=x_imag[k]*(-1.0);
 }

 for(k=0;k<N;k++)
 {
  y_real[k]=h_real[k]*x_real[k]-h_imag[k]*x_imag[k];
  y_imag[k]=h_real[k]*x_imag[k]+h_imag[k]*x_real[k];
 }

 for(n=0;n<N;n++)
 {
  y[n]=0;
  for(k=0;k<N;k++)
  {
   y[n]=y[n]+y_real[k]*cos((2*M_PI*k*n)/N)
   -y_imag[k]*sin((2*M_PI*k*n)/N);
  }
  y[n]=y[n]/N;
 }

 printf("\n The Circular of convolution is......");
 for(n=0;n<N;n++)
  printf("\n\n y[%1.0f]= %f",n,y[n]);
 getch();
}
***op***
          Circular Convolution using DFT and IDFT


 Enter the value of N=4

 Enter the sequence h(n):

 h[0]= 0

 h[1]= 1

 h[2]= 2

 h[3]= 3


 Enter the Sequence x(n):

 x[0]= 2

 x[1]= 1

 x[2]= 1

 x[3]= 2

 The Circular of convolution is......

 y[0]= 7.000000

 y[1]= 9.000000

 y[2]= 11.000000

 y[3]= 9.000000

dsp. Program for CIRCULAR CONVOLUTION of two sequences h(n) and x(n).

/* Program for CIRCULAR CONVOLUTION of two sequences h(n) and x(n).

   Inputs: 1) Length of two sequences N.
  2) Samples of two seqquences.
   Output: Circular Convolution sequence of h(n) and x(n).

   Assumptions: The given  sequence are  real and causal.
====================================================================*/

#include<stdio.h>
#include<conio.h>
#include<math.h>

void main()
{
 int N,M,n,k,m;
 float total,h[20],x[20],y[20];

 clrscr();
 printf("\n\t\t\t Circular Convolution\n\n\n");

 printf("\n Enter the value of N=");
 scanf("%d",&N);

 printf("\n\n Enter the sequence h(n):");
 for(n=0;n<N;n++)
 {
  printf("\n\n h[%d]=",n);
  scanf("%f",&h[n]);
 }

 printf("\n\n Enter the sequence x(n):");
 for(n=0;n<N;n++)
 {
  printf("\n\n x[%d]=",n);
  scanf("%f",&x[n]);
 }

 printf("\n\n\n Circular convolution is:");
 for(m=0;m<N;m++)
 {
  total=0.0;
  for(k=0;k<N;k++)
  {
   if((m-k)>=0)
    n=m-k;
   else
    n=m-k+N;
   total=total+x[k]*h[n];
  }
  y[m]=total;
  printf("\n\n y[%d]=%f",m,y[m]);
 }
 getch();
}


 
                         Circular Convolution



 Enter the value of N=4


 Enter the sequence h(n):

 h[0]=2


 h[1]=1


 h[2]=2


 h[3]=1


 Enter the sequence x(n):

 x[0]=1


 x[1]=2


 x[2]=3


 x[3]=4



 Circular convolution is:

 y[0]=14.000000

 y[1]=16.000000

 y[2]=14.000000

 y[3]=16.000000