clear all;										% 清所有的變量
clf;											% 清屏
N=1224;										% 設置數據點數
f1=200;										% 設置信號1頻率
f2=700;										% 設置信號2頻率
f3=1100;										% 設置信號3頻率
dt=1/1024;										% 設置采樣時間
sample=1/dt;									% 設置采樣頻率
randn('state',sum(100*clock));						% 產生隨機噪聲
dis=2*rand(1,N);
dis=dis-mean(dis);								% 將隨機噪聲的均值設置成0
dis=2*dis;										% 擴大隨機噪聲的幅度
for k=1:N;										% 產生待檢測的信號
    w(k)=1*sin(2*pi*f1*k*dt)+1*sin(2*pi*f2*k*dt)+1*sin(2*pi*f3*k*dt)+dis(k);
end

figure(1);										% 計算待檢測信號的均勻采樣頻譜
fw=fft(w);
s=abs(fw).^2;
f=(1:N-1)*sample/N;
plot(f,s(2:N));
xlabel('Frequency /Hz');							% 橫坐標說明
ylabel('Amplitude');								% 縱坐標說明
title('DFT of uniformly sampling signal');				% 標題說明

t=tnonunif(N,1,1);								% 調用非均勻采樣時間函數
for k=1:N;										% 產生非均勻采樣數據
    w2(k)=1*sin(2*pi*f1*t(k)*dt)+1*sin(2*pi*f2*t(k)*dt)+1*sin(2*pi*f3*t(k)*dt)+dis(k);
end

for k=1:N;										% 計算非均勻采樣數據的頻譜
    s2(k)=0;
    for n=1:N;									% 計算非均勻采樣數據的傅立葉變換
        s2(k)=s2(k)+w2(n)*(cos(2*pi*t(n)*(k-1)/N)-j*sin(2*pi*t(n)*(k-1)/N));
    end
    s2(k)=abs(s2(k)).^2;
end

figure(2);										% 畫出非均勻采樣數據的頻譜
f=(1:N-1)*sample/N;
plot(f,s2(2:N));
xlabel('Frequency /Hz');
ylabel('Amplitude');
title('DFT of additive random sampling signal');

f=(1:N-1)*sample/N;
plot(f,s2(2:N));
xlabel('Frequency /Hz');
ylabel('Amplitude');
title('DFT of additive random sampling signal');