filter.m -- output magnitude

Dear community,

I want to impose an autocorrelative structure on a random vector, so that its autocorrelation matches the one I have in another vector with real, autocorrelated data. To that end, I estimate an AR(n) model on that data vector using arburg.m (the model order is determined on the fly through partial autocorrelations), and filter the random vector with it using filter:

vec_random_autocorrelated = filter(1,ARcoefficients',vec_random);

This works quite well except for the fact that the values in vec_random_autocorrelated are now much larger then in vec_random (on average, around 3 times higher I would say). I however want the new values to have the same average 'magnitude' as the values in vec_random. I figured out that, by changing the first argument of filter (set to 1 above), the values change into the same direction (e.g. filter(.5,ARcoefficients',vec_random); yields smaller values on average). I have checked the literature, but I am still not 100% about how filter.m actually does its job. Can anyone advise on the issue? How can I use filter.m, so that it gives me values of a desired magnitude?

Help is much appreciated.

Linda

Answer 1

[Wayne King]

It sounds to me like you are describing the need to compensate for the gain of the system at DC, 0 frequency. You should be to determine that by looking at the frequency response at DC obtained from freqz()

   A = [1 -0.9]; % allpole lowpass filter
   B = 1;
   fvtool(B,A) %note that the response at DC is 20 dB (20*log10(10))
   % you can obtain that from
   [H,W] = freqz(B,A);
   H(1)

You can make H(1) equal to unity as you describe by

    [H,W] = freqz(1/H(1),A); 
    H(1)

Is your system stable? In other words, if you look at the poles (crosses) in

zplane(1,A)

are all the poles inside the unit circle?

Answer 2

[Linda]

Hello Wayne

thank you for your comment, I found it quite useful. The system seems stable (all inside the unit circle). However, I do not get the desired results when I replace B = 1 with B = 1/H (where H = freqz(1,A)) in filter.m. The resulting values are now way smaller than in the original data set. Any ideas?

[Wayne King]

Did you replace by 1/H(1), or 1/H