filter signal in time domain

Author: mohadeseh-ghafoori

DSP-HW1/DSP_HW1.m

@@ -0,0 +1,63 @@
+clear all
+close all
+clc
+a=0.95;
+teta=pi/4;
+h=[1,-2*a*cos(teta),a^2];
+Ak1=[1,0,0];
+N=3;
+
+x=audioread("C:\Users\Desktop\a2.wav");
+x=x';
+u=zeros(1,10);
+u(1)=1;
+[Bk2,r]=deconv(u,h);
+
+z1=output(x,Ak1,h,N,size(x,2));
+y1=output(z1,Ak1,Bk2,N,size(z1,2));
+e1=y1-x;
+
+sound(y1,1e4)
+
+plot(x,'r');
+hold on
+plot(y1,'b');
+plot(e1,'g');
+
+%////////////////////////////////////////////
+
+z2=output(x,Ak1,Bk2,N,size(x,2));
+y2=output(z2,Ak1,h,N,size(z2,2));
+e2=y2-x;
+
+figure
+plot(x,'r');
+hold on
+plot(y2,'b');
+plot(e2,'g');
+
+%///////////////////////////////////////////
+function y = output(x,Ak,Bk,N,M)
+y=zeros(1,M);
+for n=1:1:M
+    sum=0;
+    for m=0:1:M-1
+        if((n-m)<=0 || (n-m)>size(x,2))
+            continue;
+        elseif((m+1)>size(Bk,2))
+                break;
+        else
+        sum=sum+Bk(m+1)*x(n-m);
+        end
+    end
+    for i=1:1:N-1
+        if((n-i)<=0 || (i+1)>size(Ak,2))
+            continue;
+        else
+        sum=sum-Ak(i+1)*y(n-i);
+        end
+    end
+    y(n)=sum/Ak(1);
+end
+end
+

DSP-HW1/DSP_HW1_Report.pdf

@@ -0,0 +1,59 @@
+clc;
+clear;
+close all;
+
+b = [2,3,1,7 ];
+a = [1,0.2 ];
+
+x = randn(1000,1 );
+y = filter(b,a,x );
+y1 = filterbax(b,a,x );
+
+subplot(3,1,1);
+plot(y,'b');
+title('input signal');
+
+subplot(3,1,2);
+plot(y,'r');
+title('with filter function()');
+
+subplot(3,1,3);
+plot(y1,'g');
+title('without filter function()');
+
+
+
+function y = filterbax(b,a,x )
+% The function expects inputs to be column vectors
+
+M = length(b );
+N = length(a );
+xbuffer = [];
+ybuffer = [];
+
+for i = 1:1:length(x)
+    xbuffer = [x(i)];
+    if length(xbuffer) < M
+    index = length(xbuffer );
+    else
+    index = M ;
+    end
+
+    if length(ybuffer) < N-1
+    yindex = length(ybuffer);
+    else
+    yindex = N-1 ;
+    end
+
+    y(i) = (1/a(1)).*(sum(b(1:index).*xbuffer(1:index)) - sum(a(2:yindex+1).*ybuffer(1:yindex)));
+    ybuffer = [y(i);ybuffer ];
+end
+
+y = y(:); % Return as column vector
+end
+
+
+
+
+
+

DSP-HW1/HW1_2.m

@@ -0,0 +1,20 @@
+function output = Sigma(input1,input2,input3)
+    %intializing coefficients and input data
+    a = input1;
+    h = input2;
+    N = input3;
+    sum = 0;
+    %Processing Sigma
+    if(N>1)
+        for k=2:(N);
+            if(k<=length(a))
+                sum = sum + (a(k)*h(N-k+1));
+            end
+
+        end
+    else
+        sum=0;
+    end
+    output=sum;
+end
+

DSP-HW1/T1/Sigma.m

@@ -0,0 +1,78 @@
+%initalizing MATLAB
+    clear all
+    close all
+    clc
+%intializing Coefficients and Input Datas
+    [x Fs]=audioread('abc.wav');
+    x=x';
+    alpha = 1.1;
+    %************************
+    filter = 0 ;%if filter = 1 => filter in MATLAB , if filter = 0 => myfilter 
+    system = 1 ;%if system = 0 => h(z)=(1-alpha*e^j(tetta)*z^-1)*(1-alpha*e^-j(tetta)*z^-1)
+    theta=pi/4 ;%if system = 1 => h(z)=1/((1-alpha*e^j(tetta)*z^-1)*(1-alpha*e^-j(tetta)*z^-1))
+    %************************
+    a = [1];
+    b = [1 -2*(alpha)*cos(theta) (alpha^2)];
+%Processing 1
+    if(filter == 0 )
+        if(system == 0)
+             z=myfilter(a,b,x);
+             y=myfilter(b,a,z);
+        else
+            z=myfilter(b,a,x);
+            y=myfilter(a,b,z);
+        end
+    else
+        if(system == 0)
+             z=filter(b,a,x);
+             y=filter(a,b,z);
+        else
+            z=filter(a,b,x);
+            y=filter(b,a,z);
+        end
+    end
+
+     %sound(y,Fs);
+ %plotting
+    string_title1="\alpha = "+alpha;
+    string_title2="\theta = "+theta;
+    if(system == 0)
+        string_title3="H(z)=(1-\alphae^{+j\theta}z^{-1}).(1-\alphae^{-j\theta}z^{-1})";
+    else
+        string_title3="H(z)=1/(1-\alphae^{+j\theta}z^{-1}).(1-\alphae^{-j\theta}z^{-1})";
+    end
+    if(filter == 1 )
+        string_title4 = "filter MATLAB";
+    else
+        string_title4 = "myfilter";
+    end;
+    
+    suptitle({'','',char(string_title1),char(string_title2),char(string_title3),char(string_title4)})
+    subplot(3,2,3);
+    n = 5000:14999;
+    plot(n,x(5000:14999));
+    title('x[n]')
+    
+    subplot(3,2,4);
+    plot(n,z(5000:14999));
+    title('z[n]')
+    
+    subplot(3,2,5);
+    plot(n,y(5000:14999));
+    title('y[n]')
+    
+    subplot(3,2,6);
+    e=y-x;
+    plot(n,e(5000:14999));
+    title('e[n]')
+    
+
+
+%     ynew = filter(b,a,x);
+
+
+
+
+
+
+

DSP-HW1/T1/Untitled.m

@@ -0,0 +1,21 @@
+function output = myfilter(input1,input2,input3)
+    %intializing coefficients and input data
+    a = input1;
+    b = input2;
+    x = input3;
+    size_h=12000;
+    h=zeros(1,size_h);
+    %processing Filter
+    for i=1:size_h
+        if(i<=length(b))
+            h(i) = (b(i) - Sigma(a,h,i))/a(1);
+        else
+             h(i) = - (Sigma(a,h,i)/a(1));
+        end
+    end
+    y = conv(h,x);
+    ynew=y(1:length(x));
+    %end process
+    output = ynew;
+end
+