DSPRelated.com
Forums

Modify a FIR filter in the frequency domain a then accomplish IFFT

Started by Swift80 2 months ago8 replieslatest reply 2 months ago176 views

Hallo,

I would like to design a custom made FIR filter. I would start from an ordinary LP FIR filter designed by means of fdatool, for instance, and then modify the magnitude freq response of the filter. I would like to alter the pass-band of the FIR, passing from a flat response to an increasing slope in the pass-band.

Would it be possible to modify it in the frequency domain a the switch back to time domain by doing an IFFT and thus obtaining the filter's coefficients ?

Here is the script I used (I also modified the phase). 

It seems that I get aliasing once back in the time domain.  Why is that and how what I want to do can be accomplished ?

Thanks a lot

--------

clear all; close all; clc;

load('X_Band_sgn.mat') % Load the oriinal FIR filter (flat pass-band)

fc = 715e6;

dt = 1/fc;

duration = 10e-3;

dur = length(X_Band_sgn);

time = (0:dt:(dur-1)*dt);

df = fc/dur;

freq = (0:df:fc-df)-fc/2;

a0 = 2.25*ones(1,100);

a1 = -216*ones(1,100);

int1 = interp1([1, 201], [2.25,-216], 1:201);

ang_mod = [a0, int1, a1];

Yt = fftshift(fft(X_Band_sgn));

figure; grid minor; hold on;

plot(freq/1e6,10*log10(abs(Yt).^2))

plot(freq/1e6,unwrap(angle(Yt)))

plot(freq/1e6,ang_mod)

xlabel('Frequency (MHz)')

ylabel('Power [dB]')

s_dB = 10*log10(abs(Yt).^2);

id = find(s_dB > -1.5);

s_dB = 10*log10(abs(Yt).^2)+1;

ang = angle(Yt);


%design the rect to alterate the spectrum mag

k = 40;

dy = 0.12*k;

dx = id(end)-id(1);

y = dy/dx*(1:1:id(end)-id(1))+1;

vq = [ones(1,id(1)) y ones(1,dur-id(end))];

sp_dB = s_dB.*vq;

figure; grid minor; hold on;

plot(freq/1e6,sp_dB)

plot(freq/1e6,unwrap(angle(Yt)))

xlabel('Frequency (MHz)')

ylabel('Power [dB]')

xlim([-357.5 357.5])

ylim([-200 60])

sp_dB = sp_dB-1;

figure; grid minor; hold on;

plot(freq/1e6,sp_dB)

plot(freq/1e6,unwrap(angle(Yt)))

xlabel('Frequency (MHz)')

ylabel('Power [dB]')

xlim([-357.5 357.5])

ylim([-200 60])


% New filter

sp_X = (10.^(sp_dB/20)).*exp(1i*(wrapToPi(ang_mod)));

sp_t     = ifft(ifftshift(sp_X)); %IFFT


Yt = fftshift((fft(sp_t)));

figure; grid minor; hold on;

plot(freq/1e6,10*log10(abs(Yt).^2))

plot(freq/1e6,unwrap(angle(Yt)))

plot(freq/1e6,(angle(Yt)))

xlabel('Frequency (MHz)')

ylabel('Power [dB]')

xlim([-357.5 357.5])

ylim([-190 60])

X_Band_sgn_mod = sp_t;

save('.\X_Band_sgn_mod.mat', 'X_Band_sgn_mod');


% generate x signal 

pow = 15;

duration = 100e-3;

dur = round(duration*fc);

time = (0:dt:(dur-1)*dt);

df = fc/dur;

freq = (0:df:fc-df)-fc/2;

noise1_i = wgn(dur,1,pow,50,'dBm');

noise1_q = wgn(dur,1,pow,50,'dBm');

noise_iq_app = (noise1_i.'+1i*noise1_q.');


% filtering x signal by using the new designed filter

yX_Band = filter(X_Band_sgn_mod,1,noise_iq_app); 

%Shift the freq

y_X_Band_tot = yX_Band.*exp(-1j*2*pi*250e6*time);


Yt = fftshift(fft(y_X_Band_tot));

figure; grid minor; hold on;

plot(freq/1e6,10*log10(abs(Yt).^2))

xlabel('Frequency (MHz)')

ylabel('Power [dB]')


[ - ]
Reply by kazFebruary 13, 2024

I think you are on the wrong path.

FIR2 function can do an arbitrary response.

[ - ]
Reply by Robert WolfeFebruary 13, 2024

this recent thread discusses some of the potential perils of assuming that a FFT modification, will inverse back to the time domain as desired/assumed:

filtering using FFT/iFFT (dsprelated.com)

[ - ]
Reply by kazFebruary 13, 2024

Keep in mind that filter design(coefficients) is distinct from applying filter.

Much of above linked thread is about applying filter. For coefficient design using inverse FFT is very common practice (FIR2 is based on that)

[ - ]
Reply by Robert WolfeFebruary 13, 2024

Right.  Thanks for that clarification.

[ - ]
Reply by Swift80February 13, 2024

Thanks a lot for the answer. I'll take a look at the fir2 function to design the filter with arbitrary response.

[ - ]
Reply by Swift80February 13, 2024

Thanks for the reply Robert

[ - ]
Reply by neiroberFebruary 13, 2024

You can also design for an arbitrary response using the Matlab function firpm (Parks-McClellan synthesis).


--Neil


[ - ]
Reply by Swift80February 13, 2024

Thanks a lot Neil, I'll take a look and try