mpiir_l2() from user Matt L's PhD thesis to design IIR filters. I set the number of numerator and denominator coefficients both to the same value (between $100$ and $300$). The maximum pole raduis is $0.98$. My observation is that the result is effectively a FIR filter, with the denominator coefficients allways being $[1, 0, 0, ...]$. Basically, the result is equal to that of the
lslevin() function for least squares FIR design from the same source.
Clearly something in my filter spec makes the FIR solution have better errors than any IIR filter. The specification is the result of numerical optimization and looks pretty random to a humans eye, so I'm not sure what properties to look for.
Is there a way to get a true IIR filter from this function? I want to compare IIR anf FIR filters for the same specification and see what's better. Right now, I can't really compare when the results are basically the same.
EDIT: Here is a sample script and corresponding desired response. It seems like I have to use a third party hosting site for this, sorry.
close all; clear all; % loads the specification, weights and frequencies as three vectors % the sepcification has a relevant band that is padded with zeros to both % sides. There is a "don't care" gap between the relevant band and the % zeros. load 'desiredResponse.mat' % loads Ws, freqs, desiredResponse % search area allowedOrder = 1:1:3; % numerator and denominator order allowedShifts = 0:0.125:2; % linear phase term % other variables Fs = 48000; normalFreqs = freqs / Fs * 2 * pi; % normalized frequencies poleRadius = 0.98; % preallcoate for the search lowestError = inf; bestNum = ; bestDenom = ; errors = zeros(length(allowedOrder), length(allowedShifts)); orders = zeros(length(allowedOrder), length(allowedShifts)); shifts = zeros(length(allowedOrder), length(allowedShifts)); % search for best aproximation for nl = 1:length(allowedOrder) ord = allowedOrder(nl); for ns = 1:length(allowedShifts) orders(nl, ns) = ord; shifts(nl, ns) = allowedShifts(ns); % add the linear phase shift phaseShiftTerm = exp(-1i * shifts(nl, ns) * (ord-1)/2/Fs * 2 * pi * freqs); D = desiredResponse .* phaseShiftTerm; % design the filter (change denominator order here, if you like) [num,denom] = mpiir_l2(ord, ord, normalFreqs, D, Ws, poleRadius); % calculate filter response H = freqz(num, denom, freqs, Fs); % calculate the weighted error errors(nl, ns) = sum(Ws .* abs(H - D).^2); if errors(nl, ns) < lowestError lowestError = errors(nl, ns); bestNum = num; bestDenom = denom; end end end