%==============================================================================
%
%   Author: Imtiaz Nizami
%   Date:   11/07/2002
%
%   Time-Varying Signal Processing (EEE 598-401-B)
%   Assignment 5
%
%==============================================================================

clc;
clear all;
close all;

%=======================================
% LFM Chirp
%----------
i = sqrt(-1);
chirp_rate = 2e06;         % (Hz)^2
f_s = 102.4e03;            % set the sampling rate
T = 1/f_s;                 % sample interval
f_initial = 20e03;         % Hz initial frequency
%f_initial = f_initial/f_s;

t_final = 10e-03;          % length of time
f_final = f_initial + chirp_rate * t_final;
%f_final = f_final/f_s;
t = 0 : T : (t_final - T); % 0 to t_final secs @ f_s sample rate

signal_x = ( exp( i*2*pi * ( f_initial * t + (chirp_rate/2) * t.^2 ) ) )';
%signal_x = signal_x / norm(signal_x);
%plot(signal_x)
%pause

%=======================================
% Part 1 (i)
%-----------
[WD_signal_x, WD_t_signal_x, WD_f_signal_x] = tfrwv(signal_x);
%plot(WD_signal_x)
%pause
WD_f_signal_x = WD_f_signal_x * f_s;
imagesc(t, WD_f_signal_x, WD_signal_x)
title('Image plot - WD of LFM chirp')
xlabel('Time')
ylabel('Frequency')
pause

%=======================================
% Part 2 (i)
%-----------
[marginal_t_WD_signal_x, marginal_f_WD_signal_x, energy_WD_signal_x] = margtfr(WD_signal_x);
sum(marginal_t_WD_signal_x)
sum(marginal_f_WD_signal_x)
energy_WD_signal_x

%=======================================
% Part 3 (i)
%-----------
signal_w = noisecg(length(signal_x), 0);

%plot(signal_w)
%pause

% Part 3 (ii)
%------------
[WD_signal_w, WD_t_signal_w, WD_f_signal_w] = tfrwv(signal_w);
%plot(WD_signal_w)
%pause
WD_f_signal_w = WD_f_signal_w * f_s;
imagesc(t, WD_f_signal_w, WD_signal_w)
title('Image plot - WD of AWGN')
xlabel('Time')
ylabel('Frequency')
pause

%=======================================
% Part 4 (i)
%-----------
SNR=0;
signal_y = sigmerge(signal_x,signal_w,SNR);
Esig=mean(abs(signal_x).^2); Enoise=mean(abs(signal_y-signal_x).^2);
10*log10(Esig/Enoise)

%plot(signal_y)
%pause

%signal_temp = signal_x + signal_w;

% Part 4 (ii)
%------------
[WD_signal_y, WD_t_signal_y, WD_f_signal_y] = tfrwv(signal_y);
%plot(WD_signal_y)
%pause
WD_f_signal_y = WD_f_signal_y * f_s;
imagesc(t, WD_f_signal_y, WD_signal_y)
title('Image plot - WD of LFM chirp + AWGN')
xlabel('Time')
ylabel('Frequency')
pause

%[WD_signal_temp, WD_t_signal_temp, WD_f_signal_temp] = tfrwv(signal_temp');
%plot(WD_signal_temp)
%pause
%imagesc(WD_signal_temp)
%pause


%=======================================
% Part 5 (i)
%-----------
unitarity_signal_x_y = (abs(signal_x * signal_y'))^2;
%unitarity_WD_signal_x_y = 

% Part 5 (ii)
%------------
mesh(WD_signal_w)
title('Mesh plot - WD of AWGN, (H0)')
pause
mesh(WD_signal_y)
title('Mesh plot - WD of LFM chirp + AWGN, (H1)')
pause


% Part 5 (iii)
%-------------
matched_filter_under_H0 = sum(dot(WD_signal_x,WD_signal_w))
matched_filter_under_H1 = sum(dot(WD_signal_x,WD_signal_y))

%=======================================
% Part 6 (i)
%-----------
signal_g = chirp(t,f_initial,t_final,f_final,'linear');   % Start @ f_initial Hz, cross f_final Hz at t=t_final sec

%plot(signal_g)
%pause

% Part 6 (ii)
%------------
[WD_signal_g, WD_t_signal_g, WD_f_signal_g] = tfrwv(signal_g');
%plot(WD_signal_g)
%pause
WD_f_signal_g = WD_f_signal_g * f_s;
imagesc(t, WD_f_signal_g, WD_signal_g)
title('Image plot - WD of real chirp with 10ms duration')
xlabel('Time')
ylabel('Frequency')
pause


%=======================================
% Part 7 (i)
%-----------
t_final = 20e-03;          % length of time
f_final = f_initial + chirp_rate * t_final;
%f_final = f_final/f_s;
t = 0 : T : (t_final - T);   % 0 to t_final secs @ f_s sample rate

signal_h = chirp(t,f_initial,t_final,f_final,'linear');   % Start @ f_initial Hz, cross f_final Hz at t=t_final sec

%plot(signal_h)
%pause

% Part 7 (ii)
%------------
[WD_signal_h, WD_t_signal_h, WD_f_signal_h] = tfrwv(signal_h');
%plot(WD_signal_h)
%pause
WD_f_signal_h = WD_f_signal_h * f_s;
imagesc(t, WD_f_signal_h, WD_signal_h)
title('Image plot - WD of real chirp with 20ms duration')
xlabel('Time')
ylabel('Frequency')
pause

% Part 7 (iv)
%------------
analytic_signal_h = hilbert(signal_h);

%plot(analytic_signal_h)
%pause

[WD_analytic_signal_h, WD_t_analytic_signal_h, WD_f_analytic_signal_h] = tfrwv(analytic_signal_h');
%plot(WD_analytic_signal_h)
%pause
WD_f_analytic_signal_h = WD_f_analytic_signal_h * f_s;
imagesc(t, WD_f_analytic_signal_h, WD_analytic_signal_h)
title('Image plot - WD of analytic real chirp with 20ms duration')
xlabel('Time')
ylabel('Frequency')
pause

close all;