第8章 变换域滤波器的FPGA实现
原理:
变换域滤波:
最常用的是频域处理,即将时域信号通过快速傅里叶变换(Fast Fourier Transform,FFT)变换到频域进行滤波处理,再将频域信号通过快速傅里叶反变换(Inverse Fast Fourier Transform,IFFT)变换成时域信号,从而完成滤波处理。
DFT存在的问题:
栅栏效应和序列补零,频谱泄漏和混叠失真,频率分辨率与DFT参数的选择
DFT和FFT算法计算量:
DFT算法,乘法次数与N^2成正比
FFT算法,(N / 2)log2N次复乘运算
变换域滤波技术:
变换域技术通过合适的变换把扩频信号、热噪声以及干扰组成的混合信号映射到另一个域处理,一般将NBI映射为类似于冲激的函数,将有用信号映射为与干扰正交,具有平坦特性的波形。这样就能在较彻底去除干扰的同时对有用信号的损伤较小。
MATLAB实现:
FFT算法MATLAB仿真:
仿真FFT参数对采用FFT算法分析信号频谱的影响。产生频率分别为f 1 =2 Hz、f 2 =2.05 Hz的正弦波合成信号,采样频率f s =10 Hz。根据式(8-8),要实现分辨两个单频信号的目的,DFT的序列长度必须满足N> 400。分别仿真3种情况下的FFT变换:(1)取x (n )的128点数据,计算FFT;(2)将128点x (n )以补零的方式加长到512点,计算FFT;(3)取512点x (n ),计算FFT。
f1=2; f2=2.05; %单正弦信号的频率
fs=10; %采样频率
%对 128 点时域序列进行 FFT 分析
N=128; %FFT 分析的点数
n=0:N-1;
xn1=sin(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs); %产生 128 点时域信号序列
XK1=fft(xn1); %进行傅里叶变换,并进行归一化处理
MXK1=abs(XK1(1:N/2));
%对补零后的 512 点时域序列进行 FFT 分析
M=512;
xn2=[xn1 zeros(1,M-N)]; %在时域信号序列后补零
XK2=fft(xn2); %进行傅里叶变换,并进行归一化处理
MXK2=abs(XK2(1:M/2));
%对 512 点时域序列进行 FFT 分析
n=0:M-1;
xn3=sin(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs); %产生 128 点时域信号序列
XK3=fft(xn3); %进行傅里叶变换,并进行归一化处理
MXK3=abs(XK3(1:M/2));
%绘图
subplot(321);
x1=0:N-1;
plot(x1,xn1);xlabel('n','fontsize',8);title('128 点 x(n)','fontsize',8);
subplot(322);
k1=(0:N/2-1)*fs/N;
plot(k1,MXK1);xlabel('f(Hz)','fontsize',8);title('128 点 xn 的 FFT 变换','fontsize',8);
subplot(323);
x2=0:M-1;
plot(x2,xn2);xlabel('n','fontsize',8);title('512 点补零 x(n)','fontsize',8);
subplot(324);
k2=(0:M/2-1)*fs/M;
plot(k2,MXK2);xlabel('f(Hz)','fontsize',8);title('512 点补零 xn 的 FFT 变换','fontsize',8);
subplot(325);
plot(x2,xn3);xlabel('n','fontsize',8);title('512 点 x(n)','fontsize',8);
subplot(326);
plot(k2,MXK3);xlabel('f(Hz)','fontsize',8);title('512 点 xn 的 FFT 变换','fontsize',8);
频域滤波器MATLAB仿真:
假设在伪码码率为1 MHz的扩频信号中有一个频率为100 kHz的单音干扰,且信干比为0 dB,信噪比为10 dB。采用256点FFT/IFFT变换,窗函数为汉明窗,50%重叠,门限因子θ =20。需要仿真出滤波前后信号的频谱及时域信号情况,同时需要将滤波器的输入数据量化为12比特数据,存放在外部文本文件中,作为后续FPGA仿真测试的输入源。
Fpoint=256;%fft 变换的点数为 256
poly=[1 0 0 0 1 1 1 1]; %产生 PN 码的多项式
%根据多项式产生伪码,同时对伪码进行 4 倍采样
%E8_2_PNCodeGenerator 为自定义的函数,请在光盘中查看函数代码
pn1=E8_2_PNCodeGenerator(poly);
window1=hamming(Fpoint)'; %产生 Hamming 窗
%生成相同长度的伪码 PN 及窗函数数据 Window
pn=pn1;
window=window1;
for i=1:9
pn=[pn,pn1];
end
PN=pn(1:floor(length(pn)/Fpoint)*Fpoint);
N=length(PN);
for i=1:N/Fpoint-1
window=[window,window1];
end
%产生功率为 1 W,频率为 100 kHz 的单音干扰信号
f=100000; %信号频率为 100 kHz
Fs=4*10^6; %采样频率为 4 MHz
t=0:1/Fs:(N-1)/Fs;
jam=sqrt(2)*sin(2*pi*f*t);
S=PN+jam; %产生 SJR=0 dB 的测试信号
S2=[S(Fpoint/2+1:N),S(1:Fpoint/2)]; %对输入信号进行 50%重叠
WS=S.*window; %加窗处理
WS2=S2.*window; %加窗处理
%对两路加窗信号进行变换域滤波处理
Fout=zeros(1,N);
for n=1:N/Fpoint
FWS=fft(WS((n-1)*Fpoint+1:n*Fpoint));
FWS2=fft(WS2((n-1)*Fpoint+1:n*Fpoint));
AWS=abs(FWS).*abs(FWS);
AWS2=abs(FWS2).*abs(FWS2);
MAWS=sum(AWS)/Fpoint
MAWS2=sum(AWS2)/Fpoint;
for i=1:Fpoint
if AWS(i)>(20*MAWS)
FWS(i)=0;
end;
if AWS2(i)>(20*MAWS2)
FWS2(i)=0;
end;
end;
iS(((n-1)*Fpoint+1):n*Fpoint)=ifft(FWS);
iS2(((n-1)*Fpoint+1):n*Fpoint)=ifft(FWS2);
end
%对滤波后的两路数据进行重叠处理后输出
Fout=iS+[iS2(N-Fpoint/2+1:N),iS2(1:N-Fpoint/2)];
%绘图
Fs=fft(S,256);as=abs(Fs);
Fsout=fft(Fout,256);asout=abs(Fsout);
figure(1); %绘制滤波前后信号的频谱图
subplot(211);plot(abs(as));legend('叠加干扰的原始数据频谱(sjr=0dB)');
subplot(212);plot(asout);legend('滤波后的原始数据频谱');
figure(2); %绘制滤波前后信号的时域波形图
p=1:200;
subplot(311);plot(p,PN(1:200)>0,'-',p,S(1:200)>0,'--',p,Fout(1:200)>0,'+');
axis([0 200 -0.2 1.2]);title('判决整形后的数据');xlabel('point');
legend('原始数据','滤波前数据','滤波后数据');
subplot(312);plot(S(1:200));legend('叠加干扰的原始数据');xlabel('point');
subplot(313);plot(Fout(1:200));legend('滤波后的原始数据');xlabel('point');
%将输入信号进行量化,并写入外部文本文件中,供 FPGA 测试仿真使用
FPGA实现:
直接的变换域滤波器结构:
基于模块复用的变换域滤波器结构:
第9章 解调系统滤波器的FPGA实现
数字接收机:
Costas环: