Проект Decimators/Interpolators

В этом примере показано, как спроектировать фильтры для децимации и интерполяции. Обычно фильтры lowpass используются для децимации и для интерполяции. При десятикратном уменьшении фильтры lowpass используются, чтобы уменьшать пропускную способность сигнала до сокращения частоты дискретизации. Это сделано, чтобы минимизировать искажение из-за сокращения частоты дискретизации. При интерполяции фильтры lowpass используются, чтобы удалить спектральные изображения из сигнала с низкой ставкой. Поскольку общие сведения на создании фильтра lowpass видят пример при Разработке Низких КИХ-Фильтров Передачи

Полуполоса Interpolators и Decimators

Чтобы запуститься рассматривают изменение уровня сигнала на коэффициент 2. Полуленточные фильтры являются эффективным способом сделать это. КИХ-фильтры полуполосы реализуются в dsp.FIRHalfbandInterpolator и dsp.FIRHalfbandDecimator. Их БИХ-дубликаты, dsp.IIRHalfbandInterpolator и dsp.IIRHalfbandDecimator может быть еще более эффективный способ интерполировать/десятикратно уменьшить 2.

Fs = 1e6;
hbInterp = dsp.FIRHalfbandInterpolator('TransitionWidth',Fs/10,...
    'SampleRate',Fs);
fvtool(hbInterp) % Notice gain of 2 (6 dB) in the passband
hbDecim  = dsp.FIRHalfbandDecimator('TransitionWidth',Fs/10,...
    'SampleRate',Fs);
fvtool(hbDecim)

Фс частоты дискретизации относится к входному сигналу. В случае интерполяции фильтр сохраняет большую часть спектра от 0 до Фс/2 при ослаблении спектральных изображений. Для децимации фильтр передает приблизительно половину полосы, которая является 0 к Фс/4 и ослабляет другую половину для того, чтобы минимизировать искажение. Объем затухания может быть установлен в любое требуемое значение и для интерполяции и для децимации. Если незаданный, это принимает значение по умолчанию к 80 дБ.

Многоступенчатые проекты

Полуленточные фильтры могут быть расположены каскадом для эффективного многоступенчатого преобразования уровня. Например, интерполяция/десятикратное уменьшение 8 может быть сделана путем расположения каскадом 3 полуполос, interpolators/decimators. Вместо того, чтобы проектировать 3 фильтра вручную, dsp.SampleRateConverter спроектирует все 3, просачивается очень эффективный путь. Обратите внимание на то, что dsp.SampleRateConverter может использоваться для проектов, которые идут вне использования полуленточных фильтров. Для получения дополнительной информации об этом смотрите Эффективное Преобразование Частоты дискретизации Между Случайными факторами.

src = dsp.SampleRateConverter('InputSampleRate',10e3,...
    'OutputSampleRate',8*10e3,'Bandwidth',9e3);
info(src)
ans =

    'Overall Interpolation Factor    : 8
     Overall Decimation Factor       : 1
     Number of Filters               : 3
     Multiplications per Input Sample: 96.000000
     Number of Coefficients          : 66
     Filters:                         
        Filter 1:
        dsp.FIRInterpolator  - Interpolation Factor: 2 
        Filter 2:
        dsp.FIRInterpolator  - Interpolation Factor: 2 
        Filter 3:
        dsp.FIRInterpolator  - Interpolation Factor: 2 
     '

Простые проекты для произвольного преобразования уровня

Иногда желательно спроектировать фильтр для изменения уровня рациональным фактором независимо от фактических включенных частот дискретизации. designMultirateFIR(L,M) проектирует КИХ-фильтр для интерполяции целочисленным коэффициентом L и децимации целочисленным коэффициентом M. designMultirateFIR возвращает коэффициенты фильтра. Эти коэффициенты должны использоваться с dsp.FIRDecimator (L=1), dsp.FIRInterpolator (M=1) и dsp.FIRRateConverter (общий случай).

L         = 5;
M         = 1;
b         = designMultirateFIR(L,M);
firInterp = dsp.FIRInterpolator(L,b);
fvtool(firInterp) % Notice gain of L in the passband

Дополнительные входные параметры к designMultirateFIR допускайте более крутые переходы и лучшее затухание в полосе задерживания. Для более крутого перехода задайте большую половину многофазной длины (значение по умолчанию равняется 12).

hpl        = 20;
b2         = designMultirateFIR(L,M,hpl);
firInterp2 = dsp.FIRInterpolator(L,b2);
h = fvtool(firInterp,firInterp2);
legend(h, 'Polyphase length = 24','Polyphase length = 40')

