%% Ejemplo de canal de 2 rayos (two ray channel)
% Copyright (c) 2024 G.T.A.S. (Grupo Tratamiento Avanzado de Señal). Universidad de Cantabria.
%
% Comprobamos el efecto de un canal en el que además del rayo directo aparece una componente multicamino secundaria (por ejemplo, un eco debido a una reflexión) retardada \tau segundos y con amplitud \alpha.
% Recuerda que un canal de este tipo tendrá una respuesta impulsiva, h_c(t), y una respuesta frecuencial, |H_c(f)|^2, dadas por las expresiones:
% h_c(t)=\delta(t)+\alpha\delta(t-\tau)
% |H_c(f)|^2=(1+\alpha^2)+2\alpha\cos(2\pi f\tau)
% Vamos a observar cómo varía la respuesta en frecuencia de este canal tan simple en función del retardo y la amplitud del eco multicamino.
% En el siguiente código, modifica el valor de las variables de matlab alfa y tau a tu gusto y comprueba la respuesta frecuencial en la gráfica.

alfa = 0.2;   % Amplitud del rayo secundario |alfa|<=1
tau = 20;     % Retardo del rayo secundario

% Construimos la respuesta impulsiva del canal, asumiendo una frecuencia de muestreo f_s=1 Hertz.

hc =[];
hc(1) = 1;
hc(tau) = alfa;

% Usando Matlab calculamos la respuesta en frecuencia y la dibujamos:

[Hc,f] = freqz(hc,1,8*1024,1);
figure(1);plot(f,abs(Hc).^2);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Respuesta en frecuencia del canal');axis([0 0.5 0 4])

% Evidentemente tiene una respuesta en frecuencia no plana. Pero este canal será o no será selectivo frecuencialmente en función del ancho de banda de la señal de comunicaciones que transmitamos a través de él. Veámoslo con dos ejemplos:
% Sistema 1: Sistema con un ancho de banda muy pequeño. Por ejemplo una milésima parte de la frecuencia de muestreo:W_1=f_s/1000
% Sistema 2: Sistema con un ancho de banda relativamente grande. Por ejemplo un 10% de la frecuencia de muestreo: W_2=f_s/10
% Podemos observar con matlab en detalle la respuesta en frecuencia del canal para ambos anchos de banda (consideramos que la frecuencia central es f_c=f_s/4.

W1 = 0.0001;      % Ancho de banda: banda estrecha 
W2 = 0.1;         % Ancho de banda: banda ancha
fc = 0.25;        % Frecuencia portadora
figure(2);
subplot(121);plot(f,abs(Hc).^2);axis([fc-W1/2 fc+W1/2 0 4]);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Detalle de banda estrecha')
subplot(122);plot(f,abs(Hc).^2);axis([fc-W2/2 fc+W2/2 0 4]);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Detalle de banda ancha')

% Comprueba cómo el canal es selectivo en frecuencia o no en función del ancho de banda considerado.
%% Prueba con otros parámetros del canal
% Puedes probar a variar los valores de \alpha y \tau en el inicio del código y ejecutar de nuevo esta sección. Intenta comprender en este caso la relación entre el ancho de banda de coherencia del canal y el valor de \tau.
%
%% Canal de 2 rayos variante en el tiempo
% Repetimos el ejemplo anterior pero haciendo que los valores de  \alpha y \tau sean variables aleatorias, es decir, cambien con el tiempo (canal variante). Para ello simplemente generamos varias realizaciones de las correspondientes variables alfa y tau en el código matlab anterior.
% Ejecuta la siguiente sección de código y se desplegará una figura de Matlab en la que observar la respuesta en frecuencia de cada una de las N realizaciones del canal. Si deseas que la animación sea más lenta incrementa el argumento de la función pause (o no le pases argumento, en cuyo caso espera hasta que pulsas una tecla).

N = 200;    % Número de realizaciones del canal    
for n=1:N  
    alfa = rand(1,1);          % v.a. uniforme en [0,1]
    tau = 2+randi([0,18],1);    % v.a. uniforme en [2,20]
    hc = [];
    hc(1) = 1;
    hc(tau) = alfa;
    [Hc,f] = freqz(hc,1,8*1024,1);
    handle = figure(3);
    subplot(211);plot(f,abs(Hc).^2);axis([0 0.5 0 4]);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Respuesta en frecuencia del canal')
    subplot(223);plot(f,abs(Hc).^2);plot(f,abs(Hc).^2);axis([fc-W1/2 fc+W1/2 0 4]);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Detalle de banda estrecha')
    subplot(224);plot(f,abs(Hc).^2);axis([fc-W2/2 fc+W2/2 0 4]);grid;xlabel('frecuencia/f_s');ylabel('|H_c(f)|^2');title('Detalle de banda ancha')
    set(handle,'Visible','on')
    pause(0.2)
    %pause
end

% Observa la variabilidad del canal.
% ¿Te has fijado qué ocurre cuando nuestro sistema es banda estrecha? Un canal de dos rayos variante en el tiempo como el propuesto provoca desvanecimiento plano (flat fading).
%
%% Prueba con otras distribuciones
% Puedes probar a cambiar la estadística del canal, por ejemplo haciendo que \alpha o \tau sigan otro tipo de distribuciones. Cuidado, asegúrate que \tau>0.