Mein MATLAB Forum

% Zeitbereich
% ----------------------------------

fa = 8000; % Abtastfrequenz
fn = fa/2; % Nyquistfrequenz
N = 1024; % gewünschte FFT-Länge (N=2^x, sonst wird der DFT-Algorithmus verwendet!)
df = fa/N; % Frequenzauflösung
% Erzeugung eines Datensatzes mit N Abtastwerten
% ----------------------------------------------
t = 0 : 1/fa : (N-1)/fa; % x-Vektor
% Frequenzvorgabe in Hz als ganzzahlig Vielfaches der Frequenzauflösung der DFT/FFT:
f1 = df*100; % bei fa = 8000 Hz und N = 1024 beträgt df = 7,8125 Hz und
% f1 damit 781,25 Hz
 f1 = 784;
 f1 = df;
 phase = pi/2;

a1 = 1; % Amplitudenvorgabe
y = a1*sin(2*pi*f1*t); % y-Vektor
y = [y(1:N/2) zeros(1, N/2)];
% Graphische Darstellung
% ----------------------
% max. Amplitude zur Skalierung der graphischen Darstellung feststellen:
max_y = max(abs(y))*1.1;
fig = figure(1);
plot(y)
axis([0 N -max_y max_y])
title('Datensatz')
ylabel('Amplitude')
xlabel('N Stützstellen')
grid


% Frequenzbereich
% ----------------------------------

% Berechnung der FFT
% ------------------
H = fft(y, N);
% Berechnung des Amplitudengangs aus dem komplexen Frequenzvektor H:
amplH = abs(H);
% Amplitudenskalierung (Normierung auf N) und verschieben der Elemente des
% Amplitudenvektors, so dass die Darstellung des Amplitudengangs von -fn...0...fn
% erfolgen kann:
amplitudengang = fftshift(amplH/N);
% Graphische Darstellung
% ----------------------
% Frequenzvektoren (werden bei der graphischen Darstellung benötigt):
x_fn = 0 : df : fn-df;
x_fa = 0 : df : fa-df;
% max. Amplitude zur Skalierung der graphischen Darstellung feststellen:
%a = max([a1, a2, a3, a4, a5]); % wird später benötigt
a = a1;
fig = figure(fig+1);
stem(x_fa-fn, amplitudengang, 'b.-')
axis([-fn fn 0 a/2*1.1])
title('Amplitudengang')
ylabel('Amplitude')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid

% Ausgabe in dB
% ------------------
fig = figure(fig+1);
plot(x_fa-fn, 20*log10(amplitudengang))
%axis([-fn fn -100 20*log10(a/2)+3])
axis([-fn fn -100 3])
title('Amplitudengang')
ylabel('Amplitude in dB')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid

% Darstellung des interessierenden Frequenzbereichs des
% Amplitudengangs (0...fn) und
% daran angepasste Amplitudenskalierung (Normierung auf N/2):
amplitudengang = [amplH(1)/N amplH(2:N/2)/(N/2)]; % DC-Bin auf N normieren!
fig = figure(fig+1);
stem(x_fn, amplitudengang, 'b.-')
axis([0 fn 0 a*1.1])
title('Amplitudengang')
ylabel('Amplitude')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid

% Ausgabe in dB
% ------------------
fig = figure(fig+1);
plot(x_fn, 20*log10(amplitudengang))
axis([0 fn -100 20*log10(a)+3])
title('Amplitudengang')
ylabel('Amplitude in dB')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid


% Fensterfunktion
% ----------------------

% Anhang an die bereits erfolgte Untersuchung
% -------------------------------------------
win = hann(N)';
%y_win = y.*win; % Fensterung ohne Amplitudenkorrektur
y_win = y.*win*N/sum(win); % Fensterung mit Amplitudenkorrektur
max_y = max(abs(y_win))*1.1;

fig = figure(fig+1);
plot(y_win)
axis([0 N -max_y max_y])
title('Datensatz nach Fensterung mit Hann-Fenster')
ylabel('Amplitude')
xlabel('N Stützstellen')
grid

% Berechnung der FFT
% ------------------
H = fft(y_win, N);
% Berechnung des Amplitudengangs aus dem komplexen Frequenzvektor H:
amplH = abs(H);
% Amplitudenskalierung (Normierung auf N) und verschieben der Elemente des
% Amplitudenvektors, so dass die Darstellung des Amplitudengangs von -fn...0...fn
% erfolgen kann:
amplitudengang = fftshift(amplH/N);

% Graphische Darstellung
% ----------------------
fig = figure(fig+1);
stem(x_fa-fn, amplitudengang, 'b.-')
axis([-fn fn 0 a/2*1.1])
title('Amplitudengang nach Fensterung')
ylabel('Amplitude')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid

% Ausgabe in dB
% ------------------
fig = figure(fig+1);
plot(x_fa-fn, 20*log10(amplitudengang))
%axis([-fn fn -100 20*log10(a/2)+3])
axis([-fn fn -100 3])
title('Amplitudengang nach Fensterung')
ylabel('Amplitude in dB')
xlabel(['Auflösung: ',num2str(df),' Hz Frequenz in Hz'])
grid