WDSP: split firmin into fircore, firmin and firopt

wdsp/CMakeLists.txt

@@ -23,7 +23,9 @@ set(wdsp_SOURCES
     eq.cpp
     fcurve.cpp
     fir.cpp
+    fircore.cpp
     firmin.cpp
+    firopt.cpp
     fmd.cpp
     fmmod.cpp
     fmsq.cpp
@@ -77,7 +79,9 @@ set(wdsp_HEADERS
     eq.hpp
     fcurve.hpp
     fir.hpp
+    fircore.hpp
     firmin.hpp
+    firopt.hpp
     fmd.hpp
     fmmod.hpp
     fmsq.hpp

wdsp/RXA.cpp

@@ -51,7 +51,7 @@ warren@wpratt.com
 #include "cblock.hpp"
 #include "ssql.hpp"
 #include "iir.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "wcpAGC.hpp"
 #include "anb.hpp"
 #include "nob.hpp"

wdsp/bandpass.cpp

@@ -28,7 +28,7 @@ warren@wpratt.com
 #include "comm.hpp"
 #include "bandpass.hpp"
 #include "fir.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "RXA.hpp"
 #include "TXA.hpp"
 
@@ -40,12 +40,23 @@ namespace WDSP {
 *                                                                                                       *
 ********************************************************************************************************/
 
-BANDPASS* BANDPASS::create_bandpass (int run, int position, int size, int nc, int mp, float* in, float* out,
+BANDPASS* BANDPASS::create_bandpass (
-    float f_low, float f_high, int samplerate, int wintype, float gain)
+    int run,
+    int position,
+    int size,
+    int nc,
+    int mp,
+    float* in,
+    float* out,
+    float f_low,
+    float f_high,
+    int samplerate,
+    int wintype,
+    float gain
+)
 {
     // NOTE:  'nc' must be >= 'size'
     BANDPASS *a = new BANDPASS;
-    float* impulse;
     a->run = run;
     a->position = position;
     a->size = size;
@@ -58,7 +69,15 @@ BANDPASS* BANDPASS::create_bandpass (int run, int position, int size, int nc, in
     a->samplerate = samplerate;
     a->wintype = wintype;
     a->gain = gain;
-    impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+    float* impulse = FIR::fir_bandpass (
+        a->nc,
+        a->f_low,
+        a->f_high,
+        a->samplerate,
+        a->wintype,
+        1,
+        a->gain / (float)(2 * a->size)
+    );
     a->p = FIRCORE::create_fircore (a->size, a->in, a->out, a->nc, a->mp, impulse);
     delete[] impulse;
     return a;
@@ -92,9 +111,16 @@ void BANDPASS::setBuffers_bandpass (BANDPASS *a, float* in, float* out)
 
 void BANDPASS::setSamplerate_bandpass (BANDPASS *a, int rate)
 {
-    float* impulse;
     a->samplerate = rate;
-    impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+    float* impulse = FIR::fir_bandpass (
+        a->nc,
+        a->f_low,
+        a->f_high,
+        a->samplerate,
+        a->wintype,
+        1,
+        a->gain / (float)(2 * a->size)
+    );
     FIRCORE::setImpulse_fircore (a->p, impulse, 1);
     delete[] impulse;
 }
@@ -102,33 +128,54 @@ void BANDPASS::setSamplerate_bandpass (BANDPASS *a, int rate)
 void BANDPASS::setSize_bandpass (BANDPASS *a, int size)
 {
     // NOTE:  'size' must be <= 'nc'
-    float* impulse;
     a->size = size;
     FIRCORE::setSize_fircore (a->p, a->size);
     // recalc impulse because scale factor is a function of size
-    impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+    float* impulse = FIR::fir_bandpass (
+        a->nc,
+        a->f_low,
+        a->f_high,
+        a->samplerate,
+        a->wintype,
+        1,
+        a->gain / (float)(2 * a->size)
+    );
     FIRCORE::setImpulse_fircore (a->p, impulse, 1);
     delete[] (impulse);
 }
 
 void BANDPASS::setGain_bandpass (BANDPASS *a, float gain, int update)
 {
-    float* impulse;
     a->gain = gain;
-    impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+    float* impulse = FIR::fir_bandpass (
+        a->nc,
+        a->f_low,
+        a->f_high,
+        a->samplerate,
+        a->wintype,
+        1,
+        a->gain / (float)(2 * a->size)
+    );
     FIRCORE::setImpulse_fircore (a->p, impulse, update);
     delete[] (impulse);
 }
 
 void BANDPASS::CalcBandpassFilter (BANDPASS *a, float f_low, float f_high, float gain)
 {
-    float* impulse;
     if ((a->f_low != f_low) || (a->f_high != f_high) || (a->gain != gain))
     {
         a->f_low = f_low;
         a->f_high = f_high;
         a->gain = gain;
-        impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+        float* impulse = FIR::fir_bandpass (
+            a->nc,
+            a->f_low,
+            a->f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setImpulse_fircore (a->p, impulse, 1);
         delete[] (impulse);
     }
@@ -142,12 +189,19 @@ void BANDPASS::CalcBandpassFilter (BANDPASS *a, float f_low, float f_high, float
 
 void BANDPASS::SetBandpassFreqs (RXA& rxa, float f_low, float f_high)
 {
-    float* impulse;
     BANDPASS *a = rxa.bp1.p;
+
     if ((f_low != a->f_low) || (f_high != a->f_high))
     {
-        impulse = FIR::fir_bandpass (a->nc, f_low, f_high, a->samplerate,
+        float* impulse = FIR::fir_bandpass (
-            a->wintype, 1, a->gain / (float)(2 * a->size));
+            a->nc,
+            f_low,
+            f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setImpulse_fircore (a->p, impulse, 0);
         delete[] (impulse);
         rxa.csDSP.lock();
@@ -161,17 +215,26 @@ void BANDPASS::SetBandpassFreqs (RXA& rxa, float f_low, float f_high)
 void BANDPASS::SetBandpassNC (RXA& rxa, int nc)
 {
     // NOTE:  'nc' must be >= 'size'
-    float* impulse;
     BANDPASS *a;
     rxa.csDSP.lock();
     a = rxa.bp1.p;
+
     if (nc != a->nc)
     {
         a->nc = nc;
-        impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+        float* impulse = FIR::fir_bandpass (
+            a->nc,
+            a->f_low,
+            a->f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setNc_fircore (a->p, a->nc, impulse);
         delete[] (impulse);
     }
+
     rxa.csDSP.unlock();
 }
 
@@ -179,6 +242,7 @@ void BANDPASS::SetBandpassMP (RXA& rxa, int mp)
 {
     BANDPASS *a;
     a = rxa.bp1.p;
+
     if (mp != a->mp)
     {
         a->mp = mp;
@@ -229,33 +293,62 @@ void BANDPASS::SetBandpassMP (RXA& rxa, int mp)
 void BANDPASS::SetBandpassNC (TXA& txa, int nc)
 {
     // NOTE:  'nc' must be >= 'size'
-    float* impulse;
     BANDPASS *a;
     txa.csDSP.lock();
     a = txa.bp0.p;
+
     if (a->nc != nc)
     {
         a->nc = nc;
-        impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+        float* impulse = FIR::fir_bandpass (
+            a->nc,
+            a->f_low,
+            a->f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setNc_fircore (a->p, a->nc, impulse);
         delete[] (impulse);
     }
+
     a = txa.bp1.p;
+
     if (a->nc != nc)
     {
         a->nc = nc;
-        impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+        float* impulse = FIR::fir_bandpass (
+            a->nc,
+            a->f_low,
+            a->f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setNc_fircore (a->p, a->nc, impulse);
         delete[] (impulse);
     }
+
     a = txa.bp2.p;
+
     if (a->nc != nc)
     {
         a->nc = nc;
-        impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
+        float* impulse = FIR::fir_bandpass (
+            a->nc,
+            a->f_low,
+            a->f_high,
+            a->samplerate,
+            a->wintype,
+            1,
+            a->gain / (float)(2 * a->size)
+        );
         FIRCORE::setNc_fircore (a->p, a->nc, impulse);
         delete[] (impulse);
     }
+
     txa.csDSP.unlock();
 }
 
@@ -263,18 +356,23 @@ void BANDPASS::SetBandpassMP (TXA& txa, int mp)
 {
     BANDPASS *a;
     a = txa.bp0.p;
+
     if (mp != a->mp)
     {
         a->mp = mp;
         FIRCORE::setMp_fircore (a->p, a->mp);
     }
+
     a = txa.bp1.p;
+
     if (mp != a->mp)
     {
         a->mp = mp;
         FIRCORE::setMp_fircore (a->p, a->mp);
     }
+
     a = txa.bp2.p;
+
     if (mp != a->mp)
     {
         a->mp = mp;

wdsp/bpsnba.cpp

@@ -28,7 +28,7 @@ warren@wpratt.com
 #include "comm.hpp"
 #include "resample.hpp"
 #include "lmath.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "nbp.hpp"
 #include "amd.hpp"
 #include "anf.hpp"

wdsp/cfir.cpp

@@ -27,7 +27,7 @@ warren@wpratt.com
 #include "comm.hpp"
 #include "cfir.hpp"
 #include "fir.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "TXA.hpp"
 
 namespace WDSP {

wdsp/emph.cpp

@@ -28,7 +28,7 @@ warren@wpratt.com
 #include "comm.hpp"
 #include "emph.hpp"
 #include "fcurve.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "TXA.hpp"
 
 namespace WDSP {

wdsp/eq.cpp

@@ -27,7 +27,7 @@ warren@wpratt.com
 
 #include "comm.hpp"
 #include "eq.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "fir.hpp"
 #include "RXA.hpp"
 #include "TXA.hpp"

wdsp/fircore.cpp

@@ -0,0 +1,236 @@
+/*  firmin.c
+
+This file is part of a program that implements a Software-Defined Radio.
+
+Copyright (C) 2016 Warren Pratt, NR0V
+Copyright (C) 2024 Edouard Griffiths, F4EXB Adapted to SDRangel
+
+This program is free software; you can redistribute it and/or
+modify it under the terms of the GNU General Public License
+as published by the Free Software Foundation; either version 2
+of the License, or (at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with this program; if not, write to the Free Software
+Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
+
+The author can be reached by email at
+
+warren@wpratt.com
+
+*/
+
+#include "comm.hpp"
+#include "fir.hpp"
+#include "fircore.hpp"
+
+namespace WDSP {
+
+/********************************************************************************************************
+*                                                                                                       *
+*                                   Partitioned Overlap-Save Filter Kernel                              *
+*                                                                                                       *
+********************************************************************************************************/
+
+
+void FIRCORE::plan_fircore (FIRCORE *a)
+{
+    // must call for change in 'nc', 'size', 'out'
+    int i;
+    a->nfor = a->nc / a->size;
+    a->cset = 0;
+    a->buffidx = 0;
+    a->idxmask = a->nfor - 1;
+    a->fftin = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->fftout   = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
+    a->fmask    = new float**[2]; // (float ***) malloc0 (2 * sizeof (float **));
+    a->fmask[0] = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
+    a->fmask[1] = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
+    a->maskgen = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->pcfor = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
+    a->maskplan    = new fftwf_plan*[2]; // (fftwf_plan **) malloc0 (2 * sizeof (fftwf_plan *));
+    a->maskplan[0] = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
+    a->maskplan[1] = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
+    for (i = 0; i < a->nfor; i++)
+    {
+        a->fftout[i]   = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+        a->fmask[0][i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+        a->fmask[1][i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+        a->pcfor[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->fftin, (fftwf_complex *)a->fftout[i], FFTW_FORWARD, FFTW_PATIENT);
+        a->maskplan[0][i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[0][i], FFTW_FORWARD, FFTW_PATIENT);
+        a->maskplan[1][i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[1][i], FFTW_FORWARD, FFTW_PATIENT);
+    }
+    a->accum = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->crev = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->accum, (fftwf_complex *)a->out, FFTW_BACKWARD, FFTW_PATIENT);
+    a->masks_ready = 0;
+}
+
+void FIRCORE::calc_fircore (FIRCORE *a, int flip)
+{
+    // call for change in frequency, rate, wintype, gain
+    // must also call after a call to plan_firopt()
+    int i;
+    if (a->mp)
+        FIR::mp_imp (a->nc, a->impulse, a->imp, 16, 0);
+    else
+        memcpy (a->imp, a->impulse, a->nc * sizeof (wcomplex));
+    for (i = 0; i < a->nfor; i++)
+    {
+        // I right-justified the impulse response => take output from left side of output buff, discard right side
+        // Be careful about flipping an asymmetrical impulse response.
+        memcpy (&(a->maskgen[2 * a->size]), &(a->imp[2 * a->size * i]), a->size * sizeof(wcomplex));
+        fftwf_execute (a->maskplan[1 - a->cset][i]);
+    }
+    a->masks_ready = 1;
+    if (flip)
+    {
+        a->update.lock();
+        a->cset = 1 - a->cset;
+        a->update.unlock();
+        a->masks_ready = 0;
+    }
+}
+
+FIRCORE* FIRCORE::create_fircore (int size, float* in, float* out, int nc, int mp, float* impulse)
+{
+    FIRCORE *a = new FIRCORE;
+    a->size = size;
+    a->in = in;
+    a->out = out;
+    a->nc = nc;
+    a->mp = mp;
+    // InitializeCriticalSectionAndSpinCount (&a->update, 2500);
+    plan_fircore (a);
+    a->impulse = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
+    a->imp     = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
+    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
+    calc_fircore (a, 1);
+    return a;
+}
+
+void FIRCORE::deplan_fircore (FIRCORE *a)
+{
+    int i;
+    fftwf_destroy_plan (a->crev);
+    delete[] (a->accum);
+    for (i = 0; i < a->nfor; i++)
+    {
+        delete[] (a->fftout[i]);
+        delete[] (a->fmask[0][i]);
+        delete[] (a->fmask[1][i]);
+        fftwf_destroy_plan (a->pcfor[i]);
+        fftwf_destroy_plan (a->maskplan[0][i]);
+        fftwf_destroy_plan (a->maskplan[1][i]);
+    }
+    delete[] (a->maskplan[0]);
+    delete[] (a->maskplan[1]);
+    delete[] (a->maskplan);
+    delete[] (a->pcfor);
+    delete[] (a->maskgen);
+    delete[] (a->fmask[0]);
+    delete[] (a->fmask[1]);
+    delete[] (a->fmask);
+    delete[] (a->fftout);
+    delete[] (a->fftin);
+}
+
+void FIRCORE::destroy_fircore (FIRCORE *a)
+{
+    deplan_fircore (a);
+    delete[] (a->imp);
+    delete[] (a->impulse);
+    delete (a);
+}
+
+void FIRCORE::flush_fircore (FIRCORE *a)
+{
+    int i;
+    memset (a->fftin, 0, 2 * a->size * sizeof (wcomplex));
+    for (i = 0; i < a->nfor; i++)
+        memset (a->fftout[i], 0, 2 * a->size * sizeof (wcomplex));
+    a->buffidx = 0;
+}
+
+void FIRCORE::xfircore (FIRCORE *a)
+{
+    int i, j, k;
+    memcpy (&(a->fftin[2 * a->size]), a->in, a->size * sizeof (wcomplex));
+    fftwf_execute (a->pcfor[a->buffidx]);
+    k = a->buffidx;
+    memset (a->accum, 0, 2 * a->size * sizeof (wcomplex));
+    a->update.lock();
+    for (j = 0; j < a->nfor; j++)
+    {
+        for (i = 0; i < 2 * a->size; i++)
+        {
+            a->accum[2 * i + 0] += a->fftout[k][2 * i + 0] * a->fmask[a->cset][j][2 * i + 0] - a->fftout[k][2 * i + 1] * a->fmask[a->cset][j][2 * i + 1];
+            a->accum[2 * i + 1] += a->fftout[k][2 * i + 0] * a->fmask[a->cset][j][2 * i + 1] + a->fftout[k][2 * i + 1] * a->fmask[a->cset][j][2 * i + 0];
+        }
+        k = (k + a->idxmask) & a->idxmask;
+    }
+    a->update.unlock();
+    a->buffidx = (a->buffidx + 1) & a->idxmask;
+    fftwf_execute (a->crev);
+    memcpy (a->fftin, &(a->fftin[2 * a->size]), a->size * sizeof(wcomplex));
+}
+
+void FIRCORE::setBuffers_fircore (FIRCORE *a, float* in, float* out)
+{
+    a->in = in;
+    a->out = out;
+    deplan_fircore (a);
+    plan_fircore (a);
+    calc_fircore (a, 1);
+}
+
+void FIRCORE::setSize_fircore (FIRCORE *a, int size)
+{
+    a->size = size;
+    deplan_fircore (a);
+    plan_fircore (a);
+    calc_fircore (a, 1);
+}
+
+void FIRCORE::setImpulse_fircore (FIRCORE *a, float* impulse, int update)
+{
+    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
+    calc_fircore (a, update);
+}
+
+void FIRCORE::setNc_fircore (FIRCORE *a, int nc, float* impulse)
+{
+    // because of FFT planning, this will probably cause a glitch in audio if done during dataflow
+    deplan_fircore (a);
+    delete[] (a->impulse);
+    delete[] (a->imp);
+    a->nc = nc;
+    plan_fircore (a);
+    a->imp     = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
+    a->impulse = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
+    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
+    calc_fircore (a, 1);
+}
+
+void FIRCORE::setMp_fircore (FIRCORE *a, int mp)
+{
+    a->mp = mp;
+    calc_fircore (a, 1);
+}
+
+void FIRCORE::setUpdate_fircore (FIRCORE *a)
+{
+    if (a->masks_ready)
+    {
+        a->update.lock();
+        a->cset = 1 - a->cset;
+        a->update.unlock();
+        a->masks_ready = 0;
+    }
+}
+
+} // namespace WDSP

wdsp/fircore.hpp

@@ -0,0 +1,88 @@
+/*  firmin.h
+
+This file is part of a program that implements a Software-Defined Radio.
+
+Copyright (C) 2016 Warren Pratt, NR0V
+Copyright (C) 2024 Edouard Griffiths, F4EXB Adapted to SDRangel
+
+This program is free software; you can redistribute it and/or
+modify it under the terms of the GNU General Public License
+as published by the Free Software Foundation; either version 2
+of the License, or (at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with this program; if not, write to the Free Software
+Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
+
+The author can be reached by email at
+
+warren@wpratt.com
+
+*/
+
+/********************************************************************************************************
+*                                                                                                       *
+*                                   Partitioned Overlap-Save Filter Kernel                              *
+*                                                                                                       *
+********************************************************************************************************/
+
+#ifndef wdsp_fircore_h
+#define wdsp_fircore_h
+
+#include <QRecursiveMutex>
+#include "fftw3.h"
+#include "export.h"
+
+namespace WDSP {
+
+class WDSP_API FIRCORE
+{
+public:
+    int size;               // input/output buffer size, power of two
+    float* in;             // input buffer
+    float* out;            // output buffer, can be same as input
+    int nc;                 // number of filter coefficients, power of two, >= size
+    float* impulse;        // impulse response of filter
+    float* imp;
+    int nfor;               // number of buffers in delay line
+    float* fftin;          // fft input buffer
+    float*** fmask;        // frequency domain masks
+    float** fftout;        // fftout delay line
+    float* accum;          // frequency domain accumulator
+    int buffidx;            // fft out buffer index
+    int idxmask;            // mask for index computations
+    float* maskgen;        // input for mask generation FFT
+    fftwf_plan* pcfor;       // array of forward FFT plans
+    fftwf_plan crev;         // reverse fft plan
+    fftwf_plan** maskplan;   // plans for frequency domain masks
+    QRecursiveMutex update;
+    int cset;
+    int mp;
+    int masks_ready;
+
+    static FIRCORE* create_fircore (int size, float* in, float* out,
+        int nc, int mp, float* impulse);
+    static void xfircore (FIRCORE *a);
+    static void destroy_fircore (FIRCORE *a);
+    static void flush_fircore (FIRCORE *a);
+    static void setBuffers_fircore (FIRCORE *a, float* in, float* out);
+    static void setSize_fircore (FIRCORE *a, int size);
+    static void setImpulse_fircore (FIRCORE *a, float* impulse, int update);
+    static void setNc_fircore (FIRCORE *a, int nc, float* impulse);
+    static void setMp_fircore (FIRCORE *a, int mp);
+    static void setUpdate_fircore (FIRCORE *a);
+
+private:
+    static void plan_fircore (FIRCORE *a);
+    static void calc_fircore (FIRCORE *a, int flip);
+    static void deplan_fircore (FIRCORE *a);
+};
+
+} // namespace WDSP
+
+#endif

wdsp/firmin.cpp

@@ -127,364 +127,4 @@ void FIRMIN::setFreqs_firmin (FIRMIN *a, float f_low, float f_high)
     calc_firmin (a);
 }
 
-/********************************************************************************************************
-*                                                                                                       *
-*                               Standalone Partitioned Overlap-Save Bandpass                            *
-*                                                                                                       *
-********************************************************************************************************/
-
-void FIROPT::plan_firopt (FIROPT *a)
-{
-    // must call for change in 'nc', 'size', 'out'
-    int i;
-    a->nfor = a->nc / a->size;
-    a->buffidx = 0;
-    a->idxmask = a->nfor - 1;
-    a->fftin = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->fftout = new float*[a->nfor]; //  (float **) malloc0 (a->nfor * sizeof (float *));
-    a->fmask = new float*[a->nfor]; //  (float **) malloc0 (a->nfor * sizeof (float *));
-    a->maskgen = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->pcfor = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
-    a->maskplan = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
-    for (i = 0; i < a->nfor; i++)
-    {
-        a->fftout[i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-        a->fmask[i] = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
-        a->pcfor[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->fftin, (fftwf_complex *)a->fftout[i], FFTW_FORWARD, FFTW_PATIENT);
-        a->maskplan[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[i], FFTW_FORWARD, FFTW_PATIENT);
-    }
-    a->accum = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->crev = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->accum, (fftwf_complex *)a->out, FFTW_BACKWARD, FFTW_PATIENT);
-}
-
-void FIROPT::calc_firopt (FIROPT *a)
-{
-    // call for change in frequency, rate, wintype, gain
-    // must also call after a call to plan_firopt()
-    int i;
-    float* impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain);
-    a->buffidx = 0;
-    for (i = 0; i < a->nfor; i++)
-    {
-        // I right-justified the impulse response => take output from left side of output buff, discard right side
-        // Be careful about flipping an asymmetrical impulse response.
-        memcpy (&(a->maskgen[2 * a->size]), &(impulse[2 * a->size * i]), a->size * sizeof(wcomplex));
-        fftwf_execute (a->maskplan[i]);
-    }
-    delete[] (impulse);
-}
-
-FIROPT* FIROPT::create_firopt (int run, int position, int size, float* in, float* out,
-    int nc, float f_low, float f_high, int samplerate, int wintype, float gain)
-{
-    FIROPT *a = new FIROPT;
-    a->run = run;
-    a->position = position;
-    a->size = size;
-    a->in = in;
-    a->out = out;
-    a->nc = nc;
-    a->f_low = f_low;
-    a->f_high = f_high;
-    a->samplerate = samplerate;
-    a->wintype = wintype;
-    a->gain = gain;
-    plan_firopt (a);
-    calc_firopt (a);
-    return a;
-}
-
-void FIROPT::deplan_firopt (FIROPT *a)
-{
-    int i;
-    fftwf_destroy_plan (a->crev);
-    delete[] (a->accum);
-    for (i = 0; i < a->nfor; i++)
-    {
-        delete[] (a->fftout[i]);
-        delete[] (a->fmask[i]);
-        fftwf_destroy_plan (a->pcfor[i]);
-        fftwf_destroy_plan (a->maskplan[i]);
-    }
-    delete[] (a->maskplan);
-    delete[] (a->pcfor);
-    delete[] (a->maskgen);
-    delete[] (a->fmask);
-    delete[] (a->fftout);
-    delete[] (a->fftin);
-}
-
-void FIROPT::destroy_firopt (FIROPT *a)
-{
-    deplan_firopt (a);
-    delete (a);
-}
-
-void FIROPT::flush_firopt (FIROPT *a)
-{
-    int i;
-    memset (a->fftin, 0, 2 * a->size * sizeof (wcomplex));
-    for (i = 0; i < a->nfor; i++)
-        memset (a->fftout[i], 0, 2 * a->size * sizeof (wcomplex));
-    a->buffidx = 0;
-}
-
-void FIROPT::xfiropt (FIROPT *a, int pos)
-{
-    if (a->run && (a->position == pos))
-    {
-        int i, j, k;
-        memcpy (&(a->fftin[2 * a->size]), a->in, a->size * sizeof (wcomplex));
-        fftwf_execute (a->pcfor[a->buffidx]);
-        k = a->buffidx;
-        memset (a->accum, 0, 2 * a->size * sizeof (wcomplex));
-        for (j = 0; j < a->nfor; j++)
-        {
-            for (i = 0; i < 2 * a->size; i++)
-            {
-                a->accum[2 * i + 0] += a->fftout[k][2 * i + 0] * a->fmask[j][2 * i + 0] - a->fftout[k][2 * i + 1] * a->fmask[j][2 * i + 1];
-                a->accum[2 * i + 1] += a->fftout[k][2 * i + 0] * a->fmask[j][2 * i + 1] + a->fftout[k][2 * i + 1] * a->fmask[j][2 * i + 0];
-            }
-            k = (k + a->idxmask) & a->idxmask;
-        }
-        a->buffidx = (a->buffidx + 1) & a->idxmask;
-        fftwf_execute (a->crev);
-        memcpy (a->fftin, &(a->fftin[2 * a->size]), a->size * sizeof(wcomplex));
-    }
-    else if (a->in != a->out)
-        memcpy (a->out, a->in, a->size * sizeof (wcomplex));
-}
-
-void FIROPT::setBuffers_firopt (FIROPT *a, float* in, float* out)
-{
-    a->in = in;
-    a->out = out;
-    deplan_firopt (a);
-    plan_firopt (a);
-    calc_firopt (a);
-}
-
-void FIROPT::setSamplerate_firopt (FIROPT *a, int rate)
-{
-    a->samplerate = rate;
-    calc_firopt (a);
-}
-
-void FIROPT::setSize_firopt (FIROPT *a, int size)
-{
-    a->size = size;
-    deplan_firopt (a);
-    plan_firopt (a);
-    calc_firopt (a);
-}
-
-void FIROPT::setFreqs_firopt (FIROPT *a, float f_low, float f_high)
-{
-    a->f_low = f_low;
-    a->f_high = f_high;
-    calc_firopt (a);
-}
-
-/********************************************************************************************************
-*                                                                                                       *
-*                                   Partitioned Overlap-Save Filter Kernel                              *
-*                                                                                                       *
-********************************************************************************************************/
-
-
-void FIRCORE::plan_fircore (FIRCORE *a)
-{
-    // must call for change in 'nc', 'size', 'out'
-    int i;
-    a->nfor = a->nc / a->size;
-    a->cset = 0;
-    a->buffidx = 0;
-    a->idxmask = a->nfor - 1;
-    a->fftin = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->fftout   = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
-    a->fmask    = new float**[2]; // (float ***) malloc0 (2 * sizeof (float **));
-    a->fmask[0] = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
-    a->fmask[1] = new float*[a->nfor]; // (float **) malloc0 (a->nfor * sizeof (float *));
-    a->maskgen = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->pcfor = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
-    a->maskplan    = new fftwf_plan*[2]; // (fftwf_plan **) malloc0 (2 * sizeof (fftwf_plan *));
-    a->maskplan[0] = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
-    a->maskplan[1] = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
-    for (i = 0; i < a->nfor; i++)
-    {
-        a->fftout[i]   = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-        a->fmask[0][i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-        a->fmask[1][i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-        a->pcfor[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->fftin, (fftwf_complex *)a->fftout[i], FFTW_FORWARD, FFTW_PATIENT);
-        a->maskplan[0][i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[0][i], FFTW_FORWARD, FFTW_PATIENT);
-        a->maskplan[1][i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[1][i], FFTW_FORWARD, FFTW_PATIENT);
-    }
-    a->accum = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
-    a->crev = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->accum, (fftwf_complex *)a->out, FFTW_BACKWARD, FFTW_PATIENT);
-    a->masks_ready = 0;
-}
-
-void FIRCORE::calc_fircore (FIRCORE *a, int flip)
-{
-    // call for change in frequency, rate, wintype, gain
-    // must also call after a call to plan_firopt()
-    int i;
-    if (a->mp)
-        FIR::mp_imp (a->nc, a->impulse, a->imp, 16, 0);
-    else
-        memcpy (a->imp, a->impulse, a->nc * sizeof (wcomplex));
-    for (i = 0; i < a->nfor; i++)
-    {
-        // I right-justified the impulse response => take output from left side of output buff, discard right side
-        // Be careful about flipping an asymmetrical impulse response.
-        memcpy (&(a->maskgen[2 * a->size]), &(a->imp[2 * a->size * i]), a->size * sizeof(wcomplex));
-        fftwf_execute (a->maskplan[1 - a->cset][i]);
-    }
-    a->masks_ready = 1;
-    if (flip)
-    {
-        a->update.lock();
-        a->cset = 1 - a->cset;
-        a->update.unlock();
-        a->masks_ready = 0;
-    }
-}
-
-FIRCORE* FIRCORE::create_fircore (int size, float* in, float* out, int nc, int mp, float* impulse)
-{
-    FIRCORE *a = new FIRCORE;
-    a->size = size;
-    a->in = in;
-    a->out = out;
-    a->nc = nc;
-    a->mp = mp;
-    // InitializeCriticalSectionAndSpinCount (&a->update, 2500);
-    plan_fircore (a);
-    a->impulse = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
-    a->imp     = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
-    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
-    calc_fircore (a, 1);
-    return a;
-}
-
-void FIRCORE::deplan_fircore (FIRCORE *a)
-{
-    int i;
-    fftwf_destroy_plan (a->crev);
-    delete[] (a->accum);
-    for (i = 0; i < a->nfor; i++)
-    {
-        delete[] (a->fftout[i]);
-        delete[] (a->fmask[0][i]);
-        delete[] (a->fmask[1][i]);
-        fftwf_destroy_plan (a->pcfor[i]);
-        fftwf_destroy_plan (a->maskplan[0][i]);
-        fftwf_destroy_plan (a->maskplan[1][i]);
-    }
-    delete[] (a->maskplan[0]);
-    delete[] (a->maskplan[1]);
-    delete[] (a->maskplan);
-    delete[] (a->pcfor);
-    delete[] (a->maskgen);
-    delete[] (a->fmask[0]);
-    delete[] (a->fmask[1]);
-    delete[] (a->fmask);
-    delete[] (a->fftout);
-    delete[] (a->fftin);
-}
-
-void FIRCORE::destroy_fircore (FIRCORE *a)
-{
-    deplan_fircore (a);
-    delete[] (a->imp);
-    delete[] (a->impulse);
-    delete (a);
-}
-
-void FIRCORE::flush_fircore (FIRCORE *a)
-{
-    int i;
-    memset (a->fftin, 0, 2 * a->size * sizeof (wcomplex));
-    for (i = 0; i < a->nfor; i++)
-        memset (a->fftout[i], 0, 2 * a->size * sizeof (wcomplex));
-    a->buffidx = 0;
-}
-
-void FIRCORE::xfircore (FIRCORE *a)
-{
-    int i, j, k;
-    memcpy (&(a->fftin[2 * a->size]), a->in, a->size * sizeof (wcomplex));
-    fftwf_execute (a->pcfor[a->buffidx]);
-    k = a->buffidx;
-    memset (a->accum, 0, 2 * a->size * sizeof (wcomplex));
-    a->update.lock();
-    for (j = 0; j < a->nfor; j++)
-    {
-        for (i = 0; i < 2 * a->size; i++)
-        {
-            a->accum[2 * i + 0] += a->fftout[k][2 * i + 0] * a->fmask[a->cset][j][2 * i + 0] - a->fftout[k][2 * i + 1] * a->fmask[a->cset][j][2 * i + 1];
-            a->accum[2 * i + 1] += a->fftout[k][2 * i + 0] * a->fmask[a->cset][j][2 * i + 1] + a->fftout[k][2 * i + 1] * a->fmask[a->cset][j][2 * i + 0];
-        }
-        k = (k + a->idxmask) & a->idxmask;
-    }
-    a->update.unlock();
-    a->buffidx = (a->buffidx + 1) & a->idxmask;
-    fftwf_execute (a->crev);
-    memcpy (a->fftin, &(a->fftin[2 * a->size]), a->size * sizeof(wcomplex));
-}
-
-void FIRCORE::setBuffers_fircore (FIRCORE *a, float* in, float* out)
-{
-    a->in = in;
-    a->out = out;
-    deplan_fircore (a);
-    plan_fircore (a);
-    calc_fircore (a, 1);
-}
-
-void FIRCORE::setSize_fircore (FIRCORE *a, int size)
-{
-    a->size = size;
-    deplan_fircore (a);
-    plan_fircore (a);
-    calc_fircore (a, 1);
-}
-
-void FIRCORE::setImpulse_fircore (FIRCORE *a, float* impulse, int update)
-{
-    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
-    calc_fircore (a, update);
-}
-
-void FIRCORE::setNc_fircore (FIRCORE *a, int nc, float* impulse)
-{
-    // because of FFT planning, this will probably cause a glitch in audio if done during dataflow
-    deplan_fircore (a);
-    delete[] (a->impulse);
-    delete[] (a->imp);
-    a->nc = nc;
-    plan_fircore (a);
-    a->imp     = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
-    a->impulse = new float[a->nc * 2]; // (float *) malloc0 (a->nc * sizeof (complex));
-    memcpy (a->impulse, impulse, a->nc * sizeof (wcomplex));
-    calc_fircore (a, 1);
-}
-
-void FIRCORE::setMp_fircore (FIRCORE *a, int mp)
-{
-    a->mp = mp;
-    calc_fircore (a, 1);
-}
-
-void FIRCORE::setUpdate_fircore (FIRCORE *a)
-{
-    if (a->masks_ready)
-    {
-        a->update.lock();
-        a->cset = 1 - a->cset;
-        a->update.unlock();
-        a->masks_ready = 0;
-    }
-}
-
 } // namespace WDSP

wdsp/firmin.hpp

@@ -76,125 +76,3 @@ private:
 } // namespace WDSP
 
 #endif
-
-
-/********************************************************************************************************
-*                                                                                                       *
-*                               Standalone Partitioned Overlap-Save Bandpass                            *
-*                                                                                                       *
-********************************************************************************************************/
-
-#ifndef wdsp_firopt_h
-#define wdsp_firopt_h
-
-#include "fftw3.h"
-#include "export.h"
-
-namespace WDSP {
-
-class WDSP_API FIROPT
-{
-    int run;                // run control
-    int position;           // position at which to execute
-    int size;               // input/output buffer size, power of two
-    float* in;             // input buffer
-    float* out;            // output buffer, can be same as input
-    int nc;                 // number of filter coefficients, power of two, >= size
-    float f_low;           // low cutoff frequency
-    float f_high;          // high cutoff frequency
-    float samplerate;      // sample rate
-    int wintype;            // filter window type
-    float gain;            // filter gain
-    int nfor;               // number of buffers in delay line
-    float* fftin;          // fft input buffer
-    float** fmask;         // frequency domain masks
-    float** fftout;        // fftout delay line
-    float* accum;          // frequency domain accumulator
-    int buffidx;            // fft out buffer index
-    int idxmask;            // mask for index computations
-    float* maskgen;        // input for mask generation FFT
-    fftwf_plan* pcfor;       // array of forward FFT plans
-    fftwf_plan crev;         // reverse fft plan
-    fftwf_plan* maskplan;    // plans for frequency domain masks
-
-    static FIROPT* create_firopt (int run, int position, int size, float* in, float* out,
-        int nc, float f_low, float f_high, int samplerate, int wintype, float gain);
-    static void xfiropt (FIROPT *a, int pos);
-    static void destroy_firopt (FIROPT *a);
-    static void flush_firopt (FIROPT *a);
-    static void setBuffers_firopt (FIROPT *a, float* in, float* out);
-    static void setSamplerate_firopt (FIROPT *a, int rate);
-    static void setSize_firopt (FIROPT *a, int size);
-    static void setFreqs_firopt (FIROPT *a, float f_low, float f_high);
-
-private:
-    static void plan_firopt (FIROPT *a);
-    static void calc_firopt (FIROPT *a);
-    static void deplan_firopt (FIROPT *a);
-};
-
-} // namespace WDSP
-
-#endif
-
-
-/********************************************************************************************************
-*                                                                                                       *
-*                                   Partitioned Overlap-Save Filter Kernel                              *
-*                                                                                                       *
-********************************************************************************************************/
-
-#ifndef wdsp_fircore_h
-#define wdsp_fircore_h
-
-#include <QRecursiveMutex>
-#include "export.h"
-
-namespace WDSP {
-
-class WDSP_API FIRCORE
-{
-public:
-    int size;               // input/output buffer size, power of two
-    float* in;             // input buffer
-    float* out;            // output buffer, can be same as input
-    int nc;                 // number of filter coefficients, power of two, >= size
-    float* impulse;        // impulse response of filter
-    float* imp;
-    int nfor;               // number of buffers in delay line
-    float* fftin;          // fft input buffer
-    float*** fmask;        // frequency domain masks
-    float** fftout;        // fftout delay line
-    float* accum;          // frequency domain accumulator
-    int buffidx;            // fft out buffer index
-    int idxmask;            // mask for index computations
-    float* maskgen;        // input for mask generation FFT
-    fftwf_plan* pcfor;       // array of forward FFT plans
-    fftwf_plan crev;         // reverse fft plan
-    fftwf_plan** maskplan;   // plans for frequency domain masks
-    QRecursiveMutex update;
-    int cset;
-    int mp;
-    int masks_ready;
-
-    static FIRCORE* create_fircore (int size, float* in, float* out,
-        int nc, int mp, float* impulse);
-    static void xfircore (FIRCORE *a);
-    static void destroy_fircore (FIRCORE *a);
-    static void flush_fircore (FIRCORE *a);
-    static void setBuffers_fircore (FIRCORE *a, float* in, float* out);
-    static void setSize_fircore (FIRCORE *a, int size);
-    static void setImpulse_fircore (FIRCORE *a, float* impulse, int update);
-    static void setNc_fircore (FIRCORE *a, int nc, float* impulse);
-    static void setMp_fircore (FIRCORE *a, int mp);
-    static void setUpdate_fircore (FIRCORE *a);
-
-private:
-    static void plan_fircore (FIRCORE *a);
-    static void calc_fircore (FIRCORE *a, int flip);
-    static void deplan_fircore (FIRCORE *a);
-};
-
-} // namespace WDSP
-
-#endif

wdsp/firopt.cpp

@@ -0,0 +1,193 @@
+/*  firmin.c
+
+This file is part of a program that implements a Software-Defined Radio.
+
+Copyright (C) 2016 Warren Pratt, NR0V
+Copyright (C) 2024 Edouard Griffiths, F4EXB Adapted to SDRangel
+
+This program is free software; you can redistribute it and/or
+modify it under the terms of the GNU General Public License
+as published by the Free Software Foundation; either version 2
+of the License, or (at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with this program; if not, write to the Free Software
+Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
+
+The author can be reached by email at
+
+warren@wpratt.com
+
+*/
+
+#include "comm.hpp"
+#include "fir.hpp"
+#include "firopt.hpp"
+
+namespace WDSP {
+
+/********************************************************************************************************
+*                                                                                                       *
+*                               Standalone Partitioned Overlap-Save Bandpass                            *
+*                                                                                                       *
+********************************************************************************************************/
+
+void FIROPT::plan_firopt (FIROPT *a)
+{
+    // must call for change in 'nc', 'size', 'out'
+    int i;
+    a->nfor = a->nc / a->size;
+    a->buffidx = 0;
+    a->idxmask = a->nfor - 1;
+    a->fftin = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->fftout = new float*[a->nfor]; //  (float **) malloc0 (a->nfor * sizeof (float *));
+    a->fmask = new float*[a->nfor]; //  (float **) malloc0 (a->nfor * sizeof (float *));
+    a->maskgen = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->pcfor = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
+    a->maskplan = new fftwf_plan[a->nfor]; // (fftwf_plan *) malloc0 (a->nfor * sizeof (fftwf_plan));
+    for (i = 0; i < a->nfor; i++)
+    {
+        a->fftout[i] = new float[2 * a->size * 2]; // (float *) malloc0 (2 * a->size * sizeof (complex));
+        a->fmask[i] = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
+        a->pcfor[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->fftin, (fftwf_complex *)a->fftout[i], FFTW_FORWARD, FFTW_PATIENT);
+        a->maskplan[i] = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->maskgen, (fftwf_complex *)a->fmask[i], FFTW_FORWARD, FFTW_PATIENT);
+    }
+    a->accum = new float[2 * a->size * 2]; //  (float *) malloc0 (2 * a->size * sizeof (complex));
+    a->crev = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->accum, (fftwf_complex *)a->out, FFTW_BACKWARD, FFTW_PATIENT);
+}
+
+void FIROPT::calc_firopt (FIROPT *a)
+{
+    // call for change in frequency, rate, wintype, gain
+    // must also call after a call to plan_firopt()
+    int i;
+    float* impulse = FIR::fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain);
+    a->buffidx = 0;
+    for (i = 0; i < a->nfor; i++)
+    {
+        // I right-justified the impulse response => take output from left side of output buff, discard right side
+        // Be careful about flipping an asymmetrical impulse response.
+        memcpy (&(a->maskgen[2 * a->size]), &(impulse[2 * a->size * i]), a->size * sizeof(wcomplex));
+        fftwf_execute (a->maskplan[i]);
+    }
+    delete[] (impulse);
+}
+
+FIROPT* FIROPT::create_firopt (int run, int position, int size, float* in, float* out,
+    int nc, float f_low, float f_high, int samplerate, int wintype, float gain)
+{
+    FIROPT *a = new FIROPT;
+    a->run = run;
+    a->position = position;
+    a->size = size;
+    a->in = in;
+    a->out = out;
+    a->nc = nc;
+    a->f_low = f_low;
+    a->f_high = f_high;
+    a->samplerate = samplerate;
+    a->wintype = wintype;
+    a->gain = gain;
+    plan_firopt (a);
+    calc_firopt (a);
+    return a;
+}
+
+void FIROPT::deplan_firopt (FIROPT *a)
+{
+    int i;
+    fftwf_destroy_plan (a->crev);
+    delete[] (a->accum);
+    for (i = 0; i < a->nfor; i++)
+    {
+        delete[] (a->fftout[i]);
+        delete[] (a->fmask[i]);
+        fftwf_destroy_plan (a->pcfor[i]);
+        fftwf_destroy_plan (a->maskplan[i]);
+    }
+    delete[] (a->maskplan);
+    delete[] (a->pcfor);
+    delete[] (a->maskgen);
+    delete[] (a->fmask);
+    delete[] (a->fftout);
+    delete[] (a->fftin);
+}
+
+void FIROPT::destroy_firopt (FIROPT *a)
+{
+    deplan_firopt (a);
+    delete (a);
+}
+
+void FIROPT::flush_firopt (FIROPT *a)
+{
+    int i;
+    memset (a->fftin, 0, 2 * a->size * sizeof (wcomplex));
+    for (i = 0; i < a->nfor; i++)
+        memset (a->fftout[i], 0, 2 * a->size * sizeof (wcomplex));
+    a->buffidx = 0;
+}
+
+void FIROPT::xfiropt (FIROPT *a, int pos)
+{
+    if (a->run && (a->position == pos))
+    {
+        int i, j, k;
+        memcpy (&(a->fftin[2 * a->size]), a->in, a->size * sizeof (wcomplex));
+        fftwf_execute (a->pcfor[a->buffidx]);
+        k = a->buffidx;
+        memset (a->accum, 0, 2 * a->size * sizeof (wcomplex));
+        for (j = 0; j < a->nfor; j++)
+        {
+            for (i = 0; i < 2 * a->size; i++)
+            {
+                a->accum[2 * i + 0] += a->fftout[k][2 * i + 0] * a->fmask[j][2 * i + 0] - a->fftout[k][2 * i + 1] * a->fmask[j][2 * i + 1];
+                a->accum[2 * i + 1] += a->fftout[k][2 * i + 0] * a->fmask[j][2 * i + 1] + a->fftout[k][2 * i + 1] * a->fmask[j][2 * i + 0];
+            }
+            k = (k + a->idxmask) & a->idxmask;
+        }
+        a->buffidx = (a->buffidx + 1) & a->idxmask;
+        fftwf_execute (a->crev);
+        memcpy (a->fftin, &(a->fftin[2 * a->size]), a->size * sizeof(wcomplex));
+    }
+    else if (a->in != a->out)
+        memcpy (a->out, a->in, a->size * sizeof (wcomplex));
+}
+
+void FIROPT::setBuffers_firopt (FIROPT *a, float* in, float* out)
+{
+    a->in = in;
+    a->out = out;
+    deplan_firopt (a);
+    plan_firopt (a);
+    calc_firopt (a);
+}
+
+void FIROPT::setSamplerate_firopt (FIROPT *a, int rate)
+{
+    a->samplerate = rate;
+    calc_firopt (a);
+}
+
+void FIROPT::setSize_firopt (FIROPT *a, int size)
+{
+    a->size = size;
+    deplan_firopt (a);
+    plan_firopt (a);
+    calc_firopt (a);
+}
+
+void FIROPT::setFreqs_firopt (FIROPT *a, float f_low, float f_high)
+{
+    a->f_low = f_low;
+    a->f_high = f_high;
+    calc_firopt (a);
+}
+
+
+} // namespace WDSP

wdsp/firopt.hpp

@@ -0,0 +1,86 @@
+/*  firmin.h
+
+This file is part of a program that implements a Software-Defined Radio.
+
+Copyright (C) 2016 Warren Pratt, NR0V
+Copyright (C) 2024 Edouard Griffiths, F4EXB Adapted to SDRangel
+
+This program is free software; you can redistribute it and/or
+modify it under the terms of the GNU General Public License
+as published by the Free Software Foundation; either version 2
+of the License, or (at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with this program; if not, write to the Free Software
+Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
+
+The author can be reached by email at
+
+warren@wpratt.com
+
+*/
+
+/********************************************************************************************************
+*                                                                                                       *
+*                               Standalone Partitioned Overlap-Save Bandpass                            *
+*                                                                                                       *
+********************************************************************************************************/
+
+#ifndef wdsp_firopt_h
+#define wdsp_firopt_h
+
+#include "fftw3.h"
+#include "export.h"
+
+namespace WDSP {
+
+class WDSP_API FIROPT
+{
+    int run;                // run control
+    int position;           // position at which to execute
+    int size;               // input/output buffer size, power of two
+    float* in;             // input buffer
+    float* out;            // output buffer, can be same as input
+    int nc;                 // number of filter coefficients, power of two, >= size
+    float f_low;           // low cutoff frequency
+    float f_high;          // high cutoff frequency
+    float samplerate;      // sample rate
+    int wintype;            // filter window type
+    float gain;            // filter gain
+    int nfor;               // number of buffers in delay line
+    float* fftin;          // fft input buffer
+    float** fmask;         // frequency domain masks
+    float** fftout;        // fftout delay line
+    float* accum;          // frequency domain accumulator
+    int buffidx;            // fft out buffer index
+    int idxmask;            // mask for index computations
+    float* maskgen;        // input for mask generation FFT
+    fftwf_plan* pcfor;       // array of forward FFT plans
+    fftwf_plan crev;         // reverse fft plan
+    fftwf_plan* maskplan;    // plans for frequency domain masks
+
+    static FIROPT* create_firopt (int run, int position, int size, float* in, float* out,
+        int nc, float f_low, float f_high, int samplerate, int wintype, float gain);
+    static void xfiropt (FIROPT *a, int pos);
+    static void destroy_firopt (FIROPT *a);
+    static void flush_firopt (FIROPT *a);
+    static void setBuffers_firopt (FIROPT *a, float* in, float* out);
+    static void setSamplerate_firopt (FIROPT *a, int rate);
+    static void setSize_firopt (FIROPT *a, int size);
+    static void setFreqs_firopt (FIROPT *a, float f_low, float f_high);
+
+private:
+    static void plan_firopt (FIROPT *a);
+    static void calc_firopt (FIROPT *a);
+    static void deplan_firopt (FIROPT *a);
+};
+
+} // namespace WDSP
+
+#endif
+

wdsp/fmd.cpp

@@ -27,7 +27,7 @@ warren@wpratt.com
 
 #include "comm.hpp"
 #include "iir.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "fcurve.hpp"
 #include "fir.hpp"
 #include "wcpAGC.hpp"

wdsp/fmmod.cpp

@@ -26,7 +26,7 @@ warren@wpratt.com
 */
 
 #include "comm.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "fir.hpp"
 #include "fmmod.hpp"
 #include "TXA.hpp"

wdsp/fmsq.cpp

@@ -26,7 +26,7 @@ warren@wpratt.com
 */
 
 #include "comm.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "eq.hpp"
 #include "fmsq.hpp"
 #include "RXA.hpp"

wdsp/icfir.cpp

@@ -26,7 +26,7 @@ warren@pratt.one
 */
 
 #include "comm.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "fir.hpp"
 #include "icfir.hpp"
 

wdsp/nbp.cpp

@@ -27,7 +27,7 @@ warren@wpratt.com
 
 #include "comm.hpp"
 #include "fir.hpp"
-#include "firmin.hpp"
+#include "fircore.hpp"
 #include "bpsnba.hpp"
 #include "nbp.hpp"
 #include "RXA.hpp"

wdsp/patchpanel.cpp

@@ -32,7 +32,17 @@ warren@wpratt.com
 
 namespace WDSP {
 
-PANEL* PANEL::create_panel (int run, int size, float* in, float* out, float gain1, float gain2I, float gain2Q, int inselect, int copy)
+PANEL* PANEL::create_panel (
+    int run,
+    int size,
+    float* in,
+    float* out,
+    float gain1,
+    float gain2I,
+    float gain2Q,
+    int inselect,
+    int copy
+)
 {
     PANEL* a = new PANEL;
     a->run = run;

wdsp/patchpanel.hpp

@@ -48,7 +48,17 @@ public:
     int inselect;
     int copy;
 
-    static PANEL* create_panel (int run, int size, float* in, float* out, float gain1, float gain2I, float gain2Q, int inselect, int copy);
+    static PANEL* create_panel (
+        int run,
+        int size,
+        float* in,
+        float* out,
+        float gain1,
+        float gain2I,
+        float gain2Q,
+        int inselect,
+        int copy
+    );
     static void destroy_panel (PANEL *a);
     static void flush_panel (PANEL *a);
     static void xpanel (PANEL *a);