% Skript abtast_fir_2.m, in dem ein FIR-Tiefpassfilter
% über die Abtastung des Frequenzgangs entwickelt wird
% Im Übergangsbereich wird ein Zwischenwert eingefügt

clear;
% ------- Parameter
fs = 1000;      Ts = 1/fs;  % Abtastfrequnez, Abtastperiode
fp = 200;       % Durchlassfrequenz;
N = 64;         % Ordnung des Filters (N+1 Filter-Koeffizienten)
n = 0:N;        % Abtaststützstellen
n1 = round(1.2*fp/(fs/(N)));
Hzw = 0.5;      % Zwischenwert im Übergangsbereich
Hid_2 = [ones(1,n1-1),Hzw,zeros(1,N/2-n1+1)];  % mit Zwischenwert
Hid = [Hid_2, fliplr(Hid_2(2:end))];
% ------- Ideale begrenzte Einheitspulsantwort aus der inverse DFT
hid = (ifft(Hid));    % Einheitspulsantwort (muss reell sein)
hid = fftshift(hid);
nfft = 512;
[H,w] = freqz(hid,1,nfft, 'whole');
figure(1);    clf;
subplot(311), stem(n/(N+1), Hid);
title('Idealer abgetasteter Frequenzgang');
xlabel('f/fs');    grid on;
hold on; plot((0:nfft-1)/nfft, abs(H));   hold off;
subplot(312), stem(-N/2:N/2, hid);
title(['Einheitspulsantwort des BP-Filters über die ifft ']);
grid on; axis tight;
subplot(313), stem(n/(N+1), Hid);
title('Idealer abgetasteter Frequenzgang (Ausschnitt)');
xlabel('f/fs');    grid on;
hold on; plot((0:nfft-1)/nfft, abs(H));   hold off;
La = axis;   axis([0, 0.3, La(3:4)]);

% ------- Gewichten der symmetrischen hid mit Fenster-Funktion
hidn = hid.*hamming(N+1)';
% ------- Frequenzgang mit freqz
nfft = 512;            % nfft >> N
[Hn,w] = freqz(hidn, 1, nfft, 'whole');
figure(2);     clf;
subplot(121), stem(0:N, hidn);
title('hid gewichtet mit Hamming-Fenster');
xlabel('n');    grid on;
axis tight;
subplot(222), plot(w/(2*pi), abs(Hn));
title('Amplitudengang des FIR-TP');
xlabel('f/fs');    grid on;
subplot(224), plot(w/(2*pi), unwrap(angle(Hn)));
title('Phasengang des FIR-TP');
xlabel('f/fs');    grid on;

% -------- Darstellung der Amplitudengänge in dB
figure(3);    clf;
plot(w/(2*pi),[20*log10(abs(H)), 20*log10(abs(Hn))]);
title('Amplitudengang des FIR-Filters ohne und mit Fenster-Funktion');
xlabel('f/fs');     ylabel('dB');    grid on;
La = axis;    axis([La(1:2), -120, 10]);

% -------- Abbildung für das Buch
figure(4);  clf;
subplot(211), stem(n/(N+1), Hid);
title('Idealer abgetasteter Frequenzgang');
xlabel('f/fs');    grid on;
hold on; plot((0:nfft-1)/nfft, abs(H));   hold off;

subplot(212);
plot(w/(2*pi),[20*log10(abs(H)), 20*log10(abs(Hn))]);
title('Amplitudengang des FIR-Filters mit Zwischenwert ohne und mit Fenster-Funktion');
xlabel('f/fs');     ylabel('dB');    grid on;
La = axis;    axis([La(1:2), -120, 10]);
 
