%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% .m script of web demo "Insert title here"
% Written for GNU Octave
%
% INSTITUTE OF TELECOMMUNICATIONS
% University of Stuttgart
% www.inue-uni-stuttgart.de
% author: Marvin Geiselhart
% date: 26.05.2020
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Caching documentation:
% We save all curves for a real-valued SNR_b (Eb/sigma²) from 0 dB to 15 dB in .5 dB steps, 31 values
% we do not save the snr points!
% For each of the four error rates (Hard/Soft FER/BER) we save the error rate (0=NaN) and the number of errors collected (to estimate accuracy)
% If someone simulates a code in a region with more errors than cached, we update the cache
% Directory: "error_curves_cache/"
% Filenames: %codename%_hfer.m, %codename%_sber.m, %codename%_sber_errors.m, etc


function [] = ml_decoding_cache_viewer(code1, code2, code3, show_uncoded, show_fer, show_hd, snr_def, snr_start, snr_end, y_lim, filename)
  
  snr_b_real_db_points = 0:0.5:15 % do not change this
  
  % define formats
  landscape = "-S930,350";
  portrait  = "-S640,480";
  
  linewidth = 3;
  
  % select format
  output_format = portrait;
  
  % Create an invisible figure.
  fh = figure(1);
  set(fh, "visible", "off");
  
  % define style
  set(0, "defaultlinelinewidth", 2);
  set(0, "defaultaxesfontsize", 8);
  set(0, "defaultaxesfontname", "Helvetica");
  set(0, "defaulttextfontsize", 8);
  set(0, "defaulttextfontname", "Helvetica");
  
  % This to test and demonstrate the debug engine:
  %fprintf("This is a debug message. Checking format: %s", output_format);
  
  % swap limits if out of order
  if snr_start > snr_end
    tmp = snr_start;
    snr_start = snr_end;
    snr_end = tmp;
  endif
  
  snr_start_p = snr_start;
  snr_end_p = snr_end;
  clf;
  hold on;
  xlim([snr_start,snr_end]);
  
  % setup legend
  legendCounter = 1;
  legend_entries = {};
  if show_uncoded
    % uncoded
    x_vals = linspace(snr_start, snr_end, 31);
    switch (snr_def)
      case {"eb_real", "es_real"}
      y = Perr(x_vals);
      case {"eb_complex", "es_complex"}
      y = Perr(x_vals+dB(2));
    endswitch
    semilogy(x_vals, y, "color", 'k','linewidth',linewidth);
    legend_entries{legendCounter} = "Uncoded BPSK";
    legendCounter = legendCounter + 1;
  endif
  
  
  for i = 1:3
    % select a code
    if (i == 1)
      code = code1
      color = 'r'
      elseif (i == 2)
      code = code2
      color = [0,0.7,0]
      elseif (i == 3)
      code = code3
      color = 'b'
    endif
    
    
    use_cached = false;
    switch (code)
      case "rep"
      K=1;
      N=2;
      code_name = "Repetition Code"
      
      case "four_two"
      K=2;
      N=4;
      code_name = "(4,2) Code"
      
      case "six_three"
      K=3;
      N=6;
      code_name = "(6,3) Code"

      case "ext_hamming"
      code_name = "(8,4) Ext. Hamming Code"
      K=4;
      N=8;

      case "hamming_seven"
      code_name = "(7,4) Hamming Code"
      K=4;
      N=7;

      case "hamming_fifteen"
      code_name = "(15,11) Hamming Code"
      K=11;
      N=15;

      case "rm_one_four"
      code_name = "(16,5) Reed-Muller Code"
      K=5;
      N=16;
      case "polar"
      code_name = "(16,8) Polar Code"
      K=8;
      N=16;
      case "bluetooth"
      code_name = "(15,10) Bluetooth Code"
      K=10;
      N=15;
      case "ldpc"
      code_name = "(16,8) LDPC Code"
      K=8;
      N=16;
      case "egolay"
      code_name = "(24,12) Ext. Golay Code"
      K=12;
      N=24;
      case "bch"
      code_name = "(31,16) BCH Code"
      K = 16;
      N = 31;
      load("error_curves/bch_snr_b_real_cached.mat"); % snr_vals_db
      snr_b_real_cached = snr_vals_db
      load("error_curves/bch_sber.mat");
      load("error_curves/bch_hber.mat");
      load("error_curves/bch_sfer.mat");
      load("error_curves/bch_hfer.mat");
      use_cached = true;
      case "sixteen_eight_five"
      code_name = "(16,8) Distance-5 Code"
      K=8;
      N=16;
      case "NR"
    code_name = "(16,8,6) nonlinear NR-code"
    K = 8;
    N = 16;
    case "CCSDS"
    code_name = "(32,16) LDPC Code"
      K = 16;
      N = 32;
    otherwise
      continue
    endswitch
    
    Rc = K/N;
    
    % Map axis to correct region
    switch (snr_def)
      case "eb_real"
      offset = 0;
      case "es_real"
      offset = dB(Rc);
      case "eb_complex"
      offset = -dB(2);
      case "es_complex"
      offset = dB(Rc)-dB(2);
    endswitch
    
    x_vals = snr_b_real_db_points+offset;
    num_points = length(snr_b_real_db_points);
    start_index = max(floor((snr_start_p-offset)*2), 0)+1
    end_index = min(ceil((snr_end_p-offset)*2), 2*15)+1
    
    % PLOT
    if sum(snr_start_p>=x_vals)
      snr_start = max(x_vals(snr_start_p>=x_vals));
    else
      snr_start = snr_start_p;
    endif
    
    if sum(snr_end_p<=x_vals)
      snr_end = min(x_vals(snr_end_p<=x_vals));
    else
      snr_end = snr_end_p;
    endif
    
    if show_hd && show_fer && is_cached(code, 'hfer')
      hard_frame_error_rate = load_cached(code, 'hfer');
      semilogy(x_vals(snr_start<=x_vals & x_vals <= snr_end), hard_frame_error_rate(snr_start<=x_vals & x_vals <= snr_end) , "color", color, "linestyle", '-.',"marker", 'x', 'linewidth',linewidth, "markersize",3);
      legend_entries{legendCounter} = [code_name ', Hard Dec. FER'];
      legendCounter = legendCounter + 1;
    endif
    
    if show_hd && is_cached(code, 'hber')
      hard_bit_error_rate = load_cached(code, 'hber');
      semilogy(x_vals(snr_start<=x_vals & x_vals <= snr_end), hard_bit_error_rate(snr_start<=x_vals & x_vals <= snr_end) , "color", color, "linestyle", ':','linewidth',linewidth);
      legend_entries{legendCounter} = [code_name ', Hard Dec. BER'];
      legendCounter = legendCounter + 1;
    endif
    
    if show_fer && is_cached(code, 'sfer')
      frame_error_rate = load_cached(code, 'sfer');
      
      semilogy(x_vals(snr_start<=x_vals & x_vals <= snr_end), frame_error_rate(snr_start<=x_vals & x_vals <= snr_end) , "color", color, "linestyle", '--',"marker", 'o', 'linewidth',linewidth, "markersize",3);
      legend_entries{legendCounter} = [code_name ', Soft Dec. FER'];
      legendCounter = legendCounter + 1;
    endif
    
    if is_cached(code,'sber')
      bit_error_rate = load_cached(code, 'sber');
      
      semilogy(x_vals(snr_start<=x_vals & x_vals <= snr_end), bit_error_rate(snr_start<=x_vals & x_vals <= snr_end), "color", color,'linewidth',linewidth);
      
      if show_hd || show_fer
        legend_entries{legendCounter} = [code_name ', Soft Dec. BER'];
      else
        legend_entries{legendCounter} = code_name;
      endif
      legendCounter = legendCounter + 1;
    endif
    
  endfor
  
  switch (snr_def)
    case "eb_real"
    xlabel('E_b/N_0 in dB (real-valued channel)');
    case "es_real"
    xlabel('E_s/N_0 in dB (real-valued channel)');
    case "eb_complex"
    xlabel('E_b/N_0 in dB (complex-valued channel)');
    case "es_complex"
    xlabel('E_s/N_0 in dB (complex-valued channel)');
  endswitch
  
  if y_lim < 0
    yyy = ylim();
    ylim([10^y_lim, yyy(2)]);
  endif
  
  if show_fer
    ylabel('BER/FER');
  else
    ylabel('BER');
  endif
  grid on;
  legend(legend_entries);
  print(filename, "-dpng", output_format);
end


function error_rate = load_cached(code, simulation)
  erfname = strcat('error_curves_cache/', code, '_', simulation,'.mat')
  if exist(erfname)
    load( erfname); % error_rate
  else
    error_rate = zeros([1,length(error_rate_param)]);
  endif
endfunction

function ic = is_cached(code, simulation)
  erfname = strcat('error_curves_cache/', code, '_', simulation,'.mat')
  ic = exist(erfname)
endfunction

function y=round_down_half(x)
  y=floor(x*2)/2;
endfunction

function y=round_up_half(x)
  y=ceil(x*2)/2;
endfunction

function y=Perr(EsN0dB)
  EsN0=10.^(EsN0dB/10.0);
  y=Q(sqrt(EsN0));       % Perr
endfunction

function y=dB(x)
  y=10*log10(x);
endfunction

function y=nodB(x)
  y=10.^(x./10);
endfunction

function y = Q(x)
  y = 0.5*(1-erf(x/sqrt(2)));
endfunction
