стратегия для invfreqz в Octave или Matlab

Question

У меня есть коллекция величин и фазовых данных простейшего фильтра нижних частот. Я могу успешно подогнать фильтр на сложный частотный отклик, используя invfreqz в Octave. Однако, что, если я хочу установить тот же фильтр с более высокой частотой дискретизации, скажем, 4 раза, с теми же данными? Тем не менее, нет данных о величине и частоте, доступной до новой 4-кратной частоты Найквиста. Маг и данные фазы не могут быть собраны на более высоких частотах, поэтому его необходимо добавить каким-либо другим методом или приближением.

Что было бы самым простым способом «взломать» величину наполнителя и фазовые данные до новой частоты nyquist, чтобы invfreqz имел бы наилучшие шансы достичь хорошей подгонки для собранных данных с новой частотой дискретизации ?

Answer 1

Вот некоторая функция, которую я написал некоторое время назад, которая может использоваться для смещения частоты дискретизации низких/полосовых фильтров (я никогда не удосужился продлить ее до высокого уровня). В общем, вы должны выбрать точку, в которой ответ вашего фильтра практически не меняется, как частота «среза» (я стесняюсь сказать «отсечка», потому что это не отсечение в нормальном смысле); выше этой частоты ответ устанавливается равным нулю.

function [b, a, fNew, HNew, fOld, HOld, HCut] = ... 
    shiftFilterRate(b, a, fs1, fs2, order, fCut) 
% SHIFTFILTERRATE Shift sampling rate of filter 
% 
% [bNew, aNew] = shiftFilterRate(b, a, fs1, fs2, order, fCut) designs a 
% FIR filter of the given order to operate at sampling frequency fs2, in 
% Hz, whose frequency magnitude characteristics match those of the filter 
% whose coefficients are b, a, that operates at sampling rate fs1, in Hz. 
% The function will only try to match the filter's magnitude in the 
% region [0 fCut], in Hz. 
% 
% [bNew, aNew, fNew, HNew, fOld, HOld, HCut] = shiftFilterRate(...) 
% returns additional parameters: 
%  fNew: The frequencies at which the designed filter's transfer 
%    function was evaluated (with resolution of 1 Hz) 
%  HNew: The transfer function of the designed filter 
%  fOld: The frequencies at which the existing filter's transfer 
%    function was evaluated (with a resolution of 1 Hz) 
%  HOld: The transfer function of the existing filter. 
%  HCut: The desired frequency response of the filter used as input 
%    to the function fir2. 
% 
% FIXME: Make this work for high pass filters. 
% 

if nargin < 5, fCut = inf; end 
if nargin < 4, order = 50; end 

%% Zero padding in frequency domain 
res = 1; % Hz resolution 
N1 = fs1/res; % points at resolution before padding 

% Original freq response 
[H, f] = freqz(b, a, N1); 
nanInd = find(isnan(H)); 
% Stabilise. NOTE: Will break if nanInd contains last elem or there are 
% multiple NaNs in a row 
H(nanInd) = H(nanInd + 1); 
f = f/pi; 

% Normalise cutoff freq 
fCutNorm = fCut/(fs1/2); 

% Cut frequencies above fCut, we don't really need them and it makes the 
% FIR filter nasty 
HCut = H; 
HCut(f > fCutNorm) = 0; 

% Create new freq response 
NNew = ceil(fs2/fs1 * length(HCut)); 
fNew = linspace(0, 1, NNew)'; 
HNew = [HCut; zeros(NNew - N1, 1)]; 

% Design filter 
b = fir2(order, fNew, abs(HNew)); 
a = 1; 

HOld = H; 
fOld = f; 

if nargout > 3 
    HNew = freqz(b, a, length(fNew)); 
end 

Обратите внимание, что в вашем случае, вы можете установить его на -20 дБ на 1, так как это, как представляется, максимальное ослабление ваш фильтр обеспечивает. На этот момент: не похоже, что это очень хороший фильтр ... Есть ли причина, почему у вас есть, чтобы соответствовать этому отклику? Возможно, вам лучше сконструировать, например, Butterworth с тем же самым отсечением (вы наверняка получите большее затухание в полосе остановки и более линейную фазу).