function acqResults = acquisitionPostProcessing( acqResults, settings, ...
                                                 resultsNoncoherent )

%==============================================================================%
%                                                                              %
%   acqResults = acquisitionPostProcessing( acqResults, settings, ...          %
%                                           resultsNoncoherent )               %
%                                                                              %
%   Post-processing to calculate acquisition results                           %
%   David S. De Lorenzo                                                        %
%   28 January 2008                                                            %
%                                                                              %
%==============================================================================%
%                       GPS Software Receiver in Matlab                        %
%     Copyright (c) 2003-2017 David S. De Lorenzo / 2018-2023 Trimble Inc.     %
%==============================================================================%

%--- Number of frequency bins for the given Doppler search range ---------------
numberOfDopplerBins = length(acqResults.dopplerFreqVector);

%--- Initialize temporary storage for acquisition results output ---------------
codePhaseVector = zeros(1, numberOfDopplerBins);
peakSizeVector  = zeros(1, numberOfDopplerBins);
cppmVector      = zeros(1, numberOfDopplerBins);
cpprVector      = zeros(1, numberOfDopplerBins);
stdDevVector    = zeros(1, numberOfDopplerBins);

%--- Loop over Doppler frequency -----------------------------------------------
for dopplerSearchBin = 1:numberOfDopplerBins
    
    %--- Find the value of the maximum correlation peak ------------------------
    [peakSize,codePhaseIndex] = max(resultsNoncoherent(dopplerSearchBin,:));
    peakSizeVector(dopplerSearchBin) = peakSize;
    
    %--- Find peak with "3-point neighborhood" method --------------------------
    %      --> This is ~3x more precise than using index from max() call.
    %    The code below corrects the code-phase estimate for the peak not
    %    falling exactly on an integer sample number.  This fractional
    %    code-phase correction only really matters when the sampling
    %    frequency is close to the PRN code chipping rate, so that
    %    estimating the acquisition code-phase from the max() of the
    %    correlation function could lead to initialization problems in
    %    the DLL.  For example, if the sampling frequency is 2.5 MHz and
    %    the code chipping frequency is 1.023 Mchips/sec (meaning that
    %    there are ~2.5 samples per chip), then a maximum code-phase
    %    error of 1/5 chip could result - more if the signal is noisy.
    if (codePhaseIndex == settings.samplesPer1ms)
        peak1 = abs(resultsNoncoherent(dopplerSearchBin, end-1));
        peak2 = abs(resultsNoncoherent(dopplerSearchBin, 1));
    elseif (codePhaseIndex == 1)
        peak1 = abs(resultsNoncoherent(dopplerSearchBin, end));
        peak2 = abs(resultsNoncoherent(dopplerSearchBin, 2));
    else
        peak1 = abs(resultsNoncoherent(dopplerSearchBin, codePhaseIndex-1));
        peak2 = abs(resultsNoncoherent(dopplerSearchBin, codePhaseIndex+1));
    end
    if (peak2 > peak1)
        codePhaseInterp = codePhaseIndex + 1/2 * (peak2-peak1)/(peakSize-peak1);
    else
        codePhaseInterp = codePhaseIndex - 1/2 * (peak1-peak2)/(peakSize-peak2);
    end
    
    %--- Find the code-phase of the correlation peak ---------------------------
    codePhaseVector(dopplerSearchBin) = ...
        (codePhaseInterp - 1) * settings.codeFreq / settings.samplingFreq;
    
    %--- Find 2 chip wide code-phase exclude range around the peak -------------
    excludeRangeIndex1 = ...
        codePhaseIndex - 2 * ceil(settings.samplingFreq / settings.codeFreq);
    excludeRangeIndex2 = ...
        codePhaseIndex + 2 * ceil(settings.samplingFreq / settings.codeFreq);
    
    %--- Correct exclude range if it includes either array boundary ------------
    if (excludeRangeIndex1 < 1)
        codePhaseRange = ...
            excludeRangeIndex2 : (settings.samplesPer1ms + excludeRangeIndex1);
    elseif (excludeRangeIndex2 > settings.samplesPer1ms)
        codePhaseRange = ...
            (excludeRangeIndex2 - settings.samplesPer1ms) : excludeRangeIndex1;
    else
        codePhaseRange = ...
            [1:excludeRangeIndex1, excludeRangeIndex2:settings.samplesPer1ms];
    end
    
    %--- Find the second highest peak ------------------------------------------
    secondPeakSize = max(resultsNoncoherent(dopplerSearchBin,codePhaseRange));
    
    %--- Find the mean correlation value ---------------------------------------
    meanPeakSize = mean(resultsNoncoherent(dopplerSearchBin, codePhaseRange));
    
    %--- Find the correlation-peak ratios --------------------------------------
    %    cppm - correlation peak to mean peak ratio
    %    cppr - correlation peak to next peak ratio
    cppmVector(dopplerSearchBin) = peakSize / meanPeakSize;
    cpprVector(dopplerSearchBin) = peakSize / secondPeakSize;
    
    %--- Find the standard deviation -------------------------------------------
    %    This excludes the peak value.
    stdDevVector(dopplerSearchBin) = ...
        std(resultsNoncoherent(dopplerSearchBin, codePhaseRange));
    
end % loop over Doppler frequency

%--- Store results with the highest correlation peak to mean peak ratio --------
[~,dopplerResultBin]        = max(cppmVector);
acqResults.codePhase        = codePhaseVector(dopplerResultBin);
acqResults.dopplerFreq      = acqResults.dopplerFreqVector(dopplerResultBin);
acqResults.dopplerFreqIndex = dopplerResultBin;
acqResults.cppm             = cppmVector(dopplerResultBin);
acqResults.cppr             = cpprVector(dopplerResultBin);
acqResults.stdDev           = stdDevVector(dopplerResultBin);
acqResults.codePhaseVector  = codePhaseVector;
acqResults.peakSizeVector   = peakSizeVector;
acqResults.cppmVector       = cppmVector;
acqResults.cpprVector       = cpprVector;
acqResults.stdDevVector     = stdDevVector;
acqResults.resultsVector    = resultsNoncoherent(dopplerResultBin,:);

return
%===============================================================================