calculations_only = false;

window = 512;
overlap = 64;
fft_precision = 2048;

sample_rate = 20e6; % 20 MS/s -> 20e6 S/s

channels = [1 5 9 13];

res1 = analyzeTrace("traces/2412mhz.mat", channels(1), window, overlap, fft_precision, sample_rate);
res2 = analyzeTrace("traces/2432mhz.mat", channels(2), window, overlap, fft_precision, sample_rate);
res3 = analyzeTrace("traces/2452mhz.mat", channels(3), window, overlap, fft_precision, sample_rate);
res4 = analyzeTrace("traces/2472mhz.mat", channels(4), window, overlap, fft_precision, sample_rate);

results = [res1 res2 res3 res4];

% Calibrated color scale over all traces
all_P_dB = [
    res1.P_dB(:)
    res2.P_dB(:)
    res3.P_dB(:)
    res4.P_dB(:)
];

p_dB_min = min(all_P_dB);
p_dB_max = max(all_P_dB);

% 2
if ~calculations_only
    plotSpectrogram(res1, 1, p_dB_min, p_dB_max);
    plotSpectrogram(res2, 2, p_dB_min, p_dB_max);
    plotSpectrogram(res3, 3, p_dB_min, p_dB_max);
    plotSpectrogram(res4, 4, p_dB_min, p_dB_max);
end

%{
% 3. 
% taking the 10th percentile on the linear values or the logarithmic ones
% should not make any difference

noise_floor1_linear = prctile(mag1(:), 10);
noise_floor2_linear = prctile(mag2(:), 10);
noise_floor3_linear = prctile(mag3(:), 10);
noise_floor4_linear = prctile(mag4(:), 10);

noise_floor1 = linearTodB(noise_floor1_linear);
noise_floor2 = linearTodB(noise_floor2_linear);
noise_floor3 = linearTodB(noise_floor3_linear);
noise_floor4 = linearTodB(noise_floor4_linear);

% this is false -> linear avg on logarithmic values
% noise_floor_avg = sum([noise_floor1 noise_floor2, noise_floor3, noise_floor4]) / 4;

noise_floor_avg = linearTodB(mean([noise_floor1_linear noise_floor2_linear noise_floor3_linear noise_floor4_linear]));

% 4.
occupancy1 = sum(P1_dB(:) > noise_floor1 + 10) / numel(P1_dB) * 100;
occupancy2 = sum(P2_dB(:) > noise_floor2 + 10) / numel(P2_dB) * 100;
occupancy3 = sum(P3_dB(:) > noise_floor3 + 10) / numel(P3_dB) * 100;
occupancy4 = sum(P4_dB(:) > noise_floor4 + 10) / numel(P4_dB) * 100;
%}

% 3.
% Sum all frequency bins per time slot to obtain total channel power.
noise_floor_avg = linearTodB(mean([
    res1.noise_floor_linear
    res2.noise_floor_linear
    res3.noise_floor_linear
    res4.noise_floor_linear
    ]));

% 4.
occupancies = [res1.occupancy res2.occupancy res3.occupancy res4.occupancy];

[occupancy_max, max_idx] = max(occupancies);
busiest_channel = channels(max_idx);

if ~calculations_only
    disp("Results Task 1:" + newline)

    disp("-- Trace duration ---")
    disp("Channel  1 duration: " + res1.duration + "s")
    disp("Channel  5 duration: " + res2.duration + "s")
    disp("Channel  9 duration: " + res3.duration + "s")
    disp("Channel 13 duration: " + res4.duration + "s" + newline)

    disp("--- Noise floor ---")
    disp("Distinct:")
    disp("Channel  1 noise floor: " + res1.noise_floor + "dB")
    disp("Channel  5 noise floor: " + res2.noise_floor + "dB")
    disp("Channel  9 noise floor: " + res3.noise_floor + "dB")
    disp("Channel 13 noise floor: " + res4.noise_floor + "dB" + newline)

    disp("Combined:")
    disp("Channel all noise floor: " + noise_floor_avg + "dB" + newline)

    disp("-- Occupancy ---")
    disp("Channel  1 occupancy: " + res1.occupancy + "%")
    disp("Channel  5 occupancy: " + res2.occupancy + "%")
    disp("Channel  9 occupancy: " + res3.occupancy + "%")
    disp("Channel 13 occupancy: " + res4.occupancy + "%" + newline)

    disp("Busiest channel: " + busiest_channel + newline)
end