This is my code, the optimization is not good, I hope to find a good way to optimize the coefficients of this IIR filter? I can get good results using the yulewlaker equation,

clc,clear;

[freq,H]=freqAndHFunction();

freq=[0,freq];

magnitude=[0,H];

% 采样率

Fs= 16000;

% 进行插值

interp_freq = linspace(0, 8000, 160); % 选择插值点数

interp_magnitude = interp1(freq, magnitude, interp_freq, 'linear');

% PSO选项

options = optimoptions('particleswarm', 'SwarmSize', 150, 'MaxIterations', 2000);

% 参数界限

lb = -1e-7 * ones(1, 162); % 下界

ub = 1e-7 * ones(1, 162); % 上界

% 运行PSO

best_params = particleswarm(@objective_function, 162, lb, ub, options);

b_estimated = best_params(1:41);

a_estimated = best_params(42:82);

% 计算估计滤波器的频率响应

[H_est,~] = freqz(b_estimated, a_estimated, interp_freq, 16000);

% 可视化结果

figure;

plot(interp_freq, 20*log10(abs(H_est)), 'r'); % 估计的响应

hold on;

plot(interp_freq, 20*log10(interp_magnitude), 'b'); %原始响应

legend('Estimated Response', 'Original Response');

xlabel('Frequency (Hz)');

ylabel('Magnitude (dB)');

function error = objective_function(params)

try

b = params(1:41); % 分子系数

a = params(42:82); % 分母系数

% 计算当前IIR滤波器的频率响应

[H,~] = freqz(b, a, interp_freq, 16000);

% 计算误差

error = sum((abs(H) - interp_magnitude).^2);

% 检查滤波器的稳定性

if ~isstable(b, a)

error = inf; % 给不稳定的滤波器一个大的错误值

end

catch

error = inf; % 对于任何其他错误，返回一个大的错误值

end

end

function [f, H] = freqAndHFunction()

% Define freq array

f = [1500, 1600, 1700, 1800, 1900, 2000, 2100, 2200, 2300, 2400, 2500, 2600, 2700, 2800, 2900, 3000, 3100, 3200, 3300, 3400, 3500, 3600, 3700, 3800, 3900, 4000, 4100, 4200, 4300, 4400, 4500, 4600, 4700, 4800, 4900, 5000, 5100, 5200, 5300, 5400, 5500, 5600, 5700, 5800, 5900, 6000, 6100, 6200, 6300, 6400, 6500, 6600, 6700, 6800, 6900, 7000, 7100, 7200, 7300, 7400, 7500, 7600, 7700, 7800, 7900, 8000];

% Define tvr array

tvr = [128.5163628, 129.051837, 128.0333238, 127.200628, 128.4267629, 129.4437987, 130.9112373, 133.5186641, 135.2522284, 134.405548, 130.9848593, 128.2534974, 134.0918335, 138.610268, 140.9095005, 142.6631096, 143.7293213, 144.7608089, 145.3854305, 145.5075137, 145.554636, 144.7531332, 143.9672297, 143.410465, 142.2522647, 141.169289, 140.4608309, 139.6028821, 138.8584932, 138.4957287, 138.1202031, 137.6812736, 137.4824973, 137.3362914, 137.3640559, 137.3579424, 137.6974913, 138.0771592, 138.0983559, 137.5053153, 136.0987116, 133.593733, 131.4955087, 130.7312329, 131.25016, 131.7451198, 132.140601, 132.1775175, 131.4427206, 130.7975073, 129.799245, 128.9741335, 129.4695932, 130.020389, 129.89581, 129.3831058, 127.7675183, 124.4633058, 118.3829614, 108.3285636, 113.1397294, 115.693913, 115.7284724, 114.169323, 116.192023, 120.3573138];

% Calculate H array

H = 10.^((tvr-120)/20);

end

Answer 1

[Balaji]

Hi 飞龙 杜,

I understand that you are seeking for ways to optimize the filter coefficients of an IIR filter.

Using Particle Swarm Optimization you can improve the results by changing the lower and upper bounds to -5 to 5, which you have currently set as currently -10^-7 and 10^-7. You can further experiment with changing the number of particles and the Hybrid Function.

Refer to these links for more information on particle Swarm Optimization:

https://www.mathworks.com/help/gads/particle-swarm-options.html

https://www.mathworks.com/help/gads/particleswarm.html

You can use Yule walker method for better results.

https://in.mathworks.com/help/signal/ref/yulewalk.html

Hope this helps!

Thanks

Balaji