%% Reporting of Rank Indicator (RI) Conformance Test
% This example demonstrates how to measure the Rank Indicator (RI)
% reporting performance using the LTE Toolbox(TM) under conformance test
% conditions as defined in TS36.101 Section 9.5.1.1 [ <#9 1> ].

% Copyright 2012-2018 The MathWorks, Inc.

%% Introduction
% This example highlights the use of the 
% <matlab:doc('lteRISelect') lteRISelect> function which provides 
% estimation of the RI. The performance of the RI estimation is also
% tested. The performance requirement defined in TS36.101 Section 9.5.1.1 
% [ <#9 1> ] is as follows:
%
% * the ratio of the throughput obtained when transmitting based on UE
% reported RI and that obtained when transmitting with fixed rank 2 shall
% be >= 1 (Test 1)
%
% * the ratio of the throughput obtained when transmitting based on UE
% reported RI and that obtained when transmitting with fixed rank 1 shall
% be >= 1.05 (Test 2)
%
% This example tests that these requirements are met. 

%% Simulation Configuration
% The example is executed for a simulation length of 10 frames at an SNR of
% 0.0dB. The fixed RI value is configured to be 2, therefore simulating
% TS36.101 Section 9.5.1.1 Test 1 [ <#9 1> ]. A large number of |NFrames|
% should be used to produce meaningful results. The variable |FixedRI|
% controls which of the performance requirements described in the
% introduction are tested: |FixedRI=2| corresponds to Test 1 and
% |FixedRI=1| corresponds to Test 2. Note that a target SNR of 20.0dB
% rather than 0.0dB applies for Test 2. 

NFrames = 10;
SNRdB = 0.0;
FixedRI = 2;   

%% Propagation Channel Model Configuration
% The structure, |channel|, contains the channel model configuration
% parameters. Note that TS36.101 Section 9.5.1.1 Test 3 [ <#9 1> ] can also
% be implemented, by setting |channel.MIMOCorrelation = 'High'| (a target
% SNR of 20.0dB applies to Test 3, and the throughput ratio target is
% different).

channel.Seed = 1;                   % Channel seed
channel.NRxAnts = 2;                % 2 receive antennas
channel.DelayProfile = 'EPA';       % Delay profile
channel.DopplerFreq = 5.0;          % Doppler frequency    
channel.MIMOCorrelation = 'Low';    % Multi-antenna correlation
channel.ModelType = 'GMEDS';        % Rayleigh fading model type
channel.NormalizeTxAnts = 'On';     % Normalize for transmit antennas
channel.NormalizePathGains = 'On';  % Normalize delay profile power 
channel.InitPhase = 'Random';       % Random initial phases
channel.NTerms = 16;                % Oscillators used in fading model
    
%% Channel Estimator Configuration
% The channel estimator is configured with a structure |cec|. The variable
% |perfectChanEstimator| controls channel estimator behavior. Valid values
% are |true| or |false|. When set to |true| a perfect channel estimator is
% used otherwise an imperfect estimate of the channel is used, based on the
% values of received pilot signals. In this example, we enable the perfect
% channel estimator. 

% Configure channel estimator
cec.PilotAverage = 'UserDefined';   % Type of pilot symbol averaging
cec.FreqWindow = 9;                 % Frequency window size in REs
cec.TimeWindow = 9;                 % Time window size in REs
cec.InterpType = 'Cubic';           % 2D interpolation type
cec.InterpWindow = 'Centered';      % Interpolation window type
cec.InterpWinSize = 1;              % Interpolation window size  
        
% Channel estimator behavior
perfectChanEstimator = true;

%% Set RI and CQI/PMI Delays
% Initialize RI and CQI/PMI reporting delays in subframes; the RI delay is
% 9 subframes rather than 8 subframes in accordance with Note 5 of TS36.101
% Table 9.5.1.1-1 [ <#9 1> ]. Note that the feedback of the RI, CQI and
% PMI is assumed to be perfect, with the values being fed back in buffers
% rather than being fed back in an uplink transmissions.
    
riDelay = 9; % subframes
cqipmiDelay = 8; % subframes
        
%% Codebook Subset Restriction Bitmaps
% In order to control the PMI selection which underlies the rank selection,
% codebook subset restriction bitmaps are used as described in TS36.213
% Section 7.2 [ <#9 2> ]. |codebookSubsetRI1| configures the codebook
% subset restriction to the two precoders for rank 1, used for |FixedRI=1|.
% |codebookSubsetRI2| configures the codebook subset restriction to the
% single precoder of rank 2, used for |FixedRI=2|.
% |codebookSubsetUEReported| is the union of these codebook subset
% restrictions and allows the rank to be selected dynamically in the case
% of UE reported RI.

codebookSubsetRI1 = '000011';
codebookSubsetRI2 = '010000';
codebookSubsetUEReported = '010011';

%% System Processing 
% The main processing is split into two phases, configured via the
% |riConfig| loop variable. These phases implement the two measurements
% required in the performance test defined in TS36.101 Section 9.5.1.1
% [ <#9 1> ]:
%
% * _Fixed RI_. The first phase (|riConfig=1|) performs PDSCH transmission
% and reception where the rank is set to a fixed value (either 1 or 2,
% selected by the variable |FixedRI| above). The final throughput is
% recorded.
%
% * _UE reported RI_. In the second phase (|riConfig=2|), the number of
% transmission layers is selected on the basis of UE reported RI, with the
% RI being updated every 5 subframes. The final throughput is recorded.
% The throughput ratio (|measuredGamma|) between the fixed RI phase and the
% UE reported phase is calculated and checked against the specified
% performance requirement.
%
% Note that CQI and PMI are also updated every 5 subframes during both
% simulation phases.
%
% The processing is performed on a subframe by subframe basis using the
% following steps:
%
% * _Select RI._ The current RI is read from the oldest value in the RI
% buffer |riBuffer|. The current RI is used to configure the number of
% transmission layers and codewords for the PDSCH.
%
% * _Select CQI/PMI._ The current CQI and PMI are read from the oldest
% values in the CQI and PMI buffers |cqiBuffer| and |pmiBuffer|. The
% current PMI is used to configure the precoding for the PDSCH. Note that
% for RI=2, there are two CQI values, one for each codeword; for RI=1 there
% is one CQI value for the single codeword. 
%
% * _Select MCS according to CQI._ The Modulation and Coding Scheme (MCS)
% index corresponding to the CQI is selected for each codeword by means of
% a lookup table defined in TS36.101 Table A.4-1 [ <#9 1> ] CSI RMC RC.2
% FDD (MCS.2).
%
% * _Determine Transport Block Size and modulation order._ The MCS index
% for each codeword is passed to the <matlab:doc('lteMCS') lteMCS> function
% which calculates the corresponding Transport Block Size (TBS) index and
% modulation order; the <matlab:doc('lteTBS') lteTBS> function is then used
% to calculate the TBS for each codeword from the TBS index and the number
% of resource blocks allocated to the PDSCH.
%
% * _Transmit and receive waveform._ Transport block data is generated for
% one or two codewords as appropriate and passed to
% <matlab:doc('lteRMCDLTool') lteRMCDLTool> to create a transmitted
% downlink waveform. This waveform is then passed through a fading channel
% and AWGN noise is added. The received signal is synchronized and OFDM
% demodulated and channel estimation is performed.
%
% * _Measure PDSCH throughput._ The PDSCH and DL-SCH are decoded and the
% CRC pass/fail for each codeword is recorded to determine the data
% throughput.
%
% * _Update RI._ If an RI update is scheduled in this subframe, use the
% channel estimate to update the RI with the 
% <matlab:doc('lteRISelect') lteRISelect> function. The updated RI value is 
% recorded in the RI buffer. The codebook subset restriction bitmap ensures
% that the reported RI will remain fixed for the fixed RI phase but will be
% dynamically selected for the UE reported RI phase. If an RI update is not
% scheduled in this subframe, the previous RI value is reused.
%
% * _Update CQI/PMI._ If a CQI/PMI update is scheduled in this subframe,
% use the channel estimate to update the CQI and PMI with the
% <matlab:doc('lteCQISelect') lteCQISelect> and 
% <matlab:doc('ltePMISelect') ltePMISelect> functions. The updated CQI and
% PMI values are recorded in the CQI and PMI buffers. If a CQI/PMI update
% is not scheduled in this subframe, the previous CQI and PMI values are
% reused.

% Check that fixed RI valid is either 1 or 2 
if (~any(FixedRI==[1 2]))
    error('Fixed RI value specified must be 1 or 2.');
end     

% For each RI configuration (fixed RI and UE reported):
for riConfig = 1:2

    if (riConfig==1)
        riConfigStr = sprintf('fixed RI=%d',FixedRI);
    else    
        riConfigStr = 'UE reported RI';
    end    
    fprintf('\nSimulating with %s at %gdB SNR for %d Frame(s)\n', ...
        riConfigStr,SNRdB,NFrames);
    
    % Set up eNodeB settings for 2 codewords
    enb = struct('RC','R.3');              
    enb.CellRefP = 2;       
    enb.CFI = 3;
    enb.OCNGPDSCHEnable = 'On';
    enb.TotSubframes = 1;    
    enb.PDSCH.RVSeq = 0;
    enb.PDSCH.CSIMode = 'PUCCH 1-1';    
    enb.PDSCH.TxScheme = 'SpatialMux'; 
    enb.PDSCH.NLayers = 2;        
    enb.PDSCH.Rho = -3.0; 
    enb.PDSCH.CSI = 'On';
    ncw = 2;
    enb = lteRMCDL(enb,ncw); 
    
    % Set channel model sampling rate
    ofdmInfo = lteOFDMInfo(enb);
    channel.SamplingRate = ofdmInfo.SamplingRate;   

    % Configure appropriate codebook subset restriction and target
    % throughput ratio
    if (riConfig==1)
        if (FixedRI==1)
            enb.PDSCH.CodebookSubset = codebookSubsetRI1;
            targetGamma = 1.05;
        else
            enb.PDSCH.CodebookSubset = codebookSubsetRI2;
            targetGamma = 1.0;
        end
    else            
        enb.PDSCH.CodebookSubset = codebookSubsetUEReported;
    end   

    % Initialize variables
    totalCRC = [];      % CRC values, used for throughput calculation
    totalTBS = [];      % TBS values, used for throughput calculation
    RIReport = [];      % reported RI values    
    riXaxis = [];       % corresponding subframe numbers
    PMIReport = [];     % reported PMI values
    CQIReport = [];     % reported CQI values
    SINRReport = [];    % corresponding SINR values
    cqipmiXaxis = [];   % corresponding subframe numbers
    offsets = 0;        % Initialize frame offset value
    rng('default');     % Default random number generator seed         
        
    % Initialize RI/PMI/CQI buffers; note that the PDSCH throughput
    % resulting from these initial values is ignored, as the throughput
    % recording waits until the buffers have been filled with meaningful
    % reports    
    riBuffer = ones(1,riDelay);    
    pmiBuffer = ones(1,cqipmiDelay);
    cqiBuffer = ones(ncw,cqipmiDelay);

    % For each subframe:
    for subframeNo = 0:(NFrames*10-1)

        % Update subframe number
        enb.NSubframe = mod(subframeNo,10);
        
        % Update PMI, reading the oldest value from the PMI buffer
        cqipmiPtr = mod(subframeNo,cqipmiDelay);    
        enb.PDSCH.PMISet = pmiBuffer(cqipmiPtr+1);
        
        % Update number of transmission layers based on RI, reading the
        % oldest value from the RI buffer
        riPtr = mod(subframeNo,riDelay);
        enb.PDSCH.NLayers = riBuffer(riPtr+1);                        
        
        % Update MCS according to CQI, reading the oldest value from the
        % CQI buffer
        CQI = cqiBuffer(:,cqipmiPtr+1);
        % Select MCS according to CQI using TS36.101 Table A.4-1 CSI RMC
        % RC.2 FDD (MCS.2), which defines the relationship between CQI
        % indices and MCS indices
        IMCSTable = [-1 0 0 2 4 6 8 11 13 15 18 20 22 24 26 27];        
        IMCS = IMCSTable(CQI+1);
        
        % Determine TBS and modulation order
        [ITBS,modulation] = lteMCS(IMCS);
        if (mod(enb.NSubframe,5)==0)
            TBS = [0; 0];
        else
            TBS = double(lteTBS(size(enb.PDSCH.PRBSet,1),ITBS));
        end
        enb.PDSCH.TrBlkSizes(:,enb.NSubframe+1) = TBS;
        
        % Determine if an RI or PMI/CQI update is required in this
        % subframe, according to reporting periodicity N_pd = 5ms and
        % configuration indices cqi-pmi-ConfigurationIndex = 6 and
        % ri-ConfigurationInd = 1 from TS36.101 Table 9.5.1.1-1
        riPeriod = 5; % periodicity N_pd
        riOffset = 3; % offset deriving from configuration indices
        riUpdate = (mod(subframeNo,riPeriod)==riOffset);        
        cqipmiPeriod = 5; % periodicity N_pd
        cqipmiOffset = 4; % offset deriving from cqi-pmi-ConfigurationIndex        
        cqipmiUpdate = (mod(subframeNo,cqipmiPeriod)==cqipmiOffset);
        
        % Establish if this subframe actually needs executed for PDSCH
        % reception, RI/PMI/CQI estimation or initial timing offset
        % estimation:
        firstRxSubframe = max(cqipmiDelay+cqipmiOffset,riDelay+riOffset);
        if((any(TBS) && subframeNo>=firstRxSubframe) || ...
                (riUpdate || cqipmiUpdate) || subframeNo==0)
                                    
            % Configure data for transmission 
            % (1 codeword on 1 layer, 2 codewords on 2 layers)
            data = {randi([0 1],TBS(1),1)};
            if (enb.PDSCH.NLayers==2)                
                data = {data{1} randi([0 1],TBS(2),1)};                
            end            
            enb.PDSCH.Modulation = modulation(1:numel(data));

            % Create OFDM resource grid containing RMC transmission and
            % perform OFDM modulation
            txWaveform = lteRMCDLTool(enb,data);            

            % The initialization time for channel modeling is set each
            % subframe to simulate a continuously varying channel
            channel.InitTime = subframeNo/1000;            

            % Pass data through the fading channel model. 
            % An additional 25 samples are added to the end of the
            % waveform. These are to cover the range of delays expected
            % from the channel modeling (a combination of implementation
            % delay and channel delay spread)
            rxWaveform = lteFadingChannel(channel, ...
                            [txWaveform ; zeros(25,size(txWaveform,2))]);  

            % Calculate noise gain including compensation for downlink
            % power allocation
            SNR = 10^((SNRdB-enb.PDSCH.Rho)/20);

            % Normalize noise power to take account of sampling rate, which
            % is a function of the IFFT size used in OFDM modulation, and
            % the number of antennas
            N0 = 1/(sqrt(2.0*enb.CellRefP*double(ofdmInfo.Nfft))*SNR);

            % Create additive white Gaussian noise
            noise = N0*complex(randn(size(rxWaveform)), ...
                                randn(size(rxWaveform)));

            % Add AWGN to the received time domain waveform    
            rxWaveform = rxWaveform + noise;

            % Perform synchronization
            % An offset within the range of delays expected from the
            % channel modeling (a combination of implementation delay and
            % channel delay spread) indicates success
            if (mod(subframeNo,10)==0)
                offset = lteDLFrameOffset(enb,rxWaveform);  
                if (offset > 25)
                    offset = offsets(end);
                end
                offsets = [offsets offset]; %#ok<AGROW>
            end
            rxWaveform = rxWaveform(1+offset:end,:);

            % Perform OFDM demodulation on the received data to create the
            % received resource grid
            rxGrid = lteOFDMDemodulate(enb,rxWaveform);

            % Channel estimation
            if (perfectChanEstimator)
                chEstGrid = ...
                   lteDLPerfectChannelEstimate(enb,channel,offset);
                n = lteOFDMDemodulate(enb,noise(1+offset:end,:));
                noiseEst = var(reshape(n,numel(n),1));
            else
                [chEstGrid,noiseEst] = ... 
                    lteDLChannelEstimate(enb,enb.PDSCH, ...
                        cec,rxGrid); %#ok<UNRCH>
            end

            % If this subframe requires PDSCH reception:
            if(any(TBS) && subframeNo>=firstRxSubframe)                                        

                % Decode the PDSCH 
                ind = ltePDSCHIndices(enb,enb.PDSCH,enb.PDSCH.PRBSet);                
                pdschRx = lteExtractResources(ind,rxGrid) * ...
                            (10^(-enb.PDSCH.Rho/20));
                pdschChEst = lteExtractResources(ind,chEstGrid);
                [rxBits,rxSymbols] = ltePDSCHDecode(enb,enb.PDSCH, ...
                                        pdschRx,pdschChEst,noiseEst);

                % Decode the DL-SCH
                [decbits,crc] = lteDLSCHDecode(enb,enb.PDSCH, ...
                    enb.PDSCH.TrBlkSizes(1:numel(rxBits),enb.NSubframe+1),rxBits);

                % Record the CRC and TBS values for final throughput
                % calculation                
                totalCRC = [totalCRC crc];                  %#ok<AGROW>
                totalTBS = [totalTBS TBS(1:numel(crc)).'];  %#ok<AGROW>                

            end

            % Update RI:
            if (riUpdate)
                                
                % Update RI
                thisRI = lteRISelect(enb,enb.PDSCH,chEstGrid,noiseEst);                       
                                
                % Feed the value back to UE (in a buffer)
                riBuffer(riPtr+1) = thisRI;
                
                % Record values for plotting
                if (riConfig==2)
                    RIReport = [RIReport thisRI];           %#ok<AGROW>                    
                    riXaxis = [riXaxis subframeNo];         %#ok<AGROW>
                end

            end
                           
            % Update CQI and PMI:
            if (cqipmiUpdate)
                                                
                % Update PMI conditioned on the most recently reported RI,
                % or the fixed RI if RI reporting is not configured
                if (isempty(RIReport))
                    enb.PDSCH.NLayers = FixedRI;
                else
                    enb.PDSCH.NLayers = RIReport(end);
                end
                thisPMI = ltePMISelect(enb,enb.PDSCH,chEstGrid,noiseEst);                       
                
                % Update CQI conditioned on the most recently reported PMI
                % and RI; this includes configuring the number of codewords
                % for CQI selection based on the RI
                enb.PDSCH.PMISet = thisPMI;
                enb.PDSCH.NCodewords = enb.PDSCH.NLayers;
                [thisCQI,thisSINR] = ...
                    lteCQISelect(enb,enb.PDSCH,chEstGrid,noiseEst);   

                % Feed the values back to UE (in buffers)               
                pmiBuffer(:,cqipmiPtr+1) = thisPMI;                       
                % For CSIMode='PUCCH 1-1', the CQI for the second codeword
                % is reported as a differential from the first codeword;
                % here we convert to absolute CQI indices for both
                % codewords.
                if (numel(thisCQI)==2)
                    thisCQI(2) = thisCQI(1) - thisCQI(2); 
                end            
                % Transmit with lowest CQI=1 when CQI selection reports
                % CQI=0 (out of range)                              
                thisCQI(thisCQI==0) = 1;                
                
                cqiBuffer(:,cqipmiPtr+1) = thisCQI;              

                % Record values for plotting
                if (riConfig==2)                    
                    PMIReport = [PMIReport thisPMI];        %#ok<AGROW>
                    CQIReport = [CQIReport thisCQI(1)];     %#ok<AGROW>
                    SINRReport = [SINRReport thisSINR(1)];  %#ok<AGROW>
                    cqipmiXaxis = [cqipmiXaxis subframeNo]; %#ok<AGROW>
                end
            
            end

        end
                
        % For subframes where RI was not updated, re-use the previous value
        % in the buffer
        if (~riUpdate)
            riBuffer(riPtr+1) = riBuffer(mod(riPtr-1,riDelay)+1);
        end

        % For subframes where PMI/CQI were not updated, re-use the previous
        % values in the buffers
        if (~cqipmiUpdate)
            pmiBuffer(:,cqipmiPtr+1) = ...
                pmiBuffer(:,mod(cqipmiPtr-1,cqipmiDelay)+1);
            cqiBuffer(:,cqipmiPtr+1) = ...
                cqiBuffer(:,mod(cqipmiPtr-1,cqipmiDelay)+1);
        end
        
    end

    % Display results for the current RI configuration
    fprintf('\nResults with %s:\n',riConfigStr);
    tputTotal = sum(totalTBS);
    if (riConfig==1)
        
        % Compute and display throughput
        tputFixedRI = sum(totalTBS.*(1-totalCRC));        
        fprintf('Throughput: %d bits (%0.2f%%)\n', ...
            tputFixedRI,tputFixedRI/tputTotal*100);
        
    else
        
        % Compute and display throughput
        tputUEReported = sum(totalTBS.*(1-totalCRC));
        fprintf('Throughput: %d bits (%0.2f%%)\n', ...
            tputUEReported,tputUEReported/tputTotal*100); 
            
        % Compute and display throughput ratio            
        measuredGamma = tputUEReported/tputFixedRI;
        fprintf(['Throughput ratio (gamma): %0.3f' ...
            ' (requirement is >= %0.2f)'],measuredGamma,targetGamma);
        
    end

end    

%% Plot Results
% Two figures are produced. The first figure has two subplots: the first
% subplot shows the reported RI for each subframe; the second subplot shows
% the reported PMI for each subframe. The second figure also has two
% subplots: the first subplot shows the estimated SINR for each subframe;
% the second subplot shows the reported CQI for each subframe. These plots
% illustrate how the RI, PMI and CQI vary over time due to the fading
% channel.

figure;
subplot(2,1,1);
plot(riXaxis,RIReport,'bo-');
xlabel('Subframe number');
ylabel('Selected RI');
title('Selected RI versus subframe number');
hold on;
subplot(2,1,2);
plot(cqipmiXaxis,PMIReport,'ks-');
xlabel('Subframe number');
ylabel('Selected PMI');
title('Selected PMI versus subframe number');

figure;
subplot(2,1,1);
plot(cqipmiXaxis,SINRReport,'rx-');
xlabel('Subframe number');
ylabel('Estimated SINR (dB)');
title('Estimated SINR (1st codeword) versus subframe number');
hold on;
subplot(2,1,2);
plot(cqipmiXaxis,CQIReport,'bo-');    
xlabel('Subframe number');
ylabel('Selected CQI Index');    
title('Selected CQI Index (1st codeword) versus subframe number');

%% Selected Bibliography
% # 3GPP TS 36.101 "User Equipment (UE) radio transmission and reception"
% # 3GPP TS 36.213 "Physical layer procedures"
