From d993021880b22874713438035cdd96430fd4c575 Mon Sep 17 00:00:00 2001
From: ismagom <gomezi@tcd.ie>
Date: Sat, 27 Dec 2014 17:09:41 -0500
Subject: [PATCH] Add PDSCH and TDEC MEX test

---
 lte/phy/lib/fec/test/turbodecoder_test_mex.c | 107 ++++++++++
 lte/phy/lib/phch/test/pdsch_test_mex.c       | 196 +++++++++++++++++++
 matlab/tests/pdsch_bler.m                    |  95 +++++++++
 matlab/tests/turbodecoder_bler.m             |  50 +++++
 4 files changed, 448 insertions(+)
 create mode 100644 lte/phy/lib/fec/test/turbodecoder_test_mex.c
 create mode 100644 lte/phy/lib/phch/test/pdsch_test_mex.c
 create mode 100644 matlab/tests/pdsch_bler.m
 create mode 100644 matlab/tests/turbodecoder_bler.m

diff --git a/lte/phy/lib/fec/test/turbodecoder_test_mex.c b/lte/phy/lib/fec/test/turbodecoder_test_mex.c
new file mode 100644
index 000000000..b34d39b6f
--- /dev/null
+++ b/lte/phy/lib/fec/test/turbodecoder_test_mex.c
@@ -0,0 +1,107 @@
+/**
+ *
+ * \section COPYRIGHT
+ *
+ * Copyright 2013-2014 The libLTE Developers. See the
+ * COPYRIGHT file at the top-level directory of this distribution.
+ *
+ * \section LICENSE
+ *
+ * This file is part of the libLTE library.
+ *
+ * libLTE is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * libLTE is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Lesser General Public License for more details.
+ *
+ * A copy of the GNU Lesser General Public License can be found in
+ * the LICENSE file in the top-level directory of this distribution
+ * and at http://www.gnu.org/licenses/.
+ *
+ */
+
+#include <string.h>
+#include "liblte/phy/phy.h"
+#include "liblte/mex/mexutils.h"
+
+/** MEX function to be called from MATLAB to test the channel estimator 
+ */
+
+#define INPUT   prhs[0]
+#define NITERS  prhs[1]
+#define NOF_INPUTS 1
+
+
+void help()
+{
+  mexErrMsgTxt
+    ("[decoded_bits] = liblte_turbodecoder(input_llr, nof_iterations)\n\n");
+}
+
+/* the gateway function */
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
+{
+
+  tdec_t tdec;
+  float *input_llr;
+  uint8_t *output_data; 
+  uint32_t nof_bits; 
+  uint32_t nof_iterations; 
+  
+  if (nrhs < NOF_INPUTS) {
+    help();
+    return;
+  }
+  
+  // Read input symbols
+  uint32_t nof_symbols = mexutils_read_f(INPUT, &input_llr);
+  if (nof_symbols < 40) {
+    mexErrMsgTxt("Minimum block size is 40\n");
+    return; 
+  }
+  nof_bits = (nof_symbols-12)/3;
+  
+  if (!lte_cb_size_isvalid(nof_bits)) {
+    mexErrMsgTxt("Invalid codeblock size\n");
+    return; 
+  }
+
+
+  // read number of iterations 
+  if (nrhs > NOF_INPUTS) {
+    nof_iterations = (uint32_t) mxGetScalar(prhs[1]);
+    if (nof_iterations > 50) {
+      mexErrMsgTxt("Maximum number of iterations is 50\n");
+      return; 
+    }
+  } else {
+    nof_iterations = 5; // set the default nof iterations to 5 as in matlab 
+  }
+  
+  // allocate memory for output bits
+  output_data = vec_malloc(nof_bits * sizeof(uint8_t));
+
+  if (tdec_init(&tdec, nof_bits)) {
+    mexErrMsgTxt("Error initiating Turbo decoder\n");
+    return;
+  }
+
+  tdec_run_all(&tdec, input_llr, output_data, nof_iterations, nof_bits);
+
+  if (nlhs >= 1) { 
+    mexutils_write_uint8(output_data, &plhs[0], nof_bits, 1);  
+  }
+
+  tdec_free(&tdec);
+
+  free(input_llr);
+  free(output_data);
+
+  return;
+}
+
diff --git a/lte/phy/lib/phch/test/pdsch_test_mex.c b/lte/phy/lib/phch/test/pdsch_test_mex.c
new file mode 100644
index 000000000..0cf8429be
--- /dev/null
+++ b/lte/phy/lib/phch/test/pdsch_test_mex.c
@@ -0,0 +1,196 @@
+/**
+ *
+ * \section COPYRIGHT
+ *
+ * Copyright 2013-2014 The libLTE Developers. See the
+ * COPYRIGHT file at the top-level directory of this distribution.
+ *
+ * \section LICENSE
+ *
+ * This file is part of the libLTE library.
+ *
+ * libLTE is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+ *
+ * libLTE is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Lesser General Public License for more details.
+ *
+ * A copy of the GNU Lesser General Public License can be found in
+ * the LICENSE file in the top-level directory of this distribution
+ * and at http://www.gnu.org/licenses/.
+ *
+ */
+
+#include <string.h>
+#include "liblte/phy/phy.h"
+#include "liblte/mex/mexutils.h"
+
+/** MEX function to be called from MATLAB to test the channel estimator 
+ */
+
+#define ENBCFG  prhs[0]
+#define RNTI    prhs[1]
+#define TBS     prhs[2]
+#define INPUT   prhs[3]
+#define NOF_INPUTS 4
+
+void help()
+{
+  mexErrMsgTxt
+    ("[decoded_ok, llr, rm, bits, symbols] = liblte_pdsch(enbConfig, RNTI, rxWaveform)\n\n");
+}
+
+/* the gateway function */
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
+{
+  int i; 
+  lte_cell_t cell; 
+  pdsch_t pdsch;
+  chest_dl_t chest; 
+  lte_fft_t fft; 
+  uint32_t cfi, sf_idx; 
+  uint16_t rnti; 
+  cf_t *input_fft, *input_signal;
+  int nof_re; 
+  ra_mcs_t mcs;
+  ra_prb_t prb_alloc;
+  pdsch_harq_t harq_process;
+  uint32_t rv;
+
+  if (nrhs != NOF_INPUTS) {
+    help();
+    return;
+  }
+    
+  if (mexutils_read_cell(ENBCFG, &cell)) {
+    help();
+    return;
+  }
+  
+  if (mexutils_read_uint32_struct(ENBCFG, "CFI", &cfi)) {
+    help();
+    return;
+  }
+  if (mexutils_read_uint32_struct(ENBCFG, "NSubframe", &sf_idx)) {
+    help();
+    return;
+  }
+
+  if (chest_dl_init(&chest, cell)) {
+    fprintf(stderr, "Error initializing equalizer\n");
+    return;
+  }
+
+  if (lte_fft_init(&fft, cell.cp, cell.nof_prb)) {
+    fprintf(stderr, "Error initializing FFT\n");
+    return;
+  }
+  
+  rnti = (uint16_t) mxGetScalar(RNTI);
+    
+  nof_re = 2 * CPNORM_NSYMB * cell.nof_prb * RE_X_RB;
+
+  mcs.tbs = mxGetScalar(TBS);
+  mcs.mod = modulation;
+  
+  if (mexutils_read_uint32_struct(ENBCFG, "NSubframe", &sf_idx)) {
+    help();
+    return;
+  }
+  
+  prb_alloc.slot[0].nof_prb = cell.nof_prb;
+  for (i=0;i<prb_alloc.slot[0].nof_prb;i++) {
+    prb_alloc.slot[0].prb_idx[i] = i;
+  }
+  memcpy(&prb_alloc.slot[1], &prb_alloc.slot[0], sizeof(ra_prb_slot_t));
+
+  ra_prb_get_re_dl(&prb_alloc, cell.nof_prb, cell.nof_ports, 2, CPNORM);
+      
+  /** Allocate input buffers */
+  if (mexutils_read_cf(INPUT, &input_signal) < 0) {
+    mexErrMsgTxt("Error reading input signal\n");
+    return; 
+  }
+  input_fft = vec_malloc(SF_LEN_RE(cell.nof_prb, cell.cp) * sizeof(cf_t));
+  
+  // Set Channel estimates to 1.0 (ignore fading) 
+  cf_t *ce[MAX_PORTS];
+  for (i=0;i<cell.nof_ports;i++) {
+    ce[i] = vec_malloc(SF_LEN_RE(cell.nof_prb, cell.cp) * sizeof(cf_t));
+  }
+  
+  lte_fft_run_sf(&fft, input_signal, input_fft);
+
+  if (nrhs > NOF_INPUTS) {
+    cf_t *cearray; 
+    nof_re = mexutils_read_cf(prhs[NOF_INPUTS], &cearray);
+    for (i=0;i<cell.nof_ports;i++) {
+      for (int j=0;j<nof_re;j++) {
+        ce[i][j] = *cearray;
+        cearray++;
+      }
+    }
+  } else {
+    chest_dl_estimate(&chest, input_fft, ce, sf_idx);    
+  }
+  float noise_power;
+  if (nrhs > NOF_INPUTS + 1) {
+    noise_power = mxGetScalar(prhs[NOF_INPUTS+1]);
+  } else {
+    noise_power = chest_dl_get_noise_estimate(&chest);
+  }
+    
+  pdcch_extract_llr(&pdcch, input_fft, ce, noise_power, sf_idx, cfi);
+  
+  uint32_t nof_locations;
+  if (rnti == SIRNTI) {
+    nof_locations = pdcch_common_locations(&pdcch, locations, MAX_CANDIDATES, cfi);
+    formats = common_formats;
+    nof_formats = nof_common_formats;
+  } else {
+    nof_locations = pdcch_ue_locations(&pdcch, locations, MAX_CANDIDATES, sf_idx, cfi, rnti); 
+    formats = ue_formats; 
+    nof_formats = nof_ue_formats;
+  }
+  uint16_t crc_rem=0;   
+  dci_msg_t dci_msg; 
+  bzero(&dci_msg, sizeof(dci_msg_t));
+  
+  for (int f=0;f<nof_formats;f++) {
+    for (i=0;i<nof_locations && crc_rem != rnti;i++) {
+      if (pdcch_decode_msg(&pdcch, &dci_msg, &locations[i], formats[f], &crc_rem)) {
+        fprintf(stderr, "Error decoding DCI msg\n");
+        return;
+      }
+    }
+  }    
+  
+  if (nlhs >= 1) { 
+    plhs[0] = mxCreateLogicalScalar(crc_rem == rnti);
+  }
+  int nof_bits = (regs_pdcch_nregs(&regs, cfi) / 9) * 72;
+  if (nlhs >= 2) {
+    mexutils_write_f(pdcch.pdcch_llr, &plhs[1], nof_bits, 1);  
+  }
+  if (nlhs >= 3) {
+    mexutils_write_cf(pdcch.pdcch_symbols[0], &plhs[2], 36*pdcch.nof_cce, 1);  
+  }
+  
+  chest_dl_free(&chest);
+  lte_fft_free(&fft);
+  pdcch_free(&pdcch);
+  regs_free(&regs);
+  
+  for (i=0;i<cell.nof_ports;i++) {
+    free(ce[i]);
+  }
+  free(input_signal);
+  free(input_fft);
+  
+  return;
+}
+
diff --git a/matlab/tests/pdsch_bler.m b/matlab/tests/pdsch_bler.m
new file mode 100644
index 000000000..905c6f441
--- /dev/null
+++ b/matlab/tests/pdsch_bler.m
@@ -0,0 +1,95 @@
+%% PDSCH decoding based on RMC channels
+
+%% Cell-Wide Settings
+% A structure |enbConfig| is used to configure the eNodeB.
+clear
+
+Npackets = 25;
+SNR_values = linspace(1,6,4);
+
+%% Choose RMC 
+[waveform,rgrid,rmccFgOut] = lteRMCDLTool('R.12',[1;0;0;1]);
+waveform = sum(waveform,2);
+
+Nsf = 8; 
+
+%% Setup Fading channel model 
+cfg.Seed = 8;                  % Random channel seed
+cfg.NRxAnts = 1;               % 1 receive antenna
+cfg.DelayProfile = 'EPA';      % EVA delay spread
+cfg.DopplerFreq = 5;           % 120Hz Doppler frequency
+cfg.MIMOCorrelation = 'Low';   % Low (no) MIMO correlation
+cfg.InitTime = 0;              % Initialize at time zero
+cfg.NTerms = 16;               % Oscillators used in fading model
+cfg.ModelType = 'GMEDS';       % Rayleigh fading model type
+cfg.InitPhase = 'Random';      % Random initial phases
+cfg.NormalizePathGains = 'On'; % Normalize delay profile power 
+cfg.NormalizeTxAnts = 'On';    % Normalize for transmit antennas
+cfg.SamplingRate = rmccFgOut.SamplingRate; 
+
+% Setup channel equalizer
+cec.PilotAverage = 'UserDefined';     % Type of pilot averaging
+cec.FreqWindow = 9;                   % Frequency window size
+cec.TimeWindow = 9;                   % Time window size
+cec.InterpType = 'linear';             % 2D interpolation type
+cec.InterpWindow = 'Centered';        % Interpolation window type
+cec.InterpWinSize = 1;                % Interpolation window size
+
+addpath('../../debug/lte/phy/lib/phch/test')
+
+decoded = zeros(size(SNR_values));
+decoded_liblte = zeros(size(SNR_values));
+
+for snr_idx=1:length(SNR_values)
+    SNRdB = SNR_values(snr_idx);
+    SNR = 10^(SNRdB/10);    % Linear SNR  
+    N0  = 1/(sqrt(2.0*rmccFgOut.CellRefP*double(rmccFgOut.Nfft))*SNR);
+    for i=1:Npackets
+
+        %% Fading
+        rxWaveform = lteFadingChannel(cfg,waveform);
+        
+        %% Noise Addition
+        noise = N0*complex(randn(size(rxWaveform)), randn(size(rxWaveform)));  % Generate noise
+        rxWaveform = rxWaveform + noise; 
+       
+        %% Demodulate 
+        frame_rx = lteOFDMDemodulate(rmccFgOut, rxWaveform);
+
+        for sf_idx=0:Nsf
+            
+            subframe_rx=frame_rx(:,sf_idx*14+1:(sf_idx+1)*14);
+            rmccFgOut.NSubframe=sf_idx;
+            rmccFgOut.TotSubframes=1;
+
+            % Perform channel estimation
+            [hest, nest] = lteDLChannelEstimate(rmccFgOut, cec, subframe_rx);
+
+            [cws,symbols] = ltePDSCHDecode(rmccFgOut,rmccFgOut.PDSCH,subframe_rx,hest,nest);
+            [trblkout,blkcrc] = lteDLSCHDecode(rmccFgOut,rmccFgOut.PDSCH, ... 
+                                                    rmccFgOut.PDSCH.TrBlkSizes(sf_idx+1),cws);
+
+            decoded(snr_idx) = decoded(snr_idx) + ~blkcrc;
+
+
+            %% Same with libLTE
+            %[found_liblte, llr, pdcchSymbols2] = liblte_pdsch(rmccFgOut, ueConfig.RNTI, rxWaveform);
+            %decoded_liblte(snr_idx) = decoded_liblte(snr_idx)+found_liblte;
+        end
+    end
+    fprintf('SNR: %.1f\n',SNRdB)
+end
+
+if (length(SNR_values)>1)
+    semilogy(SNR_values,1-decoded/Npackets/(Nsf+1),'bo-',...
+             SNR_values,1-decoded_liblte/Npackets/(Nsf+1), 'ro-')
+    grid on;
+    legend('Matlab','libLTE')
+    xlabel('SNR (dB)')
+    ylabel('BLER')
+    axis([min(SNR_values) max(SNR_values) 1/Npackets/(Nsf+1) 1])
+else
+    disp(decoded)
+    disp(decoded_liblte)
+end
+
diff --git a/matlab/tests/turbodecoder_bler.m b/matlab/tests/turbodecoder_bler.m
new file mode 100644
index 000000000..c77ef8859
--- /dev/null
+++ b/matlab/tests/turbodecoder_bler.m
@@ -0,0 +1,50 @@
+
+clear
+blen=1008;
+SNR_values_db=linspace(-1,0.5,6);
+Nrealizations=10000;
+
+addpath('../../debug/lte/phy/lib/fec/test')
+
+errors1=zeros(1,length(SNR_values_db));
+errors2=zeros(1,length(SNR_values_db));
+for snr_idx=1:length(SNR_values_db)
+    SNRdB = SNR_values_db(snr_idx);             % Desired SNR in dB
+    SNR = 10^(SNRdB/20);    % Linear SNR  
+
+    for i=1:Nrealizations
+        Data = randi(2,blen,1)==1;
+        codedData = lteTurboEncode(Data);
+
+        codedsymbols = 2*double(codedData)-1;
+        
+        %% Additive Noise
+        N0 = 1/SNR;
+
+        % Create additive white Gaussian noise
+        noise = N0*randn(size(codedsymbols));   
+        
+        noisysymbols = noise + codedsymbols; 
+        
+        decodedData = lteTurboDecode(noisysymbols);
+        interleavedSymbols = reshape(reshape(noisysymbols,[],3)',1,[]);
+        [decodedData2] = liblte_turbodecoder(interleavedSymbols);
+        
+        errors1(snr_idx) = errors1(snr_idx) + any(decodedData ~= Data);
+        errors2(snr_idx) = errors2(snr_idx) + any(decodedData2 ~= Data);
+    end
+    
+    fprintf('SNR: %.2f\n', SNR_values_db(snr_idx));
+end
+
+if (length(SNR_values_db) > 1)
+    semilogy(SNR_values_db, errors1/Nrealizations, ...
+             SNR_values_db, errors2/Nrealizations)
+    grid on
+    xlabel('SNR (dB)')
+    ylabel('BLER')
+    legend('Matlab','libLTE');
+else
+    disp(errors1);
+    disp(errors2);
+end
\ No newline at end of file