Calculate PAPR in the OfdmGenerator

6 files changed, 35 insertions, 13 deletions

diff --git a/Makefile.am b/Makefile.am
index 54ae62e..1880acd 100644
--- a/Makefile.am
+++ b/Makefile.am
@@ -136,6 +136,8 @@ odr_dabmod_SOURCES  = src/DabMod.cpp \
 					  src/Utils.h \
 					  src/TII.cpp \
 					  src/TII.h \
+					  src/PAPRStats.cpp \
+					  src/PAPRStats.h \
 					  src/zmq.hpp \
 					  lib/crc.h \
 					  lib/crc.c \
diff --git a/doc/example.ini b/doc/example.ini
index ec0525c..b3faffd 100644
--- a/doc/example.ini
+++ b/doc/example.ini
@@ -141,11 +141,11 @@ rate=2048000
 enable=0
 
 ; At what amplitude the signal should be clipped
-clip=70.0
+clip=50.0
 
 ; How much to clip the error signal used to compensate the effect
 ; of clipping
-error_clip=0.05
+error_clip=0.1
 
 [firfilter]
 ; The FIR Filter can be used to create a better spectral quality.
diff --git a/src/OfdmGenerator.cpp b/src/OfdmGenerator.cpp
index 915d568..1a819ee 100644
--- a/src/OfdmGenerator.cpp
+++ b/src/OfdmGenerator.cpp
@@ -27,11 +27,8 @@
 #include "OfdmGenerator.h"
 #include "PcDebug.h"
 
-#include <complex>
-#include "fftw3.h"
 #define FFT_TYPE fftwf_complex
 
-#include <stdio.h>
 #include <string.h>
 #include <stdexcept>
 #include <assert.h>
@@ -56,7 +53,10 @@ OfdmGenerator::OfdmGenerator(size_t nbSymbols,
     myCfr(enableCfr),
     myCfrClip(cfrClip),
     myCfrErrorClip(cfrErrorClip),
-    myCfrFft(nullptr)
+    myCfrFft(nullptr),
+    // Initialise the PAPRStats to a few seconds worth of samples
+    myPaprBeforeCFR(nbSymbols * 50),
+    myPaprAfterCFR(nbSymbols * 50)
 {
     PDEBUG("OfdmGenerator::OfdmGenerator(%zu, %zu, %zu, %s) @ %p\n",
             nbSymbols, nbCarriers, spacing, inverse ? "true" : "false", this);
@@ -71,6 +71,7 @@ OfdmGenerator::OfdmGenerator(size_t nbSymbols,
     RC_ADD_PARAMETER(clip, "CFR: Clip to amplitude");
     RC_ADD_PARAMETER(errorclip, "CFR: Limit error");
     RC_ADD_PARAMETER(clip_stats, "CFR: statistics (clip ratio, errorclip ratio)");
+    RC_ADD_PARAMETER(papr, "PAPR measurements (before CFR, after CFR)");
 
     if (inverse) {
         myPosDst = (nbCarriers & 1 ? 0 : 1);
@@ -203,19 +204,25 @@ int OfdmGenerator::process(Buffer* const dataIn, Buffer* dataOut)
 
         fftwf_execute(myFftPlan); // IFFT from myFftIn to myFftOut
 
+        complexf *symbol = reinterpret_cast<complexf*>(myFftOut);
+        myPaprBeforeCFR.process_block(symbol, mySpacing);
+
         if (myCfr) {
             if (myMERCalcIndex == i) {
                 before_cfr.resize(mySpacing);
                 memcpy(before_cfr.data(), myFftOut, mySpacing * sizeof(FFT_TYPE));
             }
 
-            complexf *symbol = reinterpret_cast<complexf*>(myFftOut);
             /* cfr_one_iteration runs the myFftPlan again at the end, and
              * therefore writes the output data to myFftOut.
              */
             const auto stat = cfr_one_iteration(symbol, reference.data());
 
             // i == 0 always zero power, so the MER ends up being NaN
+            if (i > 0) {
+                myPaprAfterCFR.process_block(symbol, mySpacing);
+            }
+
             if (i > 0 and myMERCalcIndex == i) {
                 /* MER definition, ETSI ETR 290, Annex C
                  *
@@ -358,8 +365,8 @@ void OfdmGenerator::set_parameter(const std::string& parameter,
     else if (parameter == "errorclip") {
         ss >> myCfrErrorClip;
     }
-    else if (parameter == "clip_stats") {
-        throw ParameterError("Parameter 'clip_stats' is read-only");
+    else if (parameter == "clip_stats" or parameter == "papr") {
+        throw ParameterError("Parameter '" + parameter + "' is read-only");
     }
     else {
         stringstream ss_err;
@@ -406,6 +413,15 @@ const std::string OfdmGenerator::get_parameter(const std::string& parameter) con
                 "MER after CFR: " << avg_mer << " dB";
         }
     }
+    else if (parameter == "papr") {
+        const double papr_before = myPaprBeforeCFR.calculate_papr();
+        const double papr_after = myPaprAfterCFR.calculate_papr();
+
+        ss << "PAPR [dB]: " << std::fixed <<
+            (papr_before == 0 ? string("N/A") : to_string(papr_before)) <<
+            ", " <<
+            (papr_after == 0 ? string("N/A") : to_string(papr_after));
+    }
     else {
         ss << "Parameter '" << parameter <<
             "' is not exported by controllable " << get_rc_name();
diff --git a/src/OfdmGenerator.h b/src/OfdmGenerator.h
index 008d84e..8374cf9 100644
--- a/src/OfdmGenerator.h
+++ b/src/OfdmGenerator.h
@@ -30,11 +30,11 @@
 #   include "config.h"
 #endif
 
-#include "porting.h"
 #include "ModPlugin.h"
 #include "RemoteControl.h"
+#include "PAPRStats.h"
 #include "fftw3.h"
-#include <sys/types.h>
+#include <cstddef>
 #include <vector>
 #include <complex>
 
@@ -102,6 +102,10 @@ class OfdmGenerator : public ModCodec, public RemoteControllable
         std::deque<double> myClipRatios;
         std::deque<double> myErrorClipRatios;
 
+        // Measure PAPR before and after CFR
+        PAPRStats myPaprBeforeCFR;
+        PAPRStats myPaprAfterCFR;
+
         size_t myMERCalcIndex = 0;
         std::deque<double> myMERs;
 };
diff --git a/src/PAPRStats.cpp b/src/PAPRStats.cpp
index bf6acc4..2e8bca9 100644
--- a/src/PAPRStats.cpp
+++ b/src/PAPRStats.cpp
@@ -69,7 +69,7 @@ void PAPRStats::process_block(const complexf* data, size_t data_len)
     }
 }
 
-double PAPRStats::calculate_papr()
+double PAPRStats::calculate_papr() const
 {
     if (m_squared_mean.size() < m_num_blocks_to_accumulate) {
         return 0;
diff --git a/src/PAPRStats.h b/src/PAPRStats.h
index 9463a3d..92e5231 100644
--- a/src/PAPRStats.h
+++ b/src/PAPRStats.h
@@ -64,7 +64,7 @@ class PAPRStats
         /* Returns PAPR in dB if enough blocks were processed, or
          * 0 otherwise.
          */
-        double calculate_papr(void);
+        double calculate_papr(void) const;
 
     private:
         size_t m_num_blocks_to_accumulate;