%==========================================================================
%  Rechnerübung zur Signalübertragung II
%  Versuch 2: Nyquist-Kriterium
%==========================================================================

%--------------------------------------------------------------------------
%  Teil 2.1 Sendesignal und Ausgangssignal des Korrelationsfilters erzeugen
%--------------------------------------------------------------------------

clear all;
close all;
clc;

%  Einige Definitionen
%  Definition der Taktzeit: 
T = 1.0;            
%  Dargestellter Zeitbereich: 
t = 0:0.005:9;
%  Definition einer Beispielbinärfolge mit der Länge = 8
%  Sie können die Binärfolge ändern -->
a = [1 1 0 0 1 0 0 0];

%  Sendesignal m(t)
m = 0;
for i = 0:7
    m = m + a(i+1)*(heaviside(t-i)-heaviside(t-i-1));
end

%  m(t) plotten
figure(1)
plot(t,m,'LineWidth',2.5,'Color',[1 0 0]);
axis([0 8 -0.5 1.5]);
grid on
hold off
title('Abbildung 2: Verlauf des Sendesignals m(t)')
xlabel('t')
ylabel('m(t)')

%  Ausgangssignal des Korrelationsfilters g(t)
g = 0;
for i = 0:7
    g = g + a(i+1)*(((heaviside(t-i)-heaviside(t-i-1)).*(t-i))...
                  +((heaviside(t-i-1)-heaviside(t-i-2)).*(i+2-t)))  ;
end

%  g(t) plotten
figure(2)
plot(t,g,'LineWidth',2.5,'Color',[0 0 1]);
axis([0 8 -0.5 1.5]);
grid on
title('Abbildung 3: Ausgangssignal des Korrelationsfilters im Falle des ungestörten Kanals')
xlabel('t')
ylabel('g(t)')

%--------------------------------------------------------------------------
%  Teil 2.2 Erhöhung der Dauer der Autokorrelationsfunktion
%--------------------------------------------------------------------------

%  TD ist die Dauer der Autokorrelationsfunktion
%  Grundeinstellung TD = 1.0, hier kann die Dauer beispielsweise auf den
%  Wert 1.4 erhöht werden.
%  Sie können TD ändern -->
TD = 1.0;

Dif = TD - 1;

%  Ausgangssignal des Korrelationsfilters g(t)
figure(3); 
hold on;
g_TD = 0;
g_einzel = zeros(8,1801);
for i = 0:7
    g_einzel(i+1,:) = a(i+1)*(((heaviside(t-i+Dif)-heaviside(t-i-1)).*(t/TD-i/TD+(1-1/TD)))...
                +((heaviside(t-i-1)-heaviside(t-i-2-Dif)).*(i/TD-t/TD+(1+TD)/TD)));
    % die einzelnen Autokorrelationsfunktionen plotten
    switch i
        case 0
            plot(t,g_einzel(i+1,:),'b--','LineWidth',2);
        case 1
            plot(t,g_einzel(i+1,:),'r--','LineWidth',2);
        case 2
            plot(t,g_einzel(i+1,:),'c--','LineWidth',2);
        case 3
            plot(t,g_einzel(i+1,:),'m--','LineWidth',2);
            
    end
    
    g_TD = g_TD + a(i+1)*(((heaviside(t-i+Dif)-heaviside(t-i-1)).*(t/TD-i/TD+(1-1/TD)))...
                +((heaviside(t-i-1)-heaviside(t-i-2-Dif)).*(i/TD-t/TD+(1+TD)/TD)));
  
end

%  g(t) plotten

plot(t,g_TD,'g');
axis([0 8 -0.5 2]);
grid on
title('Abbildung 4: Ausgangssignal des Korrelationsfilter im Falle des ungestörten Kanals mit 4 einzeln eingezeichneten verschobenen Autokorrelationsfunktionen')
xlabel('t')
ylabel('Summe g(t)=grün, einzelne g(t-T)=blau,g(t-2T)=rot,g(t-3T)=cyan,g(t-4T)=magenta')

%  Abtastung der g_TD
%  Eingestellte Schwelle: C=0.5
C = 0.5;    
g_get(1,8) = 0;
for i = 1:8
    if g_TD(i*200+1) > C
        g_get(i) = 1;
    else g_get(i) = 0;
    end
end

%  Zum Vergleich sind hier die gesendeten Datenbits a[i] und die emfangenen 
%  Datenbits g_TD[i] untereinander dargestellt --> auf Command Window
disp('TD:')
disp(TD)
disp('die gesendeten Datenbits a[i]:')
disp(a)
disp('die empfangenen Datenbits g[i]:')
disp(g_get)

%  Anzahl der fehlerhaften Datenbits
%  Verändern Sie die Dauer der Autokorrelationsfunktion TD auf höhere
%  Werte!
error1 = 0;
for i = 1:8
    error1 = error1 + abs(g_get(i) - a(i));
end
disp('Anzahl der fehlerhaften Datenbits (ohne Rauschen):')
disp(error1)

%--------------------------------------------------------------------------
%  Teil 2.3 Der Einfluss des Rauschens
%--------------------------------------------------------------------------

%  Maximalwert des Rauschens
%  Sie können n_max ändern -->
n_max = 0.29;

%  Rauschsignal
n_e = 2*(n_max*rand(1,length(t))-n_max/2);

%  Summensignal am Ausgang des Empfängsfilters
y = g_TD + n_e;

% Verlauf einer Musterfunktion des Ausgangssignals y(t)
figure(4)
plot(t,y,'LineWidth',0.5,'Color',[1 0 0]);
axis([0 8 -0.5 2]);
grid on
title('Abbildung 5: Das verrauschte Ausgangssignal des Empfängsfilters')
xlabel('t')
ylabel('y(t)')

%  Schätzung der Datenbits aus den verrauschten Abtastwerten: a_en =
%  wenn(y_i>C,1,0)
y_get(1,8) = 0;
y_i(1,8) = 0;
for i = 1:8
    y_i(i) = y(i*200+1);
    if y(i*200+1) > C
        y_get(i) = 1;
    else y_get(i) = 0;
    end
end

%  Zum Vergleich sind hier die gesendeten Datenbits a[i], die verrauschten
%  Abtastwerte y[i] am Ausgang des Empfängers sowie die geschätzten
%  Datenbits a[en,i] angegeben --> auf Command Window
disp('n_max:')
disp(n_max)
disp('die gesendeten Datenbits a[i]:')
disp(a)
disp('die verrauschten Abtastwerte y[i]:')
disp(y_i)
disp('die geschätzten Datenbits a_en[i]:')
disp(y_get)

%  Anzahl der fehlerhaften Datenbits
%  Experimentieren Sie mit verschiedenen Rausch-Stärken(Parameter n_max)
%  und untersuchen Sie die Abtastwerte y[i].
%  Beachten Sie die Anzahl der fehlerhaften Datenbits.
error2 = 0;
for i = 1:8
    error2 = error2 + abs(y_get(i) - a(i));
end
disp('Anzahl der fehlerhaften Datenbits (mit Rauschen):')
disp(error2)