Более усовершенствованный проект Decimators

При десятикратном уменьшении пропускная способность сигнала уменьшается до соответствующего значения так, чтобы минимальное искажение произошло при сокращении частоты дискретизации. Предположим сигнал, который занимает полный интервал Найквиста (i.e. был критически произведен), имеет частоту дискретизации 800 Гц. Энергия сигнала расширяет до 400 Гц. Если мы хотели бы уменьшать частоту дискретизации на коэффициент 4 - 200 Гц, значительное искажение произойдет, если пропускная способность сигнала не будет также уменьшаться на коэффициент 4. Идеально, совершенный фильтр lowpass с сокращением на уровне 100 Гц использовался бы. На практике несколько вещей произойдут: компоненты сигнала между 0 и 100 Гц будут немного искажены неравномерностью в полосе пропускания неидеального фильтра lowpass; будет некоторое искажение из-за конечного затухания в полосе задерживания фильтра; фильтр будет иметь полосу перехода, которая исказит сигнал в такой полосе. Объемом искажения, введенного каждым из этих эффектов, можно управлять путем разработки соответствующего фильтра. В общем случае, чтобы получить лучший фильтр, более высокий порядок фильтра будет требоваться.

Давайте запустимся путем разработки простого lowpass decimator с фактором децимации 4.

M   = 4;   % Decimation factor
Fp  = 80;  % Passband-edge frequency
Fst = 100;  % Stopband-edge frequency
Ap  = 0.1; % Passband peak-to-peak ripple
Ast = 80;  % Minimum stopband attenuation
Fs  = 800; % Sampling frequency
fdDecim = fdesign.decimator(M,'lowpass',Fp,Fst,Ap,Ast,Fs) %#ok
fdDecim = 

  decimator with properties:

          MultirateType: 'Decimator'
               Response: 'Lowpass'
       DecimationFactor: 4
          Specification: 'Fp,Fst,Ap,Ast'
            Description: {4x1 cell}
    NormalizedFrequency: 0
                     Fs: 800
                  Fs_in: 800
                 Fs_out: 200
                  Fpass: 80
                  Fstop: 100
                  Apass: 0.1000
                  Astop: 80

Технические требования для фильтра решают, что полоса перехода 20 Гц приемлема между 80 и 100 Гц и что минимальное затухание для из компонентов полосы составляет 80 дБ. Также то, что максимальное искажение для компонентов интереса составляет 0,05 дБ (половина неравномерности в полосе пропускания от пика к пику). Фильтр equiripple, который соответствует этим спецификациям, может быть легко получен можно следующим образом:

decimFilter = design(fdDecim,'equiripple', 'SystemObject', true);
measure(decimFilter)
specScope1 = dsp.SpectrumAnalyzer(...
    'PlotAsTwoSidedSpectrum', false, ...
    'SpectralAverages', 50, 'OverlapPercent', 50, ...
    'Title', 'Decimator with equiripple lowpass filter',...
    'YLimits', [-50, 0], 'SampleRate', Fs/M*2);

for k = 1:1000
    inputSig = randn(500,1);            % Input
    decimatedSig = decimFilter(inputSig); % Decimator
    specScope1(upsample(decimatedSig,2)); % Spectrum
    % The upsampling is done to increase X-limits of SpectrumAnalyzer
    % beyond (1/M)*Fs/2 for better visualization
end
release(specScope1);
ans = 

Sample Rate      : 800 Hz     
Passband Edge    : 80 Hz      
3-dB Point       : 85.621 Hz  
6-dB Point       : 87.8492 Hz 
Stopband Edge    : 100 Hz     
Passband Ripple  : 0.092414 dB
Stopband Atten.  : 80.3135 dB 
Transition Width : 20 Hz      
 

Ясно из измерений, что проект соответствует спецификациям.

Используя фильтры Найквиста

Фильтры Найквиста привлекательны для децимации и интерполяции вследствие того, что часть 1/M количества коэффициентов является нулем. Полоса фильтра Найквиста обычно собирается быть равной фактору децимации, это сосредотачивает частоту среза в (1/M) *Fs/2. Обратите внимание на то, что фильтры полуполосы являются фильтрами Найквиста для случая M=2. Отметьте также тот designMultirateFIR проекты Найквист фильтруют также.

В этом примере полоса перехода сосредоточена вокруг 400/M = 100 Гц.

