Made PPM meter range-aware. Improved PPM-code. Made RtAps sensitivity-aware. Improved code to get ranges for each channel from Daq object

src/lasp/device/lasp_daq.cpp

@@ -75,14 +75,16 @@ const DataTypeDescriptor &Daq::dtypeDescr() const {
   return dtype_map.at(dataType());
 }
 
-double Daq::inputRangeForChannel(us ch) const {
-  if (!(ch < ninchannels)) {
-    throw rte("Invalid channel number");
+dvec Daq::inputRangeForEnabledChannels(const bool include_monitor) const {
+  dvec res;
+  auto chs = enabledInChannels(include_monitor);
+  for(auto& ch : chs) {
+    res.push_back(availableInputRanges.at(ch.rangeIndex));
   }
-  return availableInputRanges.at(inchannel_config[ch].rangeIndex);
+  return res;
 }
 
-us Daq::neninchannels(bool include_monitorchannel) const {
+us Daq::neninchannels(const bool include_monitorchannel) const {
   boolvec eninchannels = this->eninchannels(include_monitorchannel);
   return std::count(eninchannels.cbegin(), eninchannels.cend(), true);
 }

src/lasp/device/lasp_daq.h

@@ -136,7 +136,7 @@ public:
    *
    * @return Boolean vector
    */
-  boolvec eninchannels(bool include_monitor = true) const;
+  boolvec eninchannels(const bool include_monitor = true) const;
 
   /**
    * @brief Create an array of booleans for each enabled output channel.
@@ -153,15 +153,15 @@ public:
   double samplerate() const;
 
   /**
-   * @brief Returns the input range for each channel. Means the minimum from
-   * the absolute value of the minumum and maximum value that is allowed to
-   * pass unclipped.
-   *
-   * @param ch The channel index from the input channels.
+   * @brief Returns the input range for each channel. Maximum allowed absolute
+   * value of the signal that can pass unclipped.
    *
+   * @param include_monitor If set to true and a monitor channel is available,
+   * the first index in the array will correspond to the monitor channel.
+   * *
    * @return Maximum offset from 0 before clipping.
    */
-  double inputRangeForChannel(us ch) const;
+  dvec inputRangeForEnabledChannels(const bool include_monitor=true) const;
 
   /**
    * @brief Returns datatype (enum) corresponding to the datatype of the

src/lasp/device/lasp_daqconfig.cpp

@@ -86,9 +86,9 @@ int DaqConfiguration::getLowestEnabledOutChannel() const {
   }
   return -1;
 }
-vector<DaqChannel> DaqConfiguration::enabledInChannels() const {
+vector<DaqChannel> DaqConfiguration::enabledInChannels(const bool include_monitor) const {
   vector<DaqChannel> res;
-  if(monitorOutput) {
+  if(monitorOutput && include_monitor) {
     DaqChannel ch;
     ch.name = "Internal output monitor (loopback)";
     ch.enabled = true;

src/lasp/device/lasp_daqconfig.h

@@ -209,9 +209,9 @@ public:
    *
    * @return true if equal
    */
-  bool operator==(const DaqChannel& other) const {
-    return other.name == name && other.sensitivity == sensitivity && 
-      other.qty == qty;
+  bool operator==(const DaqChannel &other) const {
+    return other.name == name && other.sensitivity == sensitivity &&
+           other.qty == qty;
   }
 };
 
@@ -250,9 +250,13 @@ public:
   /**
    * @brief Return list of enabled input channels
    *
+   * @param include_monitor If set to true and a monitor channel is available,
+   * the first indices in the array will correspond to the monitor channel(s).
+   *
    * @return That.
    */
-  std::vector<DaqChannel> enabledInChannels() const;
+  std::vector<DaqChannel>
+  enabledInChannels(const bool include_monitor = true) const;
 
   /**
    * @brief Channel configuration for output channels
@@ -303,8 +307,8 @@ public:
   bool match(const DeviceInfo &devinfo) const;
 
   /**
-   * @brief Get the enabled highest channel number from the list of enabled input
-   * channels.
+   * @brief Get the enabled highest channel number from the list of enabled
+   * input channels.
    *
    * @return Index to the highest input channel. -1 if no input channels are
    * enabled.
@@ -342,7 +346,7 @@ public:
    * @param val true if enabled, false if disabled.
    */
   void setAllInputEnabled(bool val) {
-    for(auto& ch: inchannel_config) {
+    for (auto &ch : inchannel_config) {
       ch.enabled = val;
     }
   }
@@ -353,7 +357,7 @@ public:
    * @param val true if enabled, false if disabled.
    */
   void setAllOutputEnabled(bool val) {
-    for(auto& ch: outchannel_config) {
+    for (auto &ch : outchannel_config) {
       ch.enabled = val;
     }
   }

src/lasp/device/lasp_uldaq.cpp

@@ -182,13 +182,14 @@ public:
     boolvec eninchannels_without_mon = daq.eninchannels(false);
 
     // Initialize input, if any
+    dvec ranges = daq.inputRangeForEnabledChannels(false);
     for (us chin = 0; chin < 4; chin++) {
       if (eninchannels_without_mon[chin] == true) {
         DaqInChanDescriptor indesc;
         indesc.type = DAQI_ANALOG_SE;
         indesc.channel = chin;
 
-        double rangeval = daq.inputRangeForChannel(chin);
+        double rangeval = ranges.at(chin);
         Range rangenum;
         if (fabs(rangeval - 1.0) < 1e-8) {
           rangenum = BIP1VOLTS;

src/lasp/dsp/lasp_ppm.cpp

@@ -21,7 +21,7 @@ bool PPMHandler::inCallback_threaded(const DaqData &d) {
   Lck lck(_mtx);
 
   dmat data = d.toFloat();
-  /* data.print(); */
+
 
   const us nchannels = d.nchannels;
   assert(data.n_cols == nchannels);
@@ -31,26 +31,20 @@ bool PPMHandler::inCallback_threaded(const DaqData &d) {
     _cur_max = vd(nchannels, arma::fill::value(1e-80));
     _clip_time = vd(nchannels, arma::fill::value(-1));
   }
+
   assert(_clip_time.size() == _cur_max.size());
 
-  /// Obtain max abs values
-  vd maxabs = arma::max(arma::abs(data), 0).as_col();
-  /* maxabs.print(); */
+  /// Obtain max abs values, and scale with range for each channel
+  vd maxabs = arma::max(arma::abs(data), 0).as_col() / _max_range;
 
-  arma::uvec clip_indices = arma::find(maxabs > clip_point);
-  arma::uvec clip(nchannels, arma::fill::zeros);
-  clip.elem(clip_indices).fill(1);
+  /// Find indices for channels that have a clip
+  arma::uvec clips = maxabs > clip_point;
 
-  arma::uvec update_max_idx = arma::find(maxabs > _cur_max);
-  /* update_max_idx.print(); */
-  arma::uvec update_max(nchannels, arma::fill::zeros);
-  /* update_max.print(); */
-  update_max.elem(update_max_idx).fill(1);
-
-  assert(_cur_max.size() == _clip_time.size());
+  /// Find channels where the new maximum is higher than the previous one
+  arma::uvec update_max = maxabs > _cur_max;
 
   for (us i = 0; i < nchannels; i++) {
-    if (clip(i)) {
+    if (clips(i)) {
       /// Reset clip counter
       _clip_time(i) = 0;
     } else if (_clip_time(i) > clip_indication_time) {
@@ -90,6 +84,15 @@ void PPMHandler::reset(const Daq *daq) {
 
     _cur_max.fill(1e-80);
 
+    const us nchannels = daq->neninchannels();
+    _max_range.resize(nchannels);
+
+    dvec ranges = daq->inputRangeForEnabledChannels();
+    assert(ranges.size() == nchannels);
+    for(us i=0;i<daq->neninchannels();i++) {
+      _max_range[i] = ranges[i];
+    }
+
     _clip_time.fill(-1);
 
     const d fs = daq->samplerate();

src/lasp/dsp/lasp_ppm.h

@@ -60,6 +60,11 @@ class PPMHandler: public ThreadedInDataHandler {
    */
   vd _clip_time;
 
+  /**
+   * @brief Storage for maximum values
+   */
+  vd _max_range;
+
   public:
   /**
    * @brief Constructs Peak Programme Meter

src/lasp/dsp/lasp_rtaps.cpp

@@ -8,17 +8,17 @@ using std::endl;
 using Lck = std::scoped_lock<std::mutex>;
 
 RtAps::RtAps(StreamMgr &mgr, const Filter *freqWeightingFilter,
-             const us nfft,
-             const Window::WindowType w,
-             const d overlap_percentage, const d time_constant)
-    : ThreadedInDataHandler(mgr),
-      _ps(nfft, w, overlap_percentage, time_constant) {
+    const us nfft,
+    const Window::WindowType w,
+    const d overlap_percentage, const d time_constant)
+  : ThreadedInDataHandler(mgr),
+  _ps(nfft, w, overlap_percentage, time_constant) {
+
+    if (freqWeightingFilter != nullptr) {
+      _filterPrototype = freqWeightingFilter->clone();
+    }
 
-  if (freqWeightingFilter != nullptr) {
-    _filterPrototype = freqWeightingFilter->clone();
   }
-
-}
 RtAps::~RtAps() {
   Lck lck(_ps_mtx);
   stop();
@@ -27,9 +27,15 @@ bool RtAps::inCallback_threaded(const DaqData &data) {
 
   DEBUGTRACE_ENTER;
 
+  Lck lck(_ps_mtx);
 
   dmat fltdata = data.toFloat();
   const us nchannels = fltdata.n_cols;
+  if(nchannels != _sens.size()) {
+    cerr << "**** Error: sensitivity size does not match! *****" << endl;
+    return false;
+  }
+  fltdata.each_row() %= _sens.as_row();
 
   if (_filterPrototype) {
 
@@ -45,7 +51,7 @@ bool RtAps::inCallback_threaded(const DaqData &data) {
     }
 
     // Apply filtering
-  /* #pragma omp parallel for */
+    #pragma omp parallel for
     for (us i = 0; i < nchannels; i++) {
       vd col = fltdata.col(i);
       _freqWeightingFilters.at(i)->filter(col);
@@ -53,24 +59,32 @@ bool RtAps::inCallback_threaded(const DaqData &data) {
     }
   } // End of if(_filterPrototype)
 
-  Lck lck(_ps_mtx);
   _ps.compute(fltdata);
 
   return true;
 }
-void RtAps::reset(const Daq *daq) { // Explicitly say
-                                    // to GCC that
-                                    // the argument is
-                                    // not used.
+void RtAps::reset(const Daq *daq) {
 
   DEBUGTRACE_ENTER;
   Lck lck(_ps_mtx);
+  for (auto &filter : _freqWeightingFilters) {
+    filter->reset();
+  }
+  if(daq) {
+    _sens.resize(daq->neninchannels());
+    us i = 0;
+    for(const auto& ch: daq->enabledInChannels()) {
+      _sens[i] = ch.sensitivity;
+      i++;
+    }
+  }
+
   _ps.reset();
 }
 
 ccube RtAps::getCurrentValue() const {
 
-  /* DEBUGTRACE_ENTER; */
+  DEBUGTRACE_ENTER;
   Lck lck(_ps_mtx);
   return _ps.get_est();
 

src/lasp/dsp/lasp_rtaps.h

@@ -24,6 +24,11 @@ class RtAps : public ThreadedInDataHandler {
   std::unique_ptr<Filter> _filterPrototype;
   std::vector<std::unique_ptr<Filter>> _freqWeightingFilters;
 
+  /**
+   * @brief Storage for sensitivity values
+   */
+  vd _sens;
+
   /**
    * @brief Mutex only for _ps. Other members are only accessed in thread.
    */

src/lasp/pybind11/lasp_daqconfig.cpp

@@ -69,7 +69,9 @@ void init_daqconfiguration(py::module &m) {
       .value("UserDefined", DaqChannel::Qty::UserDefined);
   daqchannel.def_readwrite("qty", &DaqChannel::qty);
 
-  daqchannel.def("__eq__", [](const DaqChannel& a, const DaqChannel& b) { return a==b;});
+  daqchannel.def("__eq__", [](const DaqChannel &a, const DaqChannel &b) {
+    return a == b;
+  });
 
   /// DaqConfiguration
   daqconfig.def(py::init<>());
@@ -95,5 +97,6 @@ void init_daqconfiguration(py::module &m) {
                           &DaqConfiguration::outchannel_config);
   daqconfig.def("setAllInputEnabled", &DaqConfiguration::setAllInputEnabled);
   daqconfig.def("setAllOutputEnabled", &DaqConfiguration::setAllOutputEnabled);
-  daqconfig.def("enabledInChannels", &DaqConfiguration::enabledInChannels);
+  daqconfig.def("enabledInChannels", &DaqConfiguration::enabledInChannels,
+                py::arg("include_monitor") = true);
 }