How to apply FIR1 Bandpass filter ?

5 visualizaciones (últimos 30 días)
majid
majid el 28 de Jun. de 2021
Comentada: majid el 28 de Jun. de 2021
Hi,
I have noisy chirp signal which I want to apply FIR1 Bandpass filter with hanning window to final mixed signal in order to get intermediate frequency signal which I show below and will be derive after applying filtering.
and also may I ask you to explain how we should choose the bandpass filter frequency? frequencies that I used send me errors! and I read help but I cant understande where is my problem in code.
Here is my code:
%%
set(0,'defaultlinelinewidth',2)
clear all
clc
%chirp signal generation
fs=10000;
Ts=1/fs;
t = 0:1/fs:20.5e-3;
f0=20; f1=320e1; t1=20.5e-3;% lowered frequencies to show sampling problem
a=(f1-f0)/(t(end)-t(1));
f=f0+a.*t;
y1 = chirp(t,f0,t1,f1,'linear');
Y = [y1 y1(end:-1:1)]; % y1(end:-1:1) is the time mirrored version of y1, without need to call chirp function again
t1=0:1/fs:41.001e-3;
% figure; plot(t1,Y); title('Bidirectional chirp signal ');
% xlabel('Time(s)'); ylabel('Amplitude');
N=16;
% YY1=repmat(Y,1,N);
% t_YY1=0:1/fs:656.021e-3 ;
% plot(t_YY1,YY1)
% title('Bidirectional chirp signal ');
% xlabel('Time(s)'); ylabel('Amplitude');
clear t_YY1 Y y1 YY1 t1 t
%--------------------------------------------------------------------
%% Delayed signal
R=60;
delta_t=2*R/(3e8);
t_Recieved =delta_t:1/fs:20.5e-3+delta_t;
t_target_Recieved=t_Recieved(end);
y_Recieved =chirp(t_Recieved,f0,t_target_Recieved,f1,'linear');
Y_Recieved = [y_Recieved y_Recieved(end:-1:1)];
t_Y_Recieved=delta_t:1/fs:41.001e-3+delta_t;
% figure; plot(t_Y_Recieved.*1000,Y_Recieved);
YY_Recieved=repmat(Y_Recieved,1,N);
t_YY_Recieved=delta_t:1/fs:659.1e-3+delta_t ;
% plot(t_YY_Recieved,YY_Recieved)
% title('Bidirectional chirp signal ');
% xlabel('Time(s)'); ylabel('Amplitude');
clear delta_t t_Recieved t_target_Recieved t_transmit t_Y_Recieved
clear y_Recieved Y_Recieved
%% adding noise
%speckle noise
M=ones(size(YY_Recieved)); %coherence parameter
x=rand(size(YY_Recieved));
i_speckle=icdf('nbin',x,M,M./(YY_Recieved+M));
clear M x
% figure;
% plot(t_YY_Recieved.*1000,i_speckle);title(['speckle noise']);
% xlabel('Time (ms)');
% axis tight
% figure
% sig=sqrt(YY_Recieved.^2+i_speckle.^2);
% plot(t_YY_Recieved.*1000,sig);title(['signal + speckle noise']);
% axis tight
%--------------------------------------------------------------------------
%thermal noise
qe=1.602*10e-19; %electron charge
kb=1.3806504*10^(-23);
R_fb=1e6; %feedback resistance (MOhm)
To=300;%absolute temperature
B_det=1.5/(2*Ts);% IF detection band noise bandwidth
quantum_eff=0.2; %detector quantum efficiency
NEPh_dk_pc=1E+6;%Noise equivalent photon rate at the photocathode per sub-pixel from dark current
N_pix_x=14;%number of sub-pixel binned in the x direction to form a super pixel
N_pix_y=14;%number of sub-pixel binned in the y direction to form a super pixel
i_thermal=sqrt(2*qe*kb*To*B_det*N_pix_x*N_pix_y/R_fb)*randn(size(YY_Recieved));%TIA thermal noise
% figure;
% plot(t_YY_Recieved.*1000,i_thermal);
% xlabel('Time (ms)');title(['thermal noise']);
% axis tight
% figure
% sig=sqrt(i_thermal.^2+YY_Recieved.^2);
% plot(t_YY_Recieved.*1000,sig);title(['signal + thermal noise']);
% axis tight
clear kb R_fb To B_det sig
%--------------------------------------------------------------------------
%shot noise modeled by white noise
SNR=20;
i_shot=awgn(YY_Recieved,SNR);% white noise adding
% figure;
% plot(t_YY_Recieved.*1000,abs(i_shot));title(['shot noise modeled by white gaussian noise']);
% xlabel('Time (ms)');
% axis tight
% figure
% sig=sqrt(Y_Recieved.^2+i_shot.^2);
% plot(t_YY_Recieved.*1000,sig);title(['signal + gaussian noise']);
% axis tight
%--------------------------------------------------------------------------
noisy_sig1=i_thermal.^2+i_shot.^2+i_speckle.^2;
clear SNR sig i_thermal i_shot i_speckle x_8
%--------------------------------------------------------------------------
% backreflection , backscattering , atmospheric effects
P_l_pk=7.4; %laser peak power (watt)
FF=0.49; %Reciever FPA fill factor(%)
t_t=0.7; %transmitter transmitance
t_r=0.338; %reciever transmitance
D_r=0.254; %Reciever diameter (m)
rho_tgt=0.1; %target reflectance
kapa_tgt_gauss=1; %gaussian beam target return geometrical form factor
alpha=0.0002; %atmospheric extinction coefficient
rho_bg=0.3; %background reflectance
kapa_ls_gauss=1;%gaussian spot geometrical form factor
beta_vol_bs=9.95e-6;
Hs=267; %solar spectral irradiance (watt/(m^2 * um))
N_sky=3; %sky spectral irradiance (watt/(m^2 * sr*um))
delta_lambda=6.2e-9; %filter bandwidth
delta=12e-5; %reciever sub pixel size
Reciever_focal=1.016; %(m)
A_ifov=(R*delta/Reciever_focal)^2;
D_obst=0.5; %reciver central obstruction diameter
A_r=pi*D_r^2/4;%reciever aperture area
A_obst=pi*D_obst^2/4;%reciever central obstruction area
omega_rcv=(A_r-A_obst)/R^2;%solid angle subtended by reciever aperture at the target
Pb=abs((rho_bg*Hs/pi*exp(-alpha*R)+N_sky)*delta_lambda*t_r*FF*A_ifov*omega_rcv);
%power recieved after narrowband filter
%--------------------------------------------------------------------------
mt=0.5;%transmiter modulation index
mr=0.5;%reciever modulation index
em=1;%mixing efficiency
G=20;%electron gain of tube
R_max=0.243;%peak responsivity of tube(A/Watt)
gamma_mod_bg=sqrt((1-em*mr/2)^2+((em)^2*(mr^2)^2)/2);%
i_b=G*R_max*gamma_mod_bg*Pb;%%current of solar background
i_dk=G*qe*quantum_eff*gamma_mod_bg*NEPh_dk_pc*N_pix_x*N_pix_y; %dark current
i=(i_b+i_dk)*10^8;
i_back=poissrnd(i*ones(size(YY_Recieved)));
clear P_l_pk FF t_t t_r D_r rho_tgt kapa_tgt_gauss alpha rho_bg kapa_ls_gauss beta_vol_bs
clear Hs N_sky di_thermalelta_lambda delta Reciever_focal A_ifov D_obst A_r A_obst omega_rcv Pb
clear mt mr em G R_max gamma_mod_bg i_b i_dk i quantum_eff NEPh_dk_pc N_pix_x N_pix_y qe delta_lambda
% figure;
% plot(t_YY_Recieved.*1000,i_back);
% xlabel('Time (ms)');title(['backscattering noise']);
% axis tight
%
% figure
% sig=sqrt(YY_Recieved.^2+i_back.^2);
% plot(t_YY_Recieved.*1000,sig);title(['signal + atmospheric effects']);
% axis tight
%--------------------------------------------------------------------------
% noisy signal
noisy_signal=sqrt(YY_Recieved.^2+i_back.^2+noisy_sig1);
% figure;
% plot(t_YY_Recieved.*1000,noisy_signal);title(['noisy signal']);
% xlabel('Time (ms)');
% axis tight
clear i_back YY_Recieved noisy_sig1
%% LO signal
delta_t_LO=193e-9;
t_LO =delta_t_LO:1/fs:20.5e-3+delta_t_LO;
t_target2=t_LO(end);
y_LO =chirp(t_LO,f0,t_target2,f1,'linear');
Y_LO = [y_LO y_LO(end:-1:1)];
t_Y_LO=delta_t_LO:1/fs:41.001e-3+delta_t_LO;
% figure; plot(t_Y_LO.*1000,Y_LO);
YY_LO=repmat(Y_LO,1,N);
t_YY_LO=delta_t_LO:1/fs:659.1e-3+delta_t_LO ;
% plot(t_YY_LO,YY_LO)
% title('Bidirectional chirp signal ');
% xlabel('Time(s)'); ylabel('Amplitude');
clear delta_t_LO t_LO t_target2 y_LO t_Y_LO t_YY_Recieved Y_LO
%% signal mixing
mixed_signal=YY_LO.*noisy_signal;
plot(t_YY_LO,mixed_signal);
axis tight
clear YY_LO noisy_signal
  1 comentario
majid
majid el 28 de Jun. de 2021
Sorry I forgot to say,
if there is any way to do this work with low pass filter, please guide me.

Iniciar sesión para comentar.

Respuestas (0)

Categorías

Más información sobre Timing and presenting 2D and 3D stimuli en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by