TW = 20; % Transition width of 20 Hz
fdNyqDecim = fdesign.decimator(M,'nyquist',M,TW,Ast,Fs) %#ok
fdNyqDecim = 

  decimator with properties:

          MultirateType: 'Decimator'
               Response: 'Nyquist'
       DecimationFactor: 4
          Specification: 'TW,Ast'
            Description: {2x1 cell}
                   Band: 4
    NormalizedFrequency: 0
                     Fs: 800
                  Fs_in: 800
                 Fs_out: 200
        TransitionWidth: 20
                  Astop: 80

Проект окна Кайзера может быть получен прямым способом.

nyqDecim   = design(fdNyqDecim,'kaiserwin','SystemObject', true);
specScope2 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ...
    'SpectralAverages', 50, 'OverlapPercent', 50, ...
    'Title', 'Decimator with Nyquist filter',...
    'YLimits', [-50, 0],...
    'SampleRate', Fs/M*2);

for k = 1:1000
    inputSig = randn(500,1);            % Input
    decimatedSig = nyqDecim(inputSig);    % Decimator
    specScope2(upsample(decimatedSig,2)); % Spectrum
    % The upsampling is done to increase X-limits of SpectrumAnalyzer
    % beyond (1/M)*Fs/2 for better visualization
end
release(specScope2);

Искажение с фильтрами Найквиста

Предположим, что сигнал, который будет отфильтрован, имеет плоский спектр. После того, как отфильтрованный, это получает спектральную форму фильтра. После сокращения частоты дискретизации этот спектр повторяется с копиями, сосредоточенными вокруг множителей новой более низкой частоты дискретизации. Рисунок спектра подкошенного сигнала может быть найден от:

NFFT = 4096;
[H,f] = freqz(nyqDecim,NFFT,'whole',Fs);
figure;
plot(f-Fs/2,20*log10(abs(fftshift(H))))
grid on
hold on
plot(f-Fs/M,20*log10(abs(fftshift(H))),'r-')
plot(f-Fs/2-Fs/M,20*log10(abs(fftshift(H))),'k-')
legend('Baseband spectrum', ...
    'First positive replica', 'First negative replica')
title('Aliasing with Nyquist filter');
fig = gcf;
fig.Color = 'White';
hold off

Обратите внимание на то, что копии перекрываются несколько, таким образом искажение введено. Однако искажение только происходит в полосе перехода. Таким образом, значительная энергия (выше предписанных 80 дБ) от первой копии только искажает в основную полосу между 90 и 100 Гц. Поскольку фильтр переходил в этой области так или иначе, сигнал был искажен в той полосе, и искажающий там не важно.

С другой стороны, заметьте, что несмотря на то, что мы использовали ту же ширину перехода в качестве с проектом lowpass сверху, мы на самом деле сохраняем более широкую применимую полосу (90 Гц, а не 80) при сравнении этого проекта Найквиста с исходным проектом lowpass. Чтобы проиллюстрировать это, давайте выполним ту же процедуру, чтобы построить спектр подкошенного сигнала, когда проект lowpass сверху будет использоваться

[H,f] = freqz(decimFilter,NFFT,'whole',Fs);
figure;
plot(f-Fs/2,20*log10(abs(fftshift(H))))
grid on
hold on
plot(f-Fs/M,20*log10(abs(fftshift(H))),'r-')
plot(f-Fs/2-Fs/M,20*log10(abs(fftshift(H))),'k-')
legend('Baseband spectrum', ...
    'First positive replica', 'First negative replica')
title('Aliasing with lowpass filter');
fig = gcf;
fig.Color = 'White';
hold off

В этом случае нет никакого значительного перекрытия (выше 80 дБ) между копиями, однако потому что область перехода запустилась на уровне 80 Гц, получившийся подкошенный сигнал имеет меньшую применимую пропускную способность.

Интерполяция

При интерполяции сигнала основополосный ответ сигнала нужно оставить максимально неизменным. Интерполяция получена путем удаления спектральных копий, когда частота дискретизации увеличена.

Предположим, что нам произвели сигнал на уровне 48 Гц. Если это критически производится, существует значительная энергия в сигнале до 24 Гц. Если бы мы хотим интерполировать на коэффициент 4, мы идеально спроектировали бы фильтр lowpass, достигающий 192 Гц с сокращением на уровне 24 Гц. Как с децимацией, на практике приемлемая ширина перехода должна быть включена в проект фильтра lowpass, используемого для интерполяции наряду с неравномерностью в полосе пропускания и конечным затуханием в полосе задерживания. Например, рассмотрите следующие спецификации:

L   = 4;   % Interpolation factor
Fp  = 22;  % Passband-edge frequency
Fst = 24;  % Stopband-edge frequency
Ap  = 0.1; % Passband peak-to-peak ripple
Ast = 80;  % Minimum stopband attenuation
Fs  = 48;  % Sampling frequency
fdInterp = fdesign.interpolator(L,'lowpass',Fp,Fst,Ap,Ast,Fs*L) %#ok
fdInterp = 

  interpolator with properties:

          MultirateType: 'Interpolator'
               Response: 'Lowpass'
    InterpolationFactor: 4
          Specification: 'Fp,Fst,Ap,Ast'
            Description: {4x1 cell}
    NormalizedFrequency: 0
                     Fs: 192
                  Fs_in: 48
                 Fs_out: 192
                  Fpass: 22
                  Fstop: 24
                  Apass: 0.1000
                  Astop: 80

Проект equiripple, который соответствует спецификациям, может быть найден таким же образом как с decimators

EQRInterp  = design(fdInterp,'equiripple','SystemObject', true);
specScope4 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ...
    'SpectralAverages', 50, 'OverlapPercent', 50, ...
    'Title', 'Interpolator with equiripple lowpass filter', ...
    'SampleRate', Fs*L);

for k = 1:1000
    inputSig = randn(500,1);      % Input
    intrpSig = EQRInterp(inputSig); % Interpolator
    specScope4(intrpSig);           % Spectrum
end
release(specScope4);

Заметьте, что фильтр имеет усиление 6 dBm. В общих интерполяторах будет иметь усиление равным коэффициенту интерполяции. Это необходимо для сигнала, интерполируемого, чтобы обеспечить ту же область значений после интерполяции. Например,

release(EQRInterp);
sine = dsp.SineWave('Frequency', 18, 'SampleRate', Fs, ...
                    'SamplesPerFrame', 100);
intrpSig  = EQRInterp(sine());
arrayPlot = dsp.ArrayPlot('YLimits', [-2, 2], ...
                          'Title', 'Sine wave interpolated');
arrayPlot(intrpSig(200:300)) % Plot the output

Обратите внимание на то, что несмотря на то, что фильтр имеет усиление 4, интерполированный сигнал имеет ту же амплитуду как исходный сигнал.

Использование фильтров Найквиста для интерполяции

Подобно случаю децимации фильтры Найквиста привлекательны в целях интерполяции. Кроме того, учитывая, что существует коэффициент, равный нулю каждый выборки L, использование фильтров Найквиста гарантирует, что выборки от входного сигнала сохраняются неизменные при выходе. Дело обстоит не так для других фильтров lowpass, когда используется для интерполяции (с другой стороны, искажение может быть минимальным в других фильтрах, таким образом, это - не обязательно огромное соглашение).

TW = 2;
fdNyqIntrp = fdesign.interpolator(L,'nyquist',L,TW,Ast,Fs*L) %#ok
nyqInterp  = design(fdNyqIntrp,'kaiserwin', 'SystemObject', true);
specScope5 = dsp.SpectrumAnalyzer('PlotAsTwoSidedSpectrum', false, ...
    'SpectralAverages', 30, 'OverlapPercent', 50, ...
    'Title', 'Interpolator with Nyquist filter',...
    'SampleRate', Fs*L);

for k = 1:1000
    inputSig = randn(500,1);      % Input
    intrpSig = nyqInterp(inputSig); % Decimator
    specScope5(intrpSig);           % Spectrum
end
release(specScope5);
fdNyqIntrp = 

  interpolator with properties:

          MultirateType: 'Interpolator'
               Response: 'Nyquist'
    InterpolationFactor: 4
          Specification: 'TW,Ast'
            Description: {2x1 cell}
                   Band: 4
    NormalizedFrequency: 0
                     Fs: 192
                  Fs_in: 48
                 Fs_out: 192
        TransitionWidth: 2
                  Astop: 80

Аналогичным способом к децимации, когда используется для интерполяции, фильтры Найквиста позволяют определенную степень обработки изображений. Таким образом, некоторые частоты выше частоты среза не ослабляются значением Ast. Однако это происходит только в полосе перехода фильтра. С другой стороны, еще раз более широкий фрагмент основной полосы исходного сигнала обеспечен неповрежденный, когда по сравнению с lowpass фильтруют с ребром полосы задерживания в идеальной частоте среза, когда оба фильтра имеют ту же ширину перехода.