the Convolution – Part 2 of 2

01/11/2012 § 1 Comment

Part 1 was about ‘feeling’ convolution.

In this post I will write about implementing convolution.

So, lets take a look at convolution again: it’s the process of knowing the excitation present at any given instant, provided we know channel’s given impulse response of channel and resultant of past excitations.

From now on we’ll constrain ourselves to Z domain only.  I’ll get back with detailed explanations later, but here is the code.

If Z domain representation of 2nd order filter is of form:

H(z) = {a0 + a1*z^(-1) + a2*z^(-2)} / {1 + b1*z^(-1) + b2*z^(-2)}

a0, a1, a2, b1 and b2 are called filter coefficients. The roots of numerator gives the location of ‘zeros’ and roots of denominator gives location of ‘poles’.

Octave/Matlab script:

% 'out' is array that stores output
% 'in' is array that stores input
% size of 'out' and 'in' must be same

hz1=0.0; % stores last input excitation
hz2=0.0; % stores second last input excitation
hp1=0.0; % stores output of last excitation
hp2=0.0; % stores output of second last excitation

for i=1:length(x)

   % the convolution step
   out(i) = a0*in(i) + a1*hz1 + a2*hz2 - b1*hp1 - b2*hp2;

   % setup history
   hz2 = hz1; 
   hz1 = in(i); 
   hp2 = hp1; 
   hp1 = out(i);

   % you might want to scale your o/p here...
end

C/C++ code

// all data types are float (or double if you like, can be even an int)
// 'pOut' and 'pIn' are pointer to output and input buffers resp.
// 'iBufferSize' is size of buffer, which is int, of course
// a0, a1, a2, b1, b2 are filter coefficients as discribed in Z domain 
// equation above

float hz1=0.0; // stores last input excitation
float hz2=0.0; // stores second last input excitation
float hp1=0.0; // stores output of last excitation
float hp2=0.0; // stores output of second last excitation

for (int i = 0; i < iBufferSize; i++) 
{ 
   out = a0*pIn[i] + a1*hz1 + a2*hz2 - b1*hp1 - b2*hp2;
   hz2 = hz1;
   hz1 = pIn[i]; 
   hp2 = hp1;
   hp1 = out;

   // scaling the o/p
   if (out > 17000)
      out = 17000;
   if (out < -15000)
      out = -15000;

   pOut[i] = out;
}

I will upload complete test codes soon, which will run ‘right off the shelf’

Advertisements

Tagged: , , ,

§ One Response to the Convolution – Part 2 of 2

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

What’s this?

You are currently reading the Convolution – Part 2 of 2 at a tryst with DSP.

meta

%d bloggers like this: