Series Filter Loop Error

Question

S 2024-1-31

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2076931-series-filter-loop-error

评论： Mathieu NOE 2024-2-27

I am trying to make a loop that goes through different filters in cascade, so the output of one is the input to the next. I want to use the built in matlab gammatone function to do this. I have been stuck here for a bit. Thank you for your time!

fs = 20e3; 
numFilts = 32;
BW = 100; %Filter Bandwidth
filter_number = 5;
CenterFreqs = logspace(log10(50), log10(8000), numFilts);
CF1 = CenterFreqs - BW/2; %Lower cutoff frequency
CF2 = CenterFreqs + BW/2; %Upper cutoff frequency
t = linspace(0,2*pi,200);
input = sin(t) + 0.25*rand(size(t));
for ii = 1:filter_number
    output = gammatone(input, CenterFreqs, fs)
end
Unrecognized function or variable 'gammatone'.
[h{ii},f] = freqz(gammatone{ii}.Coefficients,1,4*8192,fs);

Error: (I do not have a line 32?, also I did not make gammatone so I don't think there is anything wrong with it)

Not enough input arguments.

Error in gammatone (line 32)

A = 10^((q_factor) / 20);

Error in t131 (line 13)

output = gammatone(input, CenterFreqs, fs)

15 个评论
显示 13更早的评论隐藏 13更早的评论

Mathieu NOE 2024-2-1

在 MATLAB Online 中打开

well , seems to me there is more than just one problem

according to your first post (error message : Unrecognized function or variable 'gammatone'.) it was only a problem that the gammatone function seemed to be unknow from the matlab path

when I use your gammatone code above and create that function , it works in terms of calling this function and getting some output, but what I see is not what I expected (or there is missing code )

also you must correct the line inside the for loop as CenterFreqs must be indexed

for ii = 1:filter_number

output = gammatone(input, CenterFreqs(ii), fs)

end

output

For me a gammatone should be some king of bandpass filter , which is obviously not what the function does here

if we look at the ouput , we have this , and most funny the output is even not function of the input (except i has same size)

where does this function come from and what was the supposed purpose of it ?

I suspect it has to do with generating the impulse response of the gammatone filter , but for me the function (as provided) is incomplete

Gammatone filter - Wikipedia

are you trying to combine input and output to generate the impulse response of the gammatone filter ? then the "true" filtering must happen latter in the code

fs = 20e3;

numFilts = 32; % why this ? not used

BW = 100; %Filter Bandwidth

filter_number = 5;

CenterFreqs = logspace(log10(50), log10(8000), numFilts);

CF1 = CenterFreqs - BW/2; %Lower cutoff frequency

CF2 = CenterFreqs + BW/2; %Upper cutoff frequency

t = linspace(0,2*pi,200);

input = sin(t) + 0.25*rand(size(t));

figure

hold on

for ii = 1:filter_number

output = gammatone(input, CenterFreqs(ii), fs);

plot(output)

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function output = gammatone(input, centerFrequency, fs)

% Design a gammatone filter

earQ = 9.26449;

minBW = 24.7;

erb = (centerFrequency / earQ) + minBW;

B = 1.019 .* 2 .* pi .* erb;

g = 10^(-3 / 20); % 3 dB down from peak

f = centerFrequency;

tau = 1. / (2. .* pi .* B);

% Apply the gammatone filter

t = 0:1/fs:(length(input)-1)/fs;

temp = (t - tau) > 0;

output = (t.^(4 - 1)) .* exp(-2 .* pi .* B .* (t - tau)) .* g .* cos(2 .* pi .* f .* (t - tau)) .* temp;

end

Mathieu NOE 2024-2-2

在 MATLAB Online 中打开

you must be calling another gammatone function

if you look again to the code I posted above , the gammatone function (from your own post above) has only 13 lines

so your error message at line 32 inside YOUR gammatone function clearly indicates we are not using the same function !

Error in gammatone (line 32)

A = 10^((q_factor) / 20);

NB also that there is no such code ( A = 10^((q_factor) / 20); ) below

can you please share your gammatone function ? It must have more than 32 lines as we know from now.

again , the code I posted above (according to your previous inputs) :

I highly suspect this is not the right gammatone function

fs = 20e3; 
numFilts = 32; % 
BW = 100; %Filter Bandwidth
filter_number = 5;
CenterFreqs = logspace(log10(50), log10(8000), numFilts);
CF1 = CenterFreqs - BW/2; %Lower cutoff frequency
CF2 = CenterFreqs + BW/2; %Upper cutoff frequency
t = linspace(0,2*pi,200);
input = sin(t) + 0.25*rand(size(t));
figure
hold on
for ii = 1:filter_number
    output = gammatone(input, CenterFreqs(ii), fs);
    plot(output)
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 function output = gammatone(input, centerFrequency, fs)
      % Design a gammatone filter
      earQ = 9.26449;
      minBW = 24.7;
      erb = (centerFrequency / earQ) + minBW;
      B = 1.019 .* 2 .* pi .* erb;
      g = 10^(-3 / 20); % 3 dB down from peak
      f = centerFrequency;
      tau = 1. / (2. .* pi .* B);
  
      % Apply the gammatone filter
      t = 0:1/fs:(length(input)-1)/fs;
      temp = (t - tau) > 0;
      output = (t.^(4 - 1)) .* exp(-2 .* pi .* B .* (t - tau)) .* g .* cos(2 .* pi .* f .* (t - tau)) .* temp;
  end

Mathieu NOE 2024-2-6

在 MATLAB Online 中打开

To your first comment above

I think youy are mixing some code from different horizons , I believe I have already seen this particular line (or something very similar) in one of your other post (same topic indeed) :

[h{ii},f] = freqz(gammatone{ii}.Coefficients,1,4*8192,fs);

Notive here gammatone is not a structure, it's a function that simply returns output (array) as we have declared it :

function output = gammatone(input, centerFrequency, fs)

so there is not gammatone{ii}.Coefficients field !! and you cannot do a freqz action on a non existent field

maybe you can do this (output is assumed to be the FIR coefficients as in the other post) , but I doubt that the result is matching your expectations , as the Bode plot does not look like a bank of gammatone filters

fs = 20e3;

numFilts = 32; %

BW = 100; %Filter Bandwidth

filter_number = 5;

CenterFreqs = logspace(log10(50), log10(8000), numFilts);

CF1 = CenterFreqs - BW/2; %Lower cutoff frequency

CF2 = CenterFreqs + BW/2; %Upper cutoff frequency

t = linspace(0,2*pi,200);

input = sin(t) + 0.25*rand(size(t));

figure

hold on

for ii = 1:filter_number

output = gammatone(input, CenterFreqs(ii), fs);

% plot(output)

[h{ii},f] = freqz(output,1,1024,fs);

plot(f,20*log10(abs(h{ii})))

end

set(gca,'XScale','log');

xlabel('Freq(Hz)');

ylabel('Gain (dB)');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function output = gammatone(input, centerFrequency, fs)

% Design a gammatone filter

earQ = 9.26449;

minBW = 24.7;

erb = (centerFrequency / earQ) + minBW;

B = 1.019 .* 2 .* pi .* erb;

g = 10^(-3 / 20); % 3 dB down from peak

f = centerFrequency;

tau = 1. / (2. .* pi .* B);

% Apply the gammatone filter

t = 0:1/fs:(length(input)-1)/fs;

temp = (t - tau) > 0;

output = (t.^(4 - 1)) .* exp(-2 .* pi .* B .* (t - tau)) .* g .* cos(2 .* pi .* f .* (t - tau)) .* temp;

end

Mathieu NOE 2024-2-13

在 MATLAB Online 中打开

I will transform my comment above into a answer so you can accept it (and tx for that ! )

for the g factor , we don't need it inside the gammatone function itself as in the main code we are doing it in these lines :

(up to you to choose a 0 or -3 dB peak amplitude of the filters)

    % scale the IR amplitude so that the max modulus is 0 dB
    a = 1/max(abs(tmp));
%     % or if you prefer - 3dB
%     g = 10^(-3 / 20); % 3 dB down from peak
%     a = g/max(abs(tmp));

Mathieu NOE 2024-2-13

在 MATLAB Online 中打开

FYI, the lines of code above could also be written this way if you prefer

    % scale the IR amplitude to match a required Bode plot magnitude
    g = 1; % the max modulus is 0 dB
%     g = 10^(-3 / 20);    % or if you prefer - 3dB
     a = g/max(abs(tmp));

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Mathieu NOE 2024-2-7

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2076931-series-filter-loop-error#answer_1407958

移动：Mathieu NOE 2024-2-13

在 MATLAB Online 中打开

hello again

so I wanted to learn a bit about gammatone filters , so I read this , which I suppose i also from where you started your code

PattersonÕs Ear Model (purdue.edu)

there are a few modifications I made in the code , see the original lines commented are replaced by new code

so the function output is the IR of the filter , we can use it latter on (not yet coded) to filter any time signal

IMHO, now the IR and Bode plots are correct

here we plot only the first 5 filters , even though the CF are computed for 32 filters

fs = 20e3;

numFilts = 32; %

% BW = 100; %Filter Bandwidth

filter_number = 5;

order = 4;

CenterFreqs = logspace(log10(50), log10(8000), numFilts);

% CF1 = CenterFreqs - BW/2; %Lower cutoff frequency

% CF2 = CenterFreqs + BW/2; %Upper cutoff frequency

% t = linspace(0,2*pi,200);

% input = sin(t) + 0.25*rand(size(t));

figure(1)

hold on

for ii = 1:filter_number

IR = gammatone(order, CenterFreqs(ii), fs);

[tmp,f] = freqz(IR,1,1024*2,fs);

% scale the IR amplitude so that the max modulus is 0 dB

a = 1/max(abs(tmp));

% % or if you prefer - 3dB

% g = 10^(-3 / 20); % 3 dB down from peak

% a = g/max(abs(tmp));

IR = IR*a;

[h{ii},f] = freqz(IR,1,1024*2,fs); % now store h{ii} after amplitude correction

plot(IR)

end

title('Impulse Response');

hold off

figure(2)

hold on

for ii = 1:filter_number

plot(f,20*log10(abs(h{ii})))

end

title('Bode plot');

set(gca,'XScale','log');

xlabel('Freq(Hz)');

ylabel('Gain (dB)');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function IR = gammatone(order, centerFrequency, fs)

% Design a gammatone filter

earQ = 9.26449;

minBW = 24.7;

% erb = (centerFrequency / earQ) + minBW;

erb = ((centerFrequency/earQ)^order + (minBW)^order)^(1/order);% we use the generalized

% function (depending on order)

% B = 1.019 .* 2 .* pi .* erb; % no, the 3pi factor is implemented twice in your code

B = 1.019 * erb;

% g = 10^(-3 / 20); % 3 dB down from peak % what is it for ? see main code above

f = centerFrequency;

tau = 1. / (2. .* pi .* B);

% gammatone filter IR

t = (0:1/fs:10/f);

temp = (t - tau) > 0;

% IR = (t.^(order - 1)) .* exp(-2 .* pi .* B .* (t - tau)) .* g .* cos(2*pi*f*(t - tau)) .* temp;

IR = ((t - tau).^(order - 1)) .* exp(-2*pi*B*(t - tau)).* cos(2*pi*f*(t - tau)) .* temp;

end

4 个评论
显示 2更早的评论隐藏 2更早的评论

S 2024-2-26

@Mathieu NOE I just posted it as a new question in case you wanted me to repost!

Mathieu NOE 2024-2-27

hello @S

ok , I'll have a look in your other post !

请先登录，再进行评论。

Series Filter Loop Error

15 个评论
显示 13更早的评论隐藏 13更早的评论

采纳的回答

4 个评论
显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

Series Filter Loop Error

15 个评论 显示 13更早的评论隐藏 13更早的评论

采纳的回答

4 个评论 显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

15 个评论
显示 13更早的评论隐藏 13更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论