close all hidden
clear
clc

set(groot, 'DefaultTextInterpreter', 'LaTeX');
set(groot, 'DefaultAxesTickLabelInterpreter', 'LaTeX');
set(groot, 'DefaultAxesFontName', 'LaTeX');
set(groot, 'DefaultAxesFontSize', 14);
set(groot, 'DefaultLegendInterpreter', 'LaTeX'); 
cmap = lines;

%**************************************************************************
% Basic parameters
%**************************************************************************
f_s = 200;        % Sampling freq for "continuous" signal in Hz
T_0 = 1;          % Fundamental period for signal description in s
T_max = 100;      % Signal duration in seconds

%**************************************************************************
% Calculate support variables
%**************************************************************************
N = T_max * f_s + 1;                     % Duration in samples
t = linspace(-T_max/2, T_max/2, N);      % Time scale in s
w = 2 * pi * linspace(-f_s/2, f_s/2, N); % Frequency scale in s^-1

%**************************************************************************
% Signal generation
%**************************************************************************

signalSel = 1;

switch(signalSel)
    
    case 1
        % Sin signal with heaviside function ------------------------------
        % Define auxiliary parameters
        sigma = 0.01;
        w_0 = 2 * pi * 1 / T_0;

        % Time domain signal and labels 
        v = sin(w_0 * t);
        v(t < 0) = 0;

        t_lim = 2*[-1 1];
        t_ticks = [-2 -1 0 1 2];
        t_labels = {'', '', '0', '$T_0$', ''};

        v_ticks = [-1 0 1];
        v_labels = {'-1', '0', '1'};

        % Area of conversion 
        s_re = [0 1 1 0]; 
        s_im = [1 1 -1 -1];

        s_re_labels = {'$-\infty$', '', '0', '', '$\infty$'};

        % Laplace domain signal and labels
        s_conv = sprintf('$\\mathrm{Re}\\{s\\} = %.2f $', sigma);

        w_lim = 2*w_0*[-1 1];
        w_ticks = w_0*[-2 -1 0 1 2];
        w_labels = {'', '$-\omega_0$', '0', '$\omega_0$', ''};

        V = w_0 ./ (sigma^2 + 1i * 2 * w + w_0^2 - w.^2);
        V = abs(V);
        V_lim = max(V) * [0 1];
        V_ticks = max(V) * [0 1];
        V_labels = {'0', '0.5'};
    
    case 2
        % Sin signal without heaviside function -------------------------------
        % Define auxiliary parameters
        w_0 = 2 * pi * 1 / T_0;

        % Time domain signal and labels 
        v = sin(w_0 * t);

        t_lim = 2*[-1 1];
        t_ticks = [-2 -1 0 1 2];
        t_labels = {'', '', '0', '$T_0$', ''};

        v_ticks = [-1 0 1];
        v_labels = {'-1', '0', '1'};

        % Area of conversion 
        s_re = [-0.01, -0.01 0.01 0.01]; 
        s_im = [-1 1 1 -1];

        s_re_labels = {'$-\infty$', '', '0', '', '$\infty$'};

        % Laplace domain signal and labels
        s_conv = sprintf('$\\mathrm{Re}\\{s\\} = 0 $');

        w_lim = 2*w_0*[-1 1];
        w_ticks = w_0*[-2 -1 0 1 2];
        w_labels = {'', '$-\omega_0$', '0', '$\omega_0$', ''};

        V = zeros(N, 1);
        % Find freq index closest to w_0
        [~, w_0_c] = min(abs(w - w_0));
        % Set value of found index to pi
        V(w_0_c) = 0.5;
        % Find freq index closest to -w_0
        [~, w_0_c] = min(abs(w + w_0));
        % Set value of found index to pi
        V(w_0_c) = 0.5;

        V = abs(V);

        V_lim = max(V) * [0 1];
        V_ticks = max(V) * [0 1];
        V_labels = {'0', '$\infty$'};
    
    case 3   
        % Decreasing exp with heaviside function --------------------------
        % Define auxiliary parameters
        alpha = -0.7;
        sigma = 0.6;
        T_h = log(0.5) / alpha;
        w_0 = 2 * pi * 1 / T_0;

        % Time domain signal and labels 
        v = exp(alpha*t); % sin(w_0 * t);
        v(t < 0) = 0;

        t_lim = T_h*2*[-1 1];
        t_ticks = T_h* [-2 -1 0 1 2];
        t_labels = {'', '', '0', '$T_h$', ''};

        v_ticks = [-1 0 1];
        v_labels = {'-1', '0', '1'};

        % Area of conversion 
        s_re = [-0.5, -0.5 1 1]; 
        s_im = [-1 1 1 -1];

        s_re_labels = {'$-\infty$', '$s_{\infty}$', '0', '', '$\infty$'};

        % Laplace domain signal and labels
        s_conv = sprintf('$\\mathrm{Re}\\{s\\} = %.1f $', sigma);

        w_lim = 4*pi/T_h * [-1 1];
        w_ticks = 2*pi/T_h * [-2 -1 0 1 2];
        w_labels = {'', '$-\frac{2 \pi}{T_h}$', '0', '$\frac{2 \pi}{T_h}$', ''};

        V = 1 ./ (alpha - sigma - 1i*w);
        V = abs(V);

        V_lim = 1.1 * max(V) * [0 1];
        V_ticks = 1.1 * max(V) * [0 1];
        V_labels = {'0', '$b$'};
  
    case 4
        % Increasing exp with heaviside function --------------------------
        % Define auxiliary parameters
        alpha = 0.7;
        sigma = 0.9;
        T_h = log(2) / alpha;
        w_0 = 2 * pi * 1 / T_0;

        % Time domain signal and labels 
        v = exp(alpha*t); % sin(w_0 * t);
        v(t < 0) = 0;

        t_lim = T_h*2*[-1 1];
        t_ticks = T_h* [-2 -1 0 1 2];
        t_labels = {'', '', '0', '$T_d$', ''};

        v_ticks = [-1 0 1];
        v_labels = {'-1', '0', '1'};

        % Area of conversion 
        s_re = [0.5, 0.5 1 1]; 
        s_im = [-1 1 1 -1];

        s_re_labels = {'$-\infty$', '', '0', '$s_{\infty}$', '$\infty$'};

        % Laplace domain signal and labels
        s_conv = sprintf('$\\mathrm{Re}\\{s\\} = %.1f $', sigma);

        w_lim = 4*pi/T_h * [-1 1];
        w_ticks = 2*pi/T_h * [-2 -1 0 1 2];
        w_labels = {'', '$-\frac{2 \pi}{T_d}$', '0', '$\frac{2 \pi}{T_d}$', ''};

        V = 1 ./ (alpha - sigma - 1i*w);
        V = abs(V);

        V_lim = 1.1 * max(V) * [0 1];
        V_ticks = 1.1 * max(V) * [0 1];
        V_labels = {'0', '$b$'};
   
end

%**************************************************************************
% Plots
%**************************************************************************
% Evaluate section by section to get plots after each other


figure
% set(gcf, 'Position', [2000 0 1400 1350]) % Set window size

f = subplot(3, 1, 1);
plot(t, v, 'LineWidth', 2, 'Color', cmap(1, :))
hold on
grid on
xlabel('$t$', 'FontSize', 20);
ylabel('$v(t)$', 'FontSize', 20);
title('Continuous time signal', 'FontSize', 20)
xlim(t_lim)
f.XTickLabelMode = 'manual';
f.XTick = t_ticks; 
f.XTickLabel = t_labels; 

f.YTickLabelMode = 'manual';
f.YTick = v_ticks;
f.YTickLabel = v_labels; 


%%
g = subplot(3, 1, 2);
fill(s_re, s_im, cmap(2, :), 'EdgeColor', 'none')
grid on
xlabel('$\mathrm{Re}\{s\}$', 'FontSize', 20);
ylabel('$\mathrm{Im}\{s\}$', 'FontSize', 20);
title('Area of conversion', 'FontSize', 20)
xlim([-1 1])
ylim([-1 1])

g.XTickLabelMode = 'manual';
g.XTick = [-1 -0.5 0 0.5 1];
g.XTickLabel = s_re_labels;

g.YTickLabelMode = 'manual';
g.YTick = [-1 0 1];
g.YTickLabel = {'$-\infty$', '0', '$\infty$'};


%%
h = subplot(3, 1, 3);
plot(w, V, 'LineWidth', 2, 'Color', cmap(3, :))
grid on
xlabel('$\mathrm{Im}\{s\}$', 'FontSize', 20);
ylabel('$|V(s)|$', 'FontSize', 20);
title(['Laplace spectrum (' s_conv ')'], 'FontSize', 20)
xlim(w_lim)
h.XTickLabelMode = 'manual';
h.XTick = w_ticks; 
h.XTickLabel = w_labels; 

ylim(V_lim)
h.YTickLabelMode = 'manual';
h.YTick = V_ticks;
h.YTickLabel = V_labels;

