[Convert Channel class to C implementation #44]

Continue acquisition convertion, not working.

core/c_functions/channel_l1ca.c

@@ -2,6 +2,8 @@
 #include <stdlib.h>
 #include <complex.h>
 
+#include "fft8g.h"
+
 /// CONSTANTS
 // GPS' definition of Pi
 #define PI 3.1415926535898
@@ -38,6 +40,7 @@ struct st_Channel_L1CA {
     // Satellite variables
     int  satelliteID;    // PRN number
     int* code;           // PRN code sequence (size 1023 for GPS L1CA)
+    complex double* codeFFT;
 
     // Numerically Controled Oscilator (NCO)
     double NCO_code;
@@ -53,7 +56,14 @@ struct st_Channel_L1CA {
     double codeCounter;
 
     // Acquisition parameters
-
+    complex double* cohCorrelationMap;
+    double* correlationMap;
+    int acq_nbCohIntegration;
+    int acq_nbNoncohIntegration;
+    double acq_frequencyRange;
+    double acq_frequencySteps;
+    int cohIntegrationCounter;
+    int noncohIntegrationCounter;
 
     // Tracking parameters
     double* correlatorsSpacing;
@@ -106,8 +116,45 @@ void ProcessHandler(st_Channel_L1CA channelVar){
 
 // --------------------------------------------------------------------------------------------------------------------
 
-void RunAcquisition(){
+void RunAcquisition(st_Channel_L1CA channelVar){
+
+    int nbBins = channelVar.acq_frequencyRange*2 / channelVar.acq_frequencySteps + 1;
+    complex double _correlationMap[nbBins][channelVar.size];
+
+    // Perform PCPS
+    PCPS(
+        channelVar.rfData, channelVar.codeFFT, channelVar.size, channelVar.acq_frequencyRange, 
+        channelVar.acq_frequencySteps, channelVar.samplingFrequency, _correlationMap
+    );
+
+    // Perform coherent / non-coherent integration
+    for(int i=0; i < nbBins; i++){
+        for(int j=0; j < channelVar.size; j++){
+
+            // Coherent integration
+            channelVar.cohCorrelationMap[i][j] += _correlationMap[i][j];
+
+            if (channelVar.cohIntegrationCounter+1 == channelVar.acq_nbCohIntegration){
+                // Non-coherent integration
+                channelVar.correlationMap[i][j] += abs(channelVar.cohCorrelationMap[i][j]);
+
+                // Reset variables
+                channelVar.cohCorrelationMap[i][j] = 0.0;
+                channelVar.cohIntegrationCounter = 0;
+            }
+            else{
+                channelVar.cohIntegrationCounter++;
+            }
+        }
+    }
+
+    // Perform non-coherent integration
+    if(channelVar.cohIntegrationCounter >= channelVar.acq_nbCohIntegration){
+
+    }
+
 
+    return;
 }
 
 // --------------------------------------------------------------------------------------------------------------------
@@ -158,30 +205,50 @@ void RunTracking(st_Channel_L1CA channelVar){
 
 // --------------------------------------------------------------------------------------------------------------------
 
-void PCPS(complex double* rfData, double* code, size_t size, double freqRange, double freqstep, int cohIntegration, 
-          int nonCohIntegration, double samplesPerCode, double samplingFreq){
+void PCPS(complex double* rfData, complex double* codeFFT, size_t size, double freqRange, double freqStep, 
+          double samplingFreq, complex double* r_correlationMap){
 
-    size_t size_phasePoints = cohIntegration * samplesPerCode];
-    double* phasePoints[size_phasePoints];    
+    int nbBins = freqRange*2 / freqStep + 1;
+    double phasePoints[size];
+    complex double signal[size];  
+    double freq = 0.0;
 
     // Compute the phase points
-    for(int i=0; i < cohIntegration * samplesPerCode; i++){
-        phasePoints[i] = i * 2 * PI / samplingFreq
+    for(int i=0; i < size; i++){
+        phasePoints[i] = i * 2 * PI / samplingFreq;
     }
 
     // Compute frequency bins
-    int i = 0;
-    double freq = 0.0
-    while(freq < freqRange){
-        freq = -freqRange + freqStep * i
+    for(int i = 0; i < nbBins; i++){
+        freq = -freqRange + freqStep * i;
 
-        // Generate replica
-        for(int i=0; i <size_phasePoint< i++){
+        // Generate replica and mix
+        for(int i=0; i <size; i++){
 
+            // Generate replica
+            signal[i] = cexp(I * freq * phasePoints[i]);
+
+            // Mix with carrier
+            signal[i] *= rfData[i];
         }
-    }
 
+        // Perform FFT (in-place) 
+        fft(signal, size, 0);
 
+        // Perform correlation with code
+        // TODO check if size corrects
+        for(int i=0; i < size; i++){
+            signal[i] *= codeFFT[i];
+        }
+
+        // Perform IFFT (in-place)
+        fft(signal, size, 1);
+
+        // Save results
+        r_correlationMap[i] = *signal;
+    }
+
+    return;
 }
 
 // --------------------------------------------------------------------------------------------------------------------

core/c_functions/fft8g.h

@@ -1189,3 +1189,36 @@ void rftbsub(int n, double *a, int nc, double *c)
     }
     a[m + 1] = -a[m + 1];
 }
+
+/// ===================================================================================================================
+// Added by Hans
+/*
+@brief Pad an input signal to a power of two length and perform the FFT/IFFT on it.
+@param signal Array with a real-valued signal (will be over-written with complex values).
+@param size   Size of the array.
+@param doInv  Whether to perform forward or inverse FFT.
+@return void.
+*/
+void fft(double* signal, const size_t size, const char doInv)
+{
+    // Compute the size to pad the signal to
+    size_t paddedSize = 1;
+    while (paddedSize < size)
+        paddedSize <<= 1;
+
+    // Pad the signal with zeroes
+    double paddedSignal[paddedSize];
+    for (size_t i = 0; i < size; i++)
+        paddedSignal[i] = signal[i];
+    for (size_t i = size; i < paddedSize; i++)
+        paddedSignal[i] = signal[i];
+
+    // Perform the FFT in place
+    rdft(paddedSize, doInv == 0 ? 1 : -1, paddedSignal, ip, w);
+
+    // Return the result in the same array
+    for (size_t i = 0; i < size; i++)
+        signal[i] = paddedSignal[i];
+
+    return;
+}