傅里叶变换中的假频**

所有波都可由间谐波叠加而来，傅里叶变换实质上是对某一复杂波的组成成分进行分析，即将复杂波分离为不同频率的间谐波：
这里以一简单信号作为示例：

图一信号可由图二、图三与图四叠加而来，图二、三、四均为简单间谐波，其频率分别为：20Hz、60Hz、85Hz。傅里叶变换即为叠加过程的逆过程。
在对某一信号做傅里叶变换时，需要格外注意采样间距dt,采样频率fn=1/dt，以及奈奎斯特频率（折叠频率）fs=fn/2，为了防止假频的出现，采样定理规定采样频率需大于信号最大频率的2倍，即fn>2fmax；而当采样频率 fn<2fmax时，假频可由fa=|mfn-f|求得，其中m为整数，这相当于将超出奈奎斯特频率的部分折叠回[-fs,fs]区间，导致假频的出现。下面上述信号作为例子：

上图第一张图为采样间距dt=1ms时对应频谱，此时其折叠频率fs为1/2*dt,即500Hz,此时无假频出现；
第二张图为采样间距dt=10ms时对应频谱，此时其折叠频率fs为50Hz，可观察到60Hz信号此时对应频率60-100=-40Hz,85Hz信号此时对应频率85-100=-15Hz;
第三张图为采样间距dt=20ms时对应频谱，此时其折叠频率fs为25Hz，可观察到60Hz信号此时对应频率60-50=10Hz,85Hz信号此时对应频率85-2*50=-15Hz。
以下为Matlab代码：

dt1=0.001;
dt2=0.01;
dt3=0.02;
f1=20;
f2=60;
f3=85;
t_start=-0.5;
t_end=0.5;
df=1/(t_end-t_start); %频率间隔
t1=t_start:dt1:t_end;
t2=t_start:dt2:t_end;
t3=t_start:dt3:t_end;

f1_real=zeros(1,(t_end-t_start)/dt1);
f1_image=f1_real;
f1_a=f1_real;

f2_real=zeros(1,(t_end-t_start)/dt2);
f2_image=f2_real;
f2_a=f2_real;

f3_real=zeros(1,(t_end-t_start)/dt3);
f3_image=f3_real;
f3_a=f3_real;
source1=sin(2pif1t1+pi/6)+sin(2pif2t1+pi/4)+sin(2pif3t1+pi/5);
source2=sin(2pif1t2+pi/6)+sin(2pif2t2+pi/4)+sin(2pif3t2+pi/5);
source3=sin(2pif1t3+pi/6)+sin(2pif2t3+pi/4)+sin(2pif3t3+pi/5);
% plot(t,source);
for j=1:df:1/dt1
for i=t_start:dt1:t_end-dt1
z=floor((i-t_start)/dt1+1);
f1_real(j)=f1_real(j)+cos(2piji)source1(z);
f1_image(j)=f1_image(j)+sin(2piji)*source1(z);
end
f1_a(j)=sqrt(f1_real(j)*f1_real(j)+f1_image(j)*f1_image(j));
end

for j=1:df:1/dt2
for i=t_start:dt2:t_end-dt2
z=floor((i-t_start)/dt2+1);
f2_real(j)=f2_real(j)+cos(2piji)source2(z);
f2_image(j)=f2_image(j)+sin(2pij*i)*source2(z);
end
f2_a(j)=sqrt(f2_real(j)*f2_real(j)+f2_image(j)*f2_image(j));
end

for j=1:df:1/dt3
for i=t_start:dt3:t_end-dt3
z=floor((i-t_start)/dt3+1);
f3_real(j)=f3_real(j)+cos(2piji)source3(z);
f3_image(j)=f3_image(j)+sin(2pij*i)*source3(z);
end
f3_a(j)=sqrt(f3_real(j)*f3_real(j)+f3_image(j)*f3_image(j));
end

% f2=fft(source2);
figure(1)
subplot(3,1,1);
% plot(-1/dt1/2:df:1/dt1/2-1,fftshift(f1_a));
plot(-1/dt1/2:df:1/dt1/2,abs(fftshift(fft(source1))));
subplot(3,1,2);
plot(-1/dt2/2:df:1/dt2/2,abs(fftshift(fft(source2))));
% plot(-1/dt2/2:df:1/dt2/2-1,abs(fftshift(f2_a)));
subplot(3,1,3)
plot(-1/dt3/2:df:1/dt3/2,abs(fftshift(fft(source3))));
% plot(-1/dt3/2:df:1/dt3/2-1,abs(fftshift(f3_a)));

实际上该代码并没有利用sin与cos做加法的部分，我对比了求和与直接利用fft函数求出的结果，发现使用求和得出的结果与f使用ft函数得出结果基本一致，但频谱波折很大，暂时未找到原因，就先挂在这等以后再找一找问题。