53 lines
2.5 KiB
Matlab
Executable File
53 lines
2.5 KiB
Matlab
Executable File
% ##################################################################################
|
|
% ## Loesung: Schaetzung des Leistungsdichtespektrums nach Bartlett und Welch ##
|
|
% ## -------------------------------------------------------------------------- ##
|
|
% ## Benoetigte(s) m-File(s): lwelch.m ##
|
|
% ##################################################################################
|
|
|
|
N = 2^10; N_vor=500; n = randn(1,N+N_vor); % Modellrauschprozesse
|
|
h1 = [1 -4 6 -4 +1]; x1_all = filter(h1,1,n); x1=x1_all(N_vor+1:length(x1_all));
|
|
h2 = [1 0.8]; x2_all = filter(1,h2,n); x2=x2_all(N_vor+1:length(x2_all));
|
|
NFFT = 2^10; % "Wahres" Leistungsdichtespektrum zum Vergleich
|
|
Sx1x1 = abs(fft(h1,NFFT)).^2; axis_Sx1=[0 1 0 1.5*max(Sx1x1)];
|
|
Sx2x2 = ones(1,NFFT)./abs(fft(h2,NFFT)).^2; axis_Sx2=[0 1 0 1.5*max(Sx2x2)];
|
|
fT = 0:1/NFFT:1-1/NFFT;
|
|
v_K = [1 4 16 64];
|
|
|
|
win='boxcar'; % BARTLETT fuer x1
|
|
for k=v_K
|
|
S_xx=lwelch(x1,k,win,NFFT);
|
|
figure; plot(fT,S_xx); hold on; ylabel('Sxx(exp(j*Omega))');
|
|
plot(fT,Sx1x1,'--'); axis(axis_Sx1); xlabel('Omega/2pi = f T');
|
|
title(sprintf('Bartlett-Sch. von Sx1x1(exp(j*Om.)); K=%d und L=%d/K',k, N));
|
|
text(0.73,900,'== Periodogramm');
|
|
end;
|
|
input('......zum Fortfahren: RETURN druecken');
|
|
|
|
win='hamming'; % WELCH-Verf. mit Hamming fuer x1
|
|
for k=v_K
|
|
S_xx=lwelch(x1,k,win,NFFT);
|
|
figure; plot(fT,S_xx); hold on; ylabel('Sxx(exp(j*Omega))');
|
|
plot(fT,Sx1x1,'--'); axis(axis_Sx1); xlabel('Omega/2pi = f T');
|
|
title(sprintf('Bartlett-Sch. von Sx1x1(exp(j*Om.)); K=%d und L=%d/K',k, N));
|
|
end;
|
|
input('......zum Fortfahren: RETURN druecken');
|
|
|
|
win='boxcar'; % BARTLETT fuer x2
|
|
for k=v_K
|
|
S_xx=lwelch(x2,k,win,NFFT);
|
|
figure; plot(fT,S_xx); hold on; ylabel('Sxx(exp(j*Omega))');
|
|
plot(fT,Sx2x2,'--'); axis(axis_Sx2); xlabel('Omega/2pi = f T');
|
|
title(sprintf('Bartlett-Sch. von Sx2x2(exp(j*Om.)); K=%d und L=%d/K',k, N));
|
|
text(0.73,900,'== Periodogramm');
|
|
end;
|
|
input('......zum Fortfahren: RETURN druecken');
|
|
|
|
win='hamming'; % WELCH-Verf. mit Hamming fuer x2
|
|
for k=v_K
|
|
S_xx=lwelch(x2,k,win,NFFT);
|
|
figure; plot(fT,S_xx); hold on; ylabel('Sxx(exp(j*Omega))');
|
|
plot(fT,Sx2x2,'--'); axis(axis_Sx2); xlabel('Omega/2pi = f T');
|
|
title(sprintf('Bartlett-Sch. von Sx2x2(exp(j*Om.)); K=%d und L=%d/K',k, N));
|
|
end;
|
|
% ##### EOF
|