clear all;
close all;
clc;


 fs=40/7*(10^6); % 采樣率
 CA_freq= 1.023*(10^6);
 F_if = 1.405*(10^6);

Wnn = 2*CA_freq;;
stop_width = Wnn/50;
freq_pts=2^15; 
Wp=[(F_if-Wnn/2) (F_if+Wnn/2)]*2/fs; 
Ws=[(F_if-Wnn/2-stop_width) (F_if+Wnn/2+stop_width)]*2/fs;
Rp = 0.05; Rs = 60;
[n Wn]=ellipord(Wp,Ws,Rp,Rs);
[b_nz,a_nz]=ellip(n,Rp,Rs,Wn);

[H,f,units]=freqz(b_nz,a_nz,freq_pts,fs);
freqzplot(H,f,units)
title('Noise Filter Frequency Response')
signal = sqrt(100)*randn(1,5000);
filtered_signal = filter(b_nz,a_nz,signal);
[Pxx,f]=periodogram(filtered_signal,[],1024,fs);

plot(f,10*log10(Pxx));

axis([0 f(end) -90 -30]);

grid;
xlabel('Freqency (Hz)');
ylabel('Power Spectral Density (dB/Hz)');

Noise_out = (abs(fft(filtered_signal)));
figure,

L = length(Signal_out);

semilogy([-L/2:L/2-1]/L*fs*2,fftshift((Signal_out))/max(Signal_out));
%axis([0*(10^6) 25*(10^6) -1 1])
title('noise 功率譜');
grid on