Added amplitude envelope to sweep signals, with associated methods

2024-11-01 14:17:31 +01:00 · 2024-11-01 14:17:31 +01:00 · 9caf5fe387
commit 9caf5fe387
parent a38eca47f3
4 changed files with 111 additions and 39 deletions
--- a/cpp_src/dsp/lasp_siggen.h
+++ b/cpp_src/dsp/lasp_siggen.h
@ -73,7 +73,7 @@ public:
   * @brief Mute the signal. Passes through the DC offset. No lock is hold. If
   * it just works one block later, than that is just the case.
   *
-   * @param mute if tre
+   * @param mute if true
   */
  void setMute(bool mute = true) { _muted = mute; _interruption_frame_count=0; }

--- a/cpp_src/dsp/lasp_siggen_impl.cpp
+++ b/cpp_src/dsp/lasp_siggen_impl.cpp
@ -24,60 +24,84 @@ DEBUGTRACE_VARIABLES;

 Noise::Noise(){DEBUGTRACE_ENTER}

-vd Noise::genSignalUnscaled(us nframes) {
+vd Noise::genSignalUnscaled(us nframes)
+{
  return arma::randn<vd>(nframes);
 }
 void Noise::resetImpl() {}

 Sine::Sine(const d freq) : omg(2 * arma::datum::pi * freq) { DEBUGTRACE_ENTER; }

-vd Sine::genSignalUnscaled(const us nframes) {
+vd Sine::genSignalUnscaled(const us nframes)
+{
  /* DEBUGTRACE_ENTER; */
  slock lck(_mtx);
  const d pi = arma::datum::pi;
  vd phase_vec =
      arma::linspace(phase, phase + omg * (nframes - 1) / _fs, nframes);
  phase += omg * nframes / _fs;
-  while (phase > 2 * arma::datum::pi) {
+  while (phase > 2 * arma::datum::pi)
+  {
    phase -= 2 * pi;
  }
  return arma::sin(phase_vec);
 }

-vd Periodic::genSignalUnscaled(const us nframes) {
+vd Periodic::genSignalUnscaled(const us nframes)
+{
  vd res(nframes);
  slock lck(_mtx);
-  if (_signal.size() == 0) {
+  if (_signal.size() == 0)
+  {
    throw rte("No signal defined while calling");
  }
-  for (us i = 0; i < nframes; i++) {
-    res(i) = _signal[_cur_pos];
+  if (_signal.size() != A_.size())
+  {
+    std::cout << "Seq size: " << _signal.size() << ", A size: " << A_.size() << "\n";
+    throw rte("Sequence and amplitude envelopes have different lengths");
+  }
+  for (us i = 0; i < nframes; i++)
+  {
+    res(i) = A_[_cur_pos] * _signal[_cur_pos];
    _cur_pos++;
    _cur_pos %= _signal.size();
  }
  return res;
 }

+void Periodic::setA(const vd &A)
+{
+  A_ = A;
+}
+
 Sweep::Sweep(const d fl, const d fu, const d Ts, const d Tq, const us flags)
-    : fl_(fl), fu_(fu), Ts(Ts), Tq(Tq), flags(flags) {
-  if (fl <= 0 || fu < fl || Ts <= 0) {
+    : fl_(fl), fu_(fu), Ts(Ts), Tq(Tq), flags(flags)
+{
+  if (fl <= 0 || fu < fl || Ts <= 0)
+  {
    throw rte("Invalid sweep parameters");
  }
-  if ((flags & ForwardSweep) && (flags & BackwardSweep)) {
+  if ((flags & ForwardSweep) && (flags & BackwardSweep))
+  {
    throw rte(
        "Both forward and backward sweep flag set. Please only set either one "
        "or none for a continuous sweep");
  }
-  if ((flags & LinearSweep) && (flags & LogSweep)) {
+  if ((flags & LinearSweep) && (flags & LogSweep))
+  {
    throw rte(
        "Both logsweep and linear sweep flag set. Please only set either one.");
  }
-  if (!((flags & LinearSweep) || (flags & LogSweep))) {
+  if (!((flags & LinearSweep) || (flags & LogSweep)))
+  {
    throw rte("Either LinearSweep or LogSweep should be given as flag");
  }
+
+  resetImpl();
 }

-void Sweep::resetImpl() {
+void Sweep::resetImpl()
+{
  DEBUGTRACE_ENTER;
  slock lck(_mtx);

@ -94,14 +118,18 @@ void Sweep::resetImpl() {
  const us N = Ns + Nq;

  _signal = vd(N, arma::fill::zeros);
+  fn_ = vd(N, arma::fill::zeros);
  index = 0;

  d fl, fu;
  /* Swap fl and fu for a backward sweep */
-  if (backward_sweep) {
+  if (backward_sweep)
+  {
    fu = fl_;
    fl = fu_;
-  } else {
+  }
+  else
+  {
    /* Case of continuous sweep, or forward sweep */
    fl = fl_;
    fu = fu_;
@ -110,8 +138,10 @@ void Sweep::resetImpl() {
  d phase = 0;

  /* Linear sweep */
-  if (flags & LinearSweep) {
-    if (forward_sweep || backward_sweep) {
+  if (flags & LinearSweep)
+  {
+    if (forward_sweep || backward_sweep)
+    {
      /* Forward or backward sweep */
      /* TRACE(15, "Forward or backward sweep"); */
      us K = (us)(Dt * (fl * Ns + 0.5 * (Ns - 1) * (fu - fl)));
@ -120,12 +150,16 @@ void Sweep::resetImpl() {
      /* iVARTRACE(15, K); */
      /* dVARTRACE(15, eps); */

-      for (us n = 0; n < Ns; n++) {
+      for (us n = 0; n < Ns; n++)
+      {
        _signal(n) = d_sin(phase);
        d fn = fl + ((d)n) / Ns * (fu + eps - fl);
+        fn_(n) = fn;
        phase += 2 * arma::datum::pi * Dt * fn;
      }
-    } else {
+    }
+    else
+    {
      /* Continous sweep */
      /* TRACE(15, "continuous sweep"); */

@ -150,18 +184,24 @@ void Sweep::resetImpl() {
      /* dVARTRACE(15, eps); */
      d phase = 0;

-      for (us n = 0; n <= Ns; n++) {
+      for (us n = 0; n <= Ns; n++)
+      {
        /* iVARTRACE(17, n); */
-        if (n < N) {
+        if (n < N)
+        {
          _signal[n] = d_sin(phase);
        }

        d fn;
-        if (n <= Nf) {
+        if (n <= Nf)
+        {
          fn = fl + ((d)n) / Nf * (fu - fl);
-        } else {
+        }
+        else
+        {
          fn = fu - ((d)n - Nf) / Nb * (fu + eps - fl);
        }
+        fn_(n) = fn;
        /* dbgassert(fn >= 0, "BUG"); */

        phase += 2 * number_pi * Dt * fn;
@ -169,9 +209,12 @@ void Sweep::resetImpl() {
      /* This should be a very small number!! */
      /* dVARTRACE(15, phase); */
    }
-  } else if (flags & LogSweep) {
+  }
+  else if (flags & LogSweep)
+  {
    DEBUGTRACE_PRINT("Log sweep");
-    if (forward_sweep || backward_sweep) {
+    if (forward_sweep || backward_sweep)
+    {
      /* Forward or backward sweep */
      DEBUGTRACE_PRINT("Forward or backward sweep");
      d k1 = (fu / fl);
@ -180,7 +223,8 @@ void Sweep::resetImpl() {

      /* Iterate k to the right solution */
      d E;
-      for (us iter = 0; iter < 10; iter++) {
+      for (us iter = 0; iter < 10; iter++)
+      {
        E = 1 + K / (Dt * fl) * (d_pow(k, 1.0 / Ns) - 1) - k;
        d dEdk = K / (Dt * fl) * d_pow(k, 1.0 / Ns) / (Ns * k) - 1;
        k -= E / dEdk;
@ -191,12 +235,16 @@ void Sweep::resetImpl() {
      DEBUGTRACE_PRINT(k);
      DEBUGTRACE_PRINT(E);

-      for (us n = 0; n < Ns; n++) {
+      for (us n = 0; n < Ns; n++)
+      {
        _signal[n] = d_sin(phase);
        d fn = fl * d_pow(k, ((d)n) / Ns);
+        fn_(n) = fn;
        phase += 2 * number_pi * Dt * fn;
      }
-    } else {
+    }
+    else
+    {
      DEBUGTRACE_PRINT("Continuous sweep");

      const us Nf = Ns / 2;
@ -212,7 +260,8 @@ void Sweep::resetImpl() {

      /* Newton iterations to converge k to the value such that the sweep is
       * continuous */
-      for (us iter = 0; iter < NITER_NEWTON; iter++) {
+      for (us iter = 0; iter < NITER_NEWTON; iter++)
+      {
        E = (k - 1) / (d_pow(k, 1.0 / Nf) - 1) +
            (k - 1) / (1 - d_pow(k, -1.0 / Nb)) - K / Dt / fl;
        DEBUGTRACE_PRINT(E);
@ -236,29 +285,37 @@ void Sweep::resetImpl() {
      DEBUGTRACE_PRINT(k);
      DEBUGTRACE_PRINT(E);

-      for (us n = 0; n <= Ns; n++) {
+      for (us n = 0; n <= Ns; n++)
+      {
        /* iVARTRACE(17, n); */
-        if (n < Ns) {
+        if (n < Ns)
+        {
          _signal[n] = d_sin(phase);
        }

        d fn;
-        if (n <= Nf) {
+        if (n <= Nf)
+        {
          fn = fl * d_pow(k, ((d)n) / Nf);
-        } else {
+        }
+        else
+        {
          fn = fl * k * d_pow(1 / k, ((d)n - Nf) / Nb);
        }
+        fn_(n) = fn;
        /* dbgassert(fn >= 0, "BUG"); */

        phase += 2 * number_pi * Dt * fn;
-        while (phase > 2 * number_pi) phase -= 2 * number_pi;
+        while (phase > 2 * number_pi)
+          phase -= 2 * number_pi;
        /* dVARTRACE(17, phase); */
      }
      /* This should be a very small number!! */
      DEBUGTRACE_PRINT(phase);
    }
  } // End of log sweep
-  else {
+  else
+  {
    // Either log or linear sweep had to be given as flags.
    assert(false);
  }
--- a/cpp_src/dsp/lasp_siggen_impl.h
+++ b/cpp_src/dsp/lasp_siggen_impl.h
@ -58,6 +58,7 @@ class Sine : public Siggen {
 * periodic as the frequency can be any floating point value.
 */
 class Periodic: public Siggen {
+  
  protected:
    vd _signal { 1, arma::fill::zeros};
    us _cur_pos = 0;
@ -68,8 +69,12 @@ class Periodic: public Siggen {
     * @return As stated above
     */
    vd getSequence() const { return _signal; }
+    vd A_ { 1, arma::fill::ones};
+  
+    void setA(const vd& A);

    virtual vd genSignalUnscaled(const us nframes) override final;
+
    ~Periodic() = default;

 };
@ -81,6 +86,7 @@ class Sweep : public Periodic {
  d fl_, fu_, Ts, Tq;
  us index;
  us flags;
+  vd fn_ { 1, arma::fill::zeros};  
  
  void resetImpl() override;

@ -90,6 +96,8 @@ class Sweep : public Periodic {
  static constexpr int LinearSweep = 1 << 2;
  static constexpr int LogSweep = 1 << 3;

+  vd getfn() const { return fn_; }
+
  /**
   * Create a sweep signal
   *
--- a/cpp_src/pybind11/lasp_siggen.cpp
+++ b/cpp_src/pybind11/lasp_siggen.cpp
@ -43,11 +43,18 @@ void init_siggen(py::module &m) {

  py::class_<Periodic, std::shared_ptr<Periodic>> periodic(m, "Periodic",
                                                           siggen);
+  periodic.def("setA", 
+               [](Periodic &p, const dpyarray A) {
+                p.setA(NpyToCol<d, false>(A));
+                });
+
  periodic.def("getSequence",
               [](const Sweep &s) { return ColToNpy<d>(s.getSequence()); });

  py::class_<Sweep, std::shared_ptr<Sweep>> sweep(m, "Sweep", periodic);
  sweep.def(py::init<const d, const d, const d, const d, const us>());
+  sweep.def("getfn", 
+               [](const Sweep &s) { return ColToNpy<d>(s.getfn()); });
  sweep.def_readonly_static("ForwardSweep", &Sweep::ForwardSweep);
  sweep.def_readonly_static("BackwardSweep", &Sweep::BackwardSweep);
  sweep.def_readonly_static("LinearSweep", &Sweep::LinearSweep);