% ltehdlTransmitter model Post Simulation script - ltehdlTransmitter_PostSim
% Uses workspace variables generated by the ltehdlTransmitter model

if exist('tx_samples', 'var') && length(squeeze(tx_samples.Data))>1
    % filtered output from "HDL LTE MIMO Transmitter" subsystem
    valid = squeeze(tx_valid.Data);
    samples = squeeze(tx_samples.Data);
    
    txOut{1} = samples(valid(:,1),1);
    txOut{2} = samples(valid(:,2),2);
    txOut{3} = samples(valid(:,3),3);
    txOut{4} = samples(valid(:,4),4);
    
    NoofSamples = 1920*TotalSubframes;
    
    if length(txOut{1}) >= NoofSamples
        txSamples = [];
        for portNumber = 1:Nant
            txSamples = [txSamples double(txOut{portNumber})];
        end
        
        simParams = enb;
        
        if (strcmpi(simParams.CyclicPrefix,'Normal'))
            out.noSym = 14;
        else
            out.noSym = 12;
        end
        
        samp = [1920 3840 7680 15360 30720 30720];
        % No of Samples per subframe
        out.NSampperSF = samp([6 15 25 50 75 100] == simParams.NDLRB);
        
        %% Fading channel model
        
        % Setting Delay Profile to EPA will create fading channel model as defined in
        % Annex B.2.1 of 3GPP TS 36.101 specification
        % For more details on fading channel model, execute:
        % web(fullfile(docroot, 'lte/ref/ltefadingchannel.html?container=jshelpbrowser'))
        % This channel configurations works only for default enb with a change
        % only in PHICHDuration and Ng fields.
        
        % Setting Delay Profile to off will create a fading-free channel
        
        chcfg.DelayProfile = 'EPA';               % {'off','EPA'}
        
        % The below channel configuration works for the changes only in
        % enb.PHICHDuration and enb.Ng fields. Changes in any other fields will set
        % the fading channel to fading free channel.
        if (simParams.NDLRB ~= 6) || (simParams.CellRefP ~= 4) ...
                || ~strcmp(simParams.CyclicPrefix,'Normal') ...
                || ~isequal(tx_cellids,[390 89 501 231 500])
            chcfg.DelayProfile = 'off';
        end
        
        % DO NOT MODIFY below chcfg configuration.
        if ~strcmpi(chcfg.DelayProfile,'off')
            chcfg.DopplerFreq = 5;
            chcfg.SamplingRate = 30.72e6;
            chcfg.InitTime = 0;
            chcfg.NTerms = 16;
            chcfg.ModelType = 'GMEDS';
            chcfg.NormalizePathGains = 'On';
            chcfg.InitPhase = 'Random';
            chcfg.Seed = 1;
        end
        
        chcfg.NRxAnts = 1;
        chcfg.MIMOCorrelation = 'Medium';
        chcfg.NormalizeTxAnts = 'On';
        
        rxSamples = lteFadingChannel(chcfg,txSamples);
        
        %% Cell Search using lteCellSearch
        cell_text = [];
        out.counter = 1;
        cell_text = [cell_text ...
            sprintf(['NCellID after Cell search (Number of Tx Antennas = ' num2str(Nant) '):\n'])];
        
        expCellIDs = [];
        detecCellIDs = [];
        for startIdx = 1:out.NSampperSF*10:length(rxSamples) - out.NSampperSF
            out.rx_data = rxSamples(startIdx:startIdx + out.NSampperSF - 1);
            [out.NCellID,offset(out.counter)] = lteCellSearch(simParams,double(out.rx_data)); %#ok<*SAGROW>
            
            cell_text = [cell_text sprintf('..................................\n')];
            cell_text = [cell_text sprintf('RadioFrame# %d \n',out.counter)]; %#ok<*AGROW>
            expCellID = tx_cellids(min(out.counter,length(tx_cellids)));
            detecCellID = out.NCellID;
            cell_text = [cell_text sprintf('Detected: %d Expected: %d \n',...
                detecCellID, expCellID)];
            expCellIDs = [expCellIDs expCellID];
            detecCellIDs = [detecCellIDs detecCellID];
            
            out.counter = out.counter +1;
        end
        cell_text = [cell_text newline];
        
        figCell = msgbox(cell_text,'Cell Search results');
        
        %% Compare Filtered output to LTE toolbox functions
        % Populate grid with PSS,SSS,Cell reference signals and PBCH data.
        
        sfidx = 1;
        simParams.NCellID = double(tx_cellids(sfidx));
        Frame = [];
        
        subframe = lteDLResourceGrid(simParams);
        for idx = 1:TotalSubframes
            out.cellrs = lteCellRS(simParams);
            out.cellrs = reshape(out.cellrs,simParams.NDLRB*4,length(out.cellrs)/(simParams.NDLRB*4));
            out.cellrsIdx = lteCellRSIndices(simParams);
            out.cellrsIdx = reshape(out.cellrsIdx,simParams.NDLRB*4,length(out.cellrsIdx)/(simParams.NDLRB*4));
            
            % Only PBCH Tx so other frames will be zero already.
            subframe(out.cellrsIdx) = out.cellrs;
            
            % Synchronization signals formation
            if (simParams.NSubframe == 0 || simParams.NSubframe == 5)
                out.pss = ltePSS(simParams);
                out.sss = lteSSS(simParams);
                out.pssIdx = ltePSSIndices(simParams);
                out.sssIdx = lteSSSIndices(simParams);
                
                if simParams.CellRefP == 1
                    subframe(out.pssIdx) = out.pss;
                    subframe(out.sssIdx) = out.sss;
                else
                    for i = 1:simParams.CellRefP
                        subframeSG = subframe(:,:,i);
                        subframeSG(out.pssIdx) = out.pss;
                        subframeSG(out.sssIdx) = out.sss;
                        subframe(:,:,i) = subframeSG;
                    end
                end
            end
            
            % PBCH data formation
            if (simParams.NSubframe == 0)
                if (simParams.NSubframe == 0 && mod(simParams.NFrame,4) == 0)
                    out.mib = lteMIB(simParams);
                    out.bch = lteBCH(simParams,out.mib);
                    out.pbch = ltePBCH(simParams,out.bch);
                    out.pbchIdx = ltePBCHIndices(simParams);
                else
                    out.pbchIdx = ltePBCHIndices(simParams);
                end
                pbchVal_idx = (mod(simParams.NFrame,4)*countVal+1:(mod(simParams.NFrame,4) + 1)*countVal).';
                subframe(out.pbchIdx) = out.pbch(pbchVal_idx,:);
            end
            
            Frame = [Frame subframe]; %#ok<*AGROW>
            
            if (simParams.NSubframe == 9)
                simParams.NCellID = double(tx_cellids(min(sfidx + 1,length(tx_cellids))));
                simParams.NFrame = simParams.NFrame + 1;
                simParams.NSubframe = 0;
                sfidx = sfidx + 1;
            else
                simParams.NSubframe = simParams.NSubframe + 1;
            end
            subframe (:,:) = 0;
        end
        
        out.ofdm_data = lteOFDMModulate(simParams,Frame,0);
        
        if (simParams.NDLRB == 6)
            simParams.lteBW = 1.4;
        else
            simParams.lteBW = simParams.NDLRB/5;
        end
        
        BW = simParams.lteBW;
        
        reset(hF);
        ref_data = hF(out.ofdm_data);
        
        % Show outputs of filtered samples
        TxInfo = lteOFDMInfo(simParams);
        out.Npsd = 64*1024;                    % Number of points for PSD FFT
        out.fs = TxInfo.SamplingRate/1e6;
        out.ff = (-out.Npsd/2:1:out.Npsd/2-1)/out.Npsd*out.fs; % Frequency vector for x-axis of plot
        
        figTx = figure('Visible', 'off');
        screenSize = get(groot, 'ScreenSize');
        figTx.Position(3) = screenSize(3) * (0.3+0.1*Nant);
        figTx.Position(4) = screenSize(4) * 0.5;
        movegui(figTx, 'center');
        
        for portNumber = 1:Nant
            % Plot transmitted signal PSD and time-domain waveform for '1' subframe
            txSamples = double(txOut{portNumber});
            
            
            subplot(2,Nant,portNumber);
            plot(out.ff,10*log10(abs(fftshift(pwelch(ref_data(:,portNumber),[],[],out.Npsd)))),'b');
            hold on;
            plot(out.ff,10*log10(abs(fftshift(pwelch(txSamples(1:length(ref_data(:,portNumber))) ...
                ,[],[],out.Npsd)))),'--r');
            title(sprintf('PSD of OFDM modulated LTE signal\n(BW = %s MHz)\n(port number : %d/%d)',string(BW),portNumber,Nant));
            xlabel('Frequency (MHz)');
            ylabel('PSD (dB)');
            legend('LTE Toolbox reference','HDL Implementation','Location', 'south');
            
            % Plot absolute value of transmitted signal and reference transmitted
            % signal
            subplot(2,Nant,Nant+portNumber);
            plot(1:length(ref_data(:,portNumber)),abs(ref_data(1:length(ref_data),portNumber)),'b');
            hold on;
            plot(1:length(ref_data(:,portNumber)),abs(txSamples(1:length(ref_data))),'--r');
            hold on;
            err(:,portNumber) = abs(abs(ref_data(1:length(ref_data),portNumber)) - abs(txSamples(1:length(ref_data))));
            plot(1:length(ref_data(:,portNumber)),err(:,portNumber),'k');
            axis([1 length(ref_data(:,portNumber)) -0.002 max(abs(ref_data(1:length(ref_data),portNumber)))]);
            xlabel('No of Samples');
            ylabel('Magnitude');
            legend('LTE Toolbox reference','HDL Implementation','Difference');
        end
        figTx.Visible = 'on';
        
        %% MIB Decoding
        
        mib_text = sprintf(['Cell-wide settings after MIB decoding(Number of Tx antennas = ' num2str(Nant) '):\n']);
        
        out.rx_data = rxSamples;
        enbRx.NDLRB = 6;
        enbRx.CyclicPrefix = enb.CyclicPrefix;
        RxInfo = lteOFDMInfo(enbRx);
        
        out.resample_data = resample(double(out.rx_data(offset(1)+1:end)), ...
            RxInfo.SamplingRate,round(TxInfo.SamplingRate));
        
        out.rxgrid = lteOFDMDemodulate(enbRx,double(out.resample_data));
        enbRx.NCellID = double(tx_cellids(1));
        
        out.pbchIndices = ltePBCHIndices(enbRx);
        out.pbchRx = lteExtractResources(out.pbchIndices, out.rxgrid);
        [~, ~, out.nfmod4, out.mib,~] = ltePBCHDecode(enbRx, out.pbchRx);
        enbDet = lteMIB(out.mib, enbRx);
        enbDet.NFrame = out.nfmod4 + enbDet.NFrame;
        
        mib_text = [mib_text sprintf('..................................\n')];
        mib_text = [mib_text sprintf('NDLRB: %d \nPHICH Duration: %s \nNg: %s \nSFN: %d \n',...
            enbDet.NDLRB, enbDet.PHICHDuration, enbDet.Ng, enbDet.NFrame)];
        mib_text = [mib_text newline];
        
        
        figMIB = msgbox(mib_text,'MIB Decoding Results');
    end
    
end

%% Clean up workspace variables
clear txSamples tx_valid startIdx samples valid BW i pbchVal_idx...
    mib_text cell_text  tx_samples portNumber subframeSG...
    txOut enbRx out idx TxInfo subframe offset ref_data expCellID detecCellID ...
    sfidx RxInfo Frame samp fig screenSize rxSamples chcfg simParams NoofSamples