escapadro
escapadro

Reputation: 31

How to generate frequency response given b,a coefficients of the system?

I have the following system, specified by the set of coefficients:

b = [1 2 3];
a = [1 .5 .25];

In the Z-Domain, such function will have the following transfer function:

 H(Z) = Y(Z)/X(Z)

So the frequency response will be just the unit circle, where:

 H(e^jw) = Y(e^jw)/X(e^jw)

Do I just substitute in the e^jw for 'Z' in my transfer function to obtain the frequency response of the system mathematically, on paper? Seems a bit ridiculous from my (a student's) point of view.

Upvotes: 3

Views: 1074

Answers (3)

AlessioX
AlessioX

Reputation: 3177

Have you tried freqz()? It returns the frequency response vector, h, and the corresponding angular frequency vector, w, for the digital filter with numerator and denominator polynomial coefficients stored in b and a, respectively.

In your case, simply follow the help:

[h,w]=freqz(b,a);

Upvotes: 1

SleuthEye
SleuthEye

Reputation: 14579

I would indeed be as simple as substituting exp(j*w) in your transfer function. There are of course different ways to implement this with Matlab. For the purpose of illustration, I will be assuming bs are the coefficients of the x sequence and as are the coefficients of the y sequence, such that the b are in the numerator and the as are in the denominator:

enter image description here

A direct evaluation with Matlab could be done with:

b = [1 2 3];
a = [1 .5 .25];

N = 513; % number of points at which to evaluate the transfer function
w = linspace(0,2*pi,N);
num = 0;
for i=1:length(b)
  num = num + b(i) * exp(-j*i*w);
end
den = 0;
for i=1:length(a)
  den = den + a(i) * exp(-j*i*w);
end
H = num ./ den;

This would be equivalent to the following which makes use of the builtin polyval:

N = 513; % number of points at which to evaluate the transfer function
w = linspace(0,2*pi,N);
H = polyval(fliplr(b),exp(-j*w))./polyval(fliplr(a),exp(-j*w));

Also, this is really evaluating the transfer function at discrete equally spaced angular frequencies w = 2*pi*k/N which corresponds to the Discrete Fourier Transform (DFT). As such it could also be done with:

N = 512;
H = fft(b,N) ./ fft(a,N);

Incidentally this is what freqz does, so you could also get the same result with:

N = 512;
H = freqz(b,a,N,'whole');

Upvotes: 0

hiandbaii
hiandbaii

Reputation: 1331

You do sub in e^jw for Z. This isn't ridiculous. Then you just sweep w from -pi to pi. Your freq response will be the absolute value of the result.

As Alessiox mentioned, freqz is the command you want to use in matlab.

Upvotes: 0

Related Questions