% Skript frac_delay_0.m, in dem ein 'Fractional Delay' Filter
% untersucht wird. 
clear;
% -------- Parameter der Simulation
N = 10;     % N+1 Koeffizienten des Frac-Delay-Filters
n = 0:N;
nk = 0:0.1:N;  % Dichtere Schrittweite für die Hüllen
D = 4.6;       % Sollte kleiner als N/2 sein 
Dg = floor(D);  % Ganze Zahl der Verspätung
% Einheitspulsantwort des Filters
h = sinc(n-D);
hk = sinc(nk-D);   % Hülle
hg = sinc(n-Dg);
hkg = sinc(nk-Dg); % Hülle
figure(1);    clf;%--------------
subplot(211),  z1 = stem(n, hg);   hold on;
set(z1,'LineWidth', 1.5);
plot(nk, hkg);   hold off;
title(['Einheitspulsantwort des ''Fractional Delay'' Filters',...
   ' für D = ',num2str(Dg)]);
xlabel('n');    grid on;
subplot(212),  z2 = stem(n, h);   hold on;
set(z2,'LineWidth', 1.5);
plot(nk, hk);   hold off;
title(['Einheitspulsantwort des ''Fractional Delay'' Filters',...
   ' für D = ',num2str(D)]);
xlabel('n');    grid on; 

% ------- Frequenzgang des Frac-Delay-Filters
nfft = 512;
[H,w] = freqz(h,1, nfft);
[g, w1] = grpdelay(h,nfft);

figure(2);    clf;%--------------
subplot(311),   plot(w/(2*pi), 20*log10(abs(H)));
title(['Amplitudengang des Fractional-Delay-Filters',...
    ' für D = ',num2str(D)]);
xlabel('f/fs');   grid on;
subplot(312),   plot(w/(2*pi), unwrap(angle(H)));
title('Phasengang des Fractional-Delay-Filters');
xlabel('f/fs');   grid on;
subplot(313),   plot(w1/(2*pi), g);
title('Gruppenlaufzeit als Anzahl Abtastwerte');
xlabel('f/fs');   grid on;
La = axis;   axis([La(1:2), 0, 2*D]);

% ------- Frequenzgang des Frac-Delay-Filters gewichtet mit Fensterfunktion
hw = h.*hamming(length(h))';
hwk = sinc(nk-D);   % Hülle
hwk = hwk.*hamming(length(hwk))';

[Hw,w] = freqz(hw,1, nfft);
[gw, w1] = grpdelay(hw,nfft);

figure(3);    clf;%--------------
subplot(411), z3 = stem(n, hw);   hold on;
set(z3,'LineWidth', 1.5);
plot(nk, hwk);   hold off;
title(['Einheitspulsantwort mit Fensterfunktion (D = ',num2str(D),')'])
xlabel(['n (D = ',num2str(D),')']);     grid on;
subplot(412),   plot(w/(2*pi), 20*log10(abs(Hw)));
title(['Amplitudengang des Fractional-Delay-Filters']);
xlabel('f/fs');   grid on;La = axis;   axis([La(1:3), 2]);

subplot(413),   plot(w/(2*pi), unwrap(angle(Hw)));
title('Phasengang des Fractional-Delay-Filters');
xlabel('f/fs');   grid on;
subplot(414),   plot(w1/(2*pi), gw);
title('Gruppenlaufzeit als Anzahl Abtastwerte');
xlabel(['f/fs  (D = ',num2str(D),')']);   grid on;
La = axis;   axis([La(1:2), 0, 1.5*D]);

% -------- Filter mit mehr Koeffizienten
N = 20;     % N+1 Koeffizienten des Frac-Delay-Filters
n = 0:N;
nk = 0:0.1:N;  % Dichtere Schrittweite für die Hüllen
D = 10.6;   % Sollte ungefähr kleiner als N/2 sein 
Dg = floor(D);  % Ganze Zahl der Verspätung
% Einheitspulsantwort des Filters
h = sinc(n-D);
hk = sinc(nk-D);   % Hülle
hg = sinc(n-Dg);
hkg = sinc(nk-Dg); % Hülle

% ------- Frequenzgang des Frac-Delay-Filters gewichtet mit Fensterfunktion
hw = h.*hamming(length(h))';
hwk = sinc(nk-D);   % Hülle
hwk = hwk.*hamming(length(hwk))';

[Hw,w] = freqz(hw,1, nfft);
[gw, w1] = grpdelay(hw,nfft);

figure(4);    clf;%--------------
subplot(411), z3 = stem(n, hw);   hold on;
set(z3,'LineWidth', 1.5);
plot(nk, hwk);   hold off;
title(['Einheitspulsantwort mit Fensterfunktion (D = ',num2str(D),')'])
xlabel(['n']);     grid on;
subplot(412),   plot(w/(2*pi), 20*log10(abs(Hw)));
title(['Amplitudengang des Fractional-Delay-Filters']);
xlabel('f/fs');   grid on; La = axis;   axis([La(1:3),2]);
subplot(413),   plot(w/(2*pi), unwrap(angle(Hw)));
title('Phasengang des Fractional-Delay-Filters');
xlabel('f/fs');   grid on;
subplot(414),   plot(w1/(2*pi), gw);
title('Gruppenlaufzeit als Anzahl Abtastwerte');
xlabel(['f/fs  (D = ',num2str(D),')']);   grid on;
La = axis;   axis([La(1:2), 0, 1.5*D]);





