Using for loops to implement filter function

A basic approach is the following

1. Start by using a short x, a, and b.
2. Do the calculation of the first couple of terms by hand
3. Set up the for loop using what you learned from the hand calculations

Let's do this together to see how it is done.

 x=[1,2,3,1,2,3]
 a=[1,0.1,0.2]
 b=[2,1,3]

If you type help filter in the Matlab command window you get the equation that the filter command is using. It is:

a(1)*y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb)
                      - a(2)*y(n-1) - ... - a(na+1)*y(n-na)

The first element of our a vector is a(1)=1 so this equation becomes

y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb) - a(2)*y(n-1) - ... - a(na+1)*y(n-na)

Let's plug in our numbers to calculate the first two elements of y. Any values of x and y less than 1 are zero.

y(1) = b(1)*x(1) + b(2)*0 + ... - a(2)*0 
y(1) = b(1)*x(1)

From our test case this becomes:

y(1) = 2*1 = 2

Now the second element:

y(2) = b(1)*x(2) + b(2)*x(1) - a(2)*y(1)
y(2) =2*2 +1*1 -0.1*2 = 4.8

In order for us to create the for loop implementation in Matlab we do the following:

1. Create zero padding on the input signal.
2. Create an initial output signal of all zeros.
3. Set up the equations for the first element using a for loop for b and a.
4. Set up the equations for the second element using a for loop for b and a.
5. Set up the n for loop using the equations.

1. Using the same x vector we need to add zeros onto the front of the vector. (This is called zero padding).

x_pad = [0,0,0,x] 
      = [0,0,0,1,2,3,1,2,3]

2. Initialize the output signal.

y=zeros(size(x_pad))

3. The first real number of our input (i.e. after the zeros) is the fourth element. The first real element of the filtered signal is given by

y(4) = b(1)*x_pad(4) 

We want to set up the for loop such that for the n=4 element of the output signal the matrix element of the x_pad is 4. This means that the index of the x_pad is n-k+1 = 4-1+1=4

    n=4
    bsum=0
    for k=1:length(b)
        bsum=bsum+b(k)*x_pad(n-k+1)
    ...

4. The second real element of the filtered signal is given by

y(5) = b(1)*x_pad(5) + b(2)*x_pad(4) - a(2)*y(4) 

We want to set up the for loop such that for the n=5 element of the output signal

    n=5
    bsum=0
    for k=1:length(b)
        bsum=bsum+b(k)*x_pad(n-k+1)
    ...
    
    asum=0
    for k=2:length(b)
        asum=asum+a(k)*y(n-k+1)
    ...   

5. Using the two equations we can now set up the for loop for the signal

for n=4:length(x_pad)
    bsum=0
    for k=1:length(b)
        bsum=bsum+b(k)*x_pad(n-k+1)
    ...
    
    asum=0
    for k=2:length(b)
        asum=asum+a(k)*y(n-k+1)
    ...
    
    y(n)=bsum-asum
    ...