1
$\begingroup$

I'm currently at a loss for why my code seems to not be working. I'm doing a sine sweep on a filter to calculate its freq and phase response. The frequency response checks out, however the phase response doesn't look right when compared to MATLABs freqz function. Plus since I know the filter is a comb, I can tell its wrong as well.

Here's what I've got so far.

%   Filter Coefficients
b = [1];
a = [1, 0, 0, 0, 0, 0, 0, 0, 0, -0.96];

%  sweep is my sine wave that goes from 0->fs/2
filtered_sweep = filter(b, a, sweep);

% This section correctly Grabs and plots the Freq Response
fft_sweep = fft(filtered_sweep);
freq_response = abs(fft_sweep);
freq_response = 20*log10(freq_response/1);      % Convert to dB scale
freq_response = freq_response(1:end/2);         % Cut to pi fs/2    
figure;
plot(freq_response);

And Now here is where I try and get the Phase two different ways with the same wrong result. From what I understand to get the phase you just need the angle of the complex values. In MATLAB I thought angle() should work, but just to be safe I tried to manually compute the angle myself using atan with the imag and real components.

% Pre Allocate Vector
% Manual Angle Extraction without using Angle()
phase = zeros(1, length(fft_sweep));
for i=1:length(fft_sweep)
    phase(i) = atan2(imag(fft_sweep(i)), real(fft_sweep(i)));
end

% Using angle()
phase_response = angle(fft_sweep);
phase_response = phase_response(1:end/2);
figure;
plot(unwrap(phase_response));
grid on;

Both methods give the same incorrect phase response. Any help would be appreciated.

$\endgroup$
  • $\begingroup$ You're determining the phase of the system plus the phase of the sweep. $\endgroup$ – Matt L. Feb 26 '15 at 17:44
2
$\begingroup$

You missing one, very important step. In order to get the impulse response of your filter while using the sweep sine, you must get the output of your system $h[n]$ to excitation with sweep-sine signal $x[n]$:

$y[n]=h[n]\star x[n]$

In your case that's the filtering with filter function. Now in order to get the impulse response you must convolve your signal with the inverse filter $f[n]$. Generally that's the time-inverted sweep signal (although special kind of pre-conditioning, amplitude scaling, etc. might be done as well):

$h[n]=y[n]\star f[n]$

Now you have your impulse response that you can use for further analysis. Please take a look at this answer tackling very similar problem.

For more theory please see:

A. Farina - Advancements in Impulse Response Measurements by Sine Sweeps

A. Farina - Simultaneous Measurement of Impulse Response and Distortion with a Swept-Sine Technique

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.