1

我的目标是获得特定中心频率的 1/3 倍频程频带频率信号。我想使用 Matlab 的 octaveFilter 功能,但我希望在 1000 Hz 处获得单个 1/3 倍频程频带峰值,但计算出 1000 Hz 的最左边和最右边的非常正的声压级。

见输出

我做错了什么?

    Fs = 48000;                         % Sampling rate
    t = 0:1/Fs:1-1/Fs;                  % Time vector of 1 second
    f = 1000;                           % Frequency of signal
    dpres = sin(2*pi*f*t);              % Signal in [Pa]

    frCen = [100 300 600 800 1000 1200 1300 1600];

    [Spl,frCen] = CompOctSplFreq(dpres,frCen)

    figure()
    semilogx(frCen,Spl,'ro-')
    title('1/3-Octave Filtered SPL over Frequency')
    xlabel('Center Frequency of Octave Band Filter [Hz]')
    ylabel('SPL [dB]')        

            function [Spl,frCen] = CompOctSplFreq(dpres,frCen);

                % Reads singal in pressure over time
                freqNum = length(frCen);
                Fs = length(dpres);

                % Setting up the 1/3-octave filters for all center frequencies
                for f = 1:freqNum
                    octaveFilterBank{f} = octaveFilter(frCen(f),'1/3 octave','SampleRate',Fs,'FilterOrder',12);
                    % Filtering the signal with the corresponding filters
                    dpresFiltered(:,f) = octaveFilterBank{f}(dpres');
                    % Getting the average for each filter frequency band
                    drms(f) = sqrt(sum(dpresFiltered(:,f).^2)/length(dpresFiltered(:,f)));
                end

                % Converting the root mean square pressure to SPL            
                pRef = 20e-06;                      % Reference pressure            
                Spl = 20*log10(drms/pRef); 


            end
4

1 回答 1

1

主要问题是您正在计算信号的功率,在某些情况下包括显着的瞬态响应。为了说明这一点,下面将推导预期的稳态响应,并比较计算得到的响应,包括和不包括瞬态部分。

基于octaveFilter 的算法描述,可以构建滤波器组中每个滤波器的理想稳态响应。下图中给出了这些响应:

过滤响应

请注意,响应确实在其指定范围之外显着下降,但它们并没有完全消失。因此,频率超出指定频带的输入信号仍会在滤波器组中的所有滤波器上产生一些输出(尽管很小)。对于 1000Hz 的正弦输入信号,您可以通过在上图中用 1000Hz 垂直线查找相应滤波器响应的截距来获得每个滤波器输​​出的稳态响应。这应该为您提供以下图表,了解滤波器组对 1000Hz 时 95dB 输入的响应:

对 1000Hz 输入的预期稳态响应

现在您可能会注意到,对于以 1000Hz 为中心的滤波器,您的结果非常接近这些稳态响应,但随着我们远离 1000Hz,差异会增加。如果我们查看时域中的一个滤波信号,我们可能会注意到在信号的初始瞬态期间确实可以找到很大一部分能量:

时域瞬态响应

简单地摆脱瞬态应该会让您更接近预期的稳态响应。更棘手的部分是弄清楚需要丢弃多少样本。作为一般经验法则,我通常使用大约 5 个时间常数,其中一个时间常数大约是通带滤波器半带宽的倒数。这可以通过以下方式计算:

G = 10^(3/10);
b = 3; % for 1/3-octave
fhi  = frCen(f)*G^(0.5/b);
flow = frCen(f)*G^(-0.5/b);
M = 10/(fhi-flow); % 5 time-constants (one time-constant is ~ 1/((fhi-flow)/2) )
drms(f) = sqrt(sum(dpresFiltered(M:end,f).^2)/length(dpresFiltered(M:end,f)));

将不包括与瞬态对应的初始样本的结果响应与您之前的结果结果进行比较,将表明这确实对图表产生了重大影响,并且现在的结果与我们之前的理想稳态响应图表更加一致:

比较

于 2018-08-11T15:10:07.407 回答