1
$\begingroup$

I am trying to build a filter with different frequency bands in Matlab. I would like the filter to be as invisible as possible, but after the filtering, a great amount of frequencies are being attenuated and I believe it is where the filters overlap, but also one can see that the low frequencies are also attenuated though there is not any overlap there.

The following images show the results:

Original signal

enter image description here

Filtered signal

enter image description here

Could you help me understand what I am doing wrong?

This is the code I am using:

[x,Fs] = audioread('drum-loop.wav');

% define constants  - - - - - - - - - - - - - - - - - - - - - - - - - - - - 

ntaps = 255;                                                               % number of filter-taps

% sampling frequencies for each branch
fs_A = 48000;                   %
fs_B = 24000;                   %
fs_C = 6000;                    %
fs_D = 1500;                    %

% decimation for each branch
dec_A = 1;                      %
dec_B = 2;                      %
dec_C = 4;                      %
dec_D = 4;                      %

% decimation filters  - - - - - - - - - - - - - - - - - - - - - - - - - - -
x_B = decimate(x, dec_B);
x_C = decimate(x_B, dec_C);
x_D = decimate(x_C, dec_D);

% interpolation filters   - - - - - - - - - - - - - - - - - - - - - - - - -
% ref: https://es.mathworks.com/help/dsp/ref/mfilt.html
intfilter_B = mfilt.firinterp(dec_B);                                      % decimate signal for branch B
intfilter_C = mfilt.firinterp(dec_C);                                      % decimate signal for branch C
intfilter_D = mfilt.firinterp(dec_D);                                      % decimate signal for branch D


% define filters  - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
ovlp1 = 2000;
ovlp2 = 500;
ovlp3 = 250;

% branch A: band pass filter.
bp_A = designfilt('bandpassfir','FilterOrder',ntaps-1, ...
         'CutoffFrequency1',fs_B/2-ovlp1,'CutoffFrequency2',(fs_A/2-1), ...
         'SampleRate',fs_A);
% branch B: band pass filter.
bp_B = designfilt('bandpassfir','FilterOrder',ntaps-1, ...
         'CutoffFrequency1',fs_C/2-ovlp2,'CutoffFrequency2',(fs_B/2-1), ...
         'SampleRate',fs_B);
% branch C: band pass filter.
bp_C = designfilt('bandpassfir','FilterOrder',ntaps-1, ...
         'CutoffFrequency1',fs_D/2-ovlp3,'CutoffFrequency2',(fs_C/2-1), ...
         'SampleRate',fs_C);
% branch D: band pass filter.
bp_D = designfilt('lowpassfir','FilterOrder',ntaps-1, ...
         'CutoffFrequency',(fs_D/2-1), 'SampleRate',fs_D);

% apply filters
x_A_f = filtfilt(bp_A,x);
x_B_f = filtfilt(bp_B,x_B);
x_C_f = filtfilt(bp_C,x_C);
x_D_f = filtfilt(bp_D,x_D);

% summation filter
y = x_A_f + interp(x_B_f,dec_B) + interp(interp(x_C_f,dec_C),dec_B) + ...
    interp(interp(interp(x_D_f,dec_D),dec_C),dec_B);
y = y/max(y);                                                              % normalize

% write audio
audiowrite('output.wav',y,48000);

% plot figures - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
% compute fft
figure;
subplot(3,1,1);plot(x);title('original signal')
NFFT = 2^nextpow2(length(x));                                              % Next power of 2 from length of y (in samples)
Y = fft(x,NFFT)/Fs;     
f = Fs/2*linspace(0,1,NFFT/2+1);
subplot(3,1,2);spectrogram(x,blackman(128),60,128,1e3,'yaxis');title('spectrogram')
% subplot(3,1,3);plot(f,2*angle(Y(1:NFFT/2+1)));title('phase')
subplot(3,1,3);plot(f,abs(Y(1:NFFT/2+1)));title('fft') 

% compute fft
figure;
subplot(3,1,1);plot(y);title('filtered signal')
NFFT = 2^nextpow2(length(y));                                           % Next power of 2 from length of y (in samples)
Y = fft(y,NFFT)/Fs;     
f = Fs/2*linspace(0,1,NFFT/2+1);
subplot(3,1,2);spectrogram(y,blackman(128),60,128,1e3,'yaxis');title('spectrogram')
% subplot(3,1,3);plot(f,2*angle(Y(1:NFFT/2+1)));title('phase')
subplot(3,1,3);plot(f,abs(Y(1:NFFT/2+1)));title('fft')
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.