1
0
mirror of https://github.com/f4exb/sdrangel.git synced 2024-11-21 15:51:47 -05:00

Compare commits

...

9 Commits

Author SHA1 Message Date
Edouard Griffiths
9e1033d68f
Merge de756413e8 into fcd43df711 2024-08-03 12:00:25 +00:00
f4exb
de756413e8 WDSP: RXA and TXA rework 2024-08-03 13:54:42 +02:00
f4exb
d6159067a8 WDSP: Sonar lint fixes (2) 2024-08-03 11:05:12 +02:00
Edouard Griffiths
fcd43df711
Merge pull request #2225 from james-pcdr/patch-2
Update scanner.py: fix print to stderr
2024-08-03 03:00:40 +02:00
Edouard Griffiths
ff7c06311b
Merge pull request #2226 from jfdelnero/master
Fix audio glitches in the DAB plugin
2024-08-02 18:56:01 +02:00
Jean-François DEL NERO
5888645957 Don't decimate the audio signal if not needed. 2024-08-02 12:23:27 +02:00
Jean-François DEL NERO
2fddaff6d2 Fix audio glitches in the DAB plugin.
The glitches were generated by an int16 integer overflow.

The issue appeared when the audio was near or at the saturation level.
When the input audio signal is saturated, the polyphase filter based interpolation/decimation functions tend to increase the samples values and then make them pass the int16 limits. The int16 sample scale() parameter defeat the min/max limitation.

This fix removes the intermediate int16 type conversion by using the Complex Real type.

fixes f4exb/sdrangel#1978
2024-08-02 12:17:00 +02:00
f4exb
8941835466 WDSP: Sonar lint fixes (1) 2024-08-02 08:01:46 +02:00
james-pcdr
4cf2c0b7c7
Update scanner.py: fix print to stderr 2024-07-31 13:33:28 -04:00
68 changed files with 2976 additions and 3028 deletions

1
.gitignore vendored
View File

@ -28,6 +28,7 @@ obj-x86_64-linux-gnu/*
**/venv*/
*.pyc
.DS_Store
.sonarlint
### Go ###
# Binaries for programs and plugins

View File

@ -366,7 +366,7 @@ void DABDemodSink::tii(int tii)
}
}
static int16_t scale(int16_t sample, float factor)
static int16_t scale(Real sample, float factor)
{
int32_t prod = (int32_t)(((int32_t)sample) * factor);
prod = std::min(prod, 32767);
@ -403,7 +403,12 @@ void DABDemodSink::audio(int16_t *buffer, int size, int samplerate, bool stereo)
ci.real(0.0f);
ci.imag(0.0f);
}
if (m_audioInterpolatorDistance < 1.0f) // interpolate
if (m_audioInterpolatorDistance == 1.0f)
{
processOneAudioSample(ci);
}
else if (m_audioInterpolatorDistance < 1.0f) // interpolate
{
while (!m_audioInterpolator.interpolate(&m_audioInterpolatorDistanceRemain, ci, &ca))
{

View File

@ -125,14 +125,14 @@ WDSPRxSink::WDSPRxSink() :
m_sPeak = 0.0;
m_sCount = m_wdspBufSize;
m_rxa = WDSP::RXA::create_rxa(
m_rxa = new WDSP::RXA(
m_wdspSampleRate, // input samplerate
m_wdspSampleRate, // output samplerate
m_wdspSampleRate, // sample rate for mainstream dsp processing (dsp)
m_wdspBufSize // number complex samples processed per buffer in mainstream dsp processing
);
m_rxa->setSpectrumProbe(&m_spectrumProbe);
WDSP::RXA::SetPassband(*m_rxa, 0, m_Bandwidth);
m_rxa->setPassband(0, m_Bandwidth);
applyChannelSettings(m_channelSampleRate, m_channelFrequencyOffset, true);
applySettings(m_settings, true);
@ -140,7 +140,7 @@ WDSPRxSink::WDSPRxSink() :
WDSPRxSink::~WDSPRxSink()
{
WDSP::RXA::destroy_rxa(m_rxa);
delete m_rxa;
}
void WDSPRxSink::feed(const SampleVector::const_iterator& begin, const SampleVector::const_iterator& end)
@ -189,7 +189,7 @@ void WDSPRxSink::processOneSample(Complex &ci)
if (++m_inCount == m_rxa->get_insize())
{
WDSP::RXA::xrxa(m_rxa);
m_rxa->execute();
m_sCount = m_wdspBufSize;
m_sAvg = m_rxa->smeter->getMeter(WDSP::RXA::RXA_S_AV);
@ -306,7 +306,7 @@ void WDSPRxSink::applyAudioSampleRate(int sampleRate)
m_interpolatorDistanceRemain = 0;
m_interpolatorDistance = (Real) m_channelSampleRate / (Real) m_wdspSampleRate;
WDSP::RXA::setOutputSamplerate(m_rxa, sampleRate);
m_rxa->setOutputSamplerate(sampleRate);
m_audioFifo.setSize(sampleRate);
m_audioSampleRate = sampleRate;
@ -446,24 +446,24 @@ void WDSPRxSink::applySettings(const WDSPRxSettings& settings, bool force)
m_interpolatorDistanceRemain = 0;
m_interpolatorDistance = (Real) m_channelSampleRate / (Real) m_audioSampleRate;
WDSP::RXA::SetPassband(*m_rxa, fLow, fHigh);
WDSP::RXA::NBPSetWindow(*m_rxa, m_settings.m_profiles[m_settings.m_profileIndex].m_fftWindow);
m_rxa->setPassband(fLow, fHigh);
m_rxa->nbpSetWindow(m_settings.m_profiles[m_settings.m_profileIndex].m_fftWindow);
if (settings.m_demod == WDSPRxProfile::DemodSSB)
{
if (dsb) {
WDSP::RXA::SetMode(*m_rxa, WDSP::RXA::RXA_DSB);
m_rxa->setMode(WDSP::RXA::RXA_DSB);
} else {
WDSP::RXA::SetMode(*m_rxa, usb ? WDSP::RXA::RXA_USB : WDSP::RXA::RXA_LSB);
m_rxa->setMode(usb ? WDSP::RXA::RXA_USB : WDSP::RXA::RXA_LSB);
}
}
else if (settings.m_demod == WDSPRxProfile::DemodAM)
{
WDSP::RXA::SetMode(*m_rxa, WDSP::RXA::RXA_AM);
m_rxa->setMode(WDSP::RXA::RXA_AM);
}
else if (settings.m_demod == WDSPRxProfile::DemodSAM)
{
WDSP::RXA::SetMode(*m_rxa, WDSP::RXA::RXA_SAM);
m_rxa->setMode(WDSP::RXA::RXA_SAM);
if (dsb) {
m_rxa->amd->setSBMode(0);
@ -473,7 +473,7 @@ void WDSPRxSink::applySettings(const WDSPRxSettings& settings, bool force)
}
else if (settings.m_demod == WDSPRxProfile::DemodFMN)
{
WDSP::RXA::SetMode(*m_rxa, WDSP::RXA::RXA_FM);
m_rxa->setMode(WDSP::RXA::RXA_FM);
}
}
@ -486,18 +486,18 @@ void WDSPRxSink::applySettings(const WDSPRxSettings& settings, bool force)
if ((m_settings.m_dnr != settings.m_dnr)
|| (m_settings.m_nrScheme != settings.m_nrScheme) || force)
{
WDSP::RXA::SetANRRun(*m_rxa, 0);
WDSP::RXA::SetEMNRRun(*m_rxa, 0);
m_rxa->setANRRun(0);
m_rxa->setEMNRRun(0);
if (settings.m_dnr)
{
switch (settings.m_nrScheme)
{
case WDSPRxProfile::NRSchemeNR:
WDSP::RXA::SetANRRun(*m_rxa, 1);
m_rxa->setANRRun(1);
break;
case WDSPRxProfile::NRSchemeNR2:
WDSP::RXA::SetEMNRRun(*m_rxa, 1);
m_rxa->setEMNRRun(1);
break;
default:
break;
@ -510,12 +510,12 @@ void WDSPRxSink::applySettings(const WDSPRxSettings& settings, bool force)
switch (settings.m_nrPosition)
{
case WDSPRxProfile::NRPositionPreAGC:
WDSP::RXA::SetANRPosition(*m_rxa, 0);
WDSP::RXA::SetEMNRPosition(*m_rxa, 0);
m_rxa->setANRPosition(0);
m_rxa->setEMNRPosition(0);
break;
case WDSPRxProfile::NRPositionPostAGC:
WDSP::RXA::SetANRPosition(*m_rxa, 1);
WDSP::RXA::SetEMNRPosition(*m_rxa, 1);
m_rxa->setANRPosition(1);
m_rxa->setEMNRPosition(1);
break;
default:
break;
@ -560,12 +560,12 @@ void WDSPRxSink::applySettings(const WDSPRxSettings& settings, bool force)
}
if ((m_settings.m_anf != settings.m_anf) || force) {
WDSP::RXA::SetANFRun(*m_rxa, settings.m_anf ? 1 : 0);
m_rxa->setANFRun(settings.m_anf ? 1 : 0);
}
// Caution: Causes corruption
if ((m_settings.m_snb != settings.m_snb) || force) {
WDSP::RXA::SetSNBARun(*m_rxa, settings.m_snb ? 1 : 0);
m_rxa->setSNBARun(settings.m_snb ? 1 : 0);
}
// CW Peaking

View File

@ -373,7 +373,7 @@ def main():
pass
except Exception as ex:
tb = traceback.format_exc()
print >> sys.stderr, tb
print(tb, file=sys.stderr)
if __name__ == "__main__":

View File

@ -62,6 +62,7 @@ set(wdsp_SOURCES
sphp.cpp
ssql.cpp
TXA.cpp
unit.cpp
varsamp.cpp
wcpAGC.cpp
)
@ -129,6 +130,7 @@ set(wdsp_HEADERS
sphp.hpp
ssql.hpp
TXA.hpp
unit.hpp
varsamp.hpp
wcpAGC.hpp
)

File diff suppressed because it is too large Load Diff

View File

@ -128,70 +128,68 @@ public:
PANEL *panel;
RESAMPLE *rsmpout;
static RXA* create_rxa (
RXA(
int in_rate, // input samplerate
int out_rate, // output samplerate
int dsp_rate, // sample rate for mainstream dsp processing
int dsp_size // number complex samples processed per buffer in mainstream dsp processing
);
static void destroy_rxa (RXA *rxa);
static void flush_rxa (RXA *rxa);
static void xrxa (RXA *rxa);
int get_insize() const { return dsp_insize; }
int get_outsize() const { return dsp_outsize; }
float *get_inbuff() { return inbuff; }
float *get_outbuff() { return outbuff; }
RXA(const RXA&) = delete;
RXA& operator=(const RXA& other) = delete;
~RXA();
void flush();
void execute();
void setInputSamplerate(int _in_rate);
void setOutputSamplerate(int _out_rate);
void setDSPSamplerate(int _dsp_rate);
void setDSPBuffsize(int _dsp_size);
int get_insize() const { return Unit::dsp_insize; }
int get_outsize() const { return Unit::dsp_outsize; }
float *get_inbuff() { return Unit::inbuff; }
float *get_outbuff() { return Unit::outbuff; }
void setSpectrumProbe(BufferProbe *_spectrumProbe);
static void setInputSamplerate (RXA *rxa, int in_rate);
static void setOutputSamplerate (RXA *rxa, int out_rate);
static void setDSPSamplerate (RXA *rxa, int dsp_rate);
static void setDSPBuffsize (RXA *rxa, int dsp_size);
// RXA Properties
static void SetMode (RXA& rxa, int mode);
static void ResCheck (RXA& rxa);
static void bp1Check (RXA& rxa, int amd_run, int snba_run, int emnr_run, int anf_run, int anr_run);
static void bp1Set (RXA& rxa);
static void bpsnbaCheck (RXA& rxa, int mode, int notch_run);
static void bpsnbaSet (RXA& rxa);
void setMode (int mode);
void resCheck ();
void bp1Check (int amd_run, int snba_run, int emnr_run, int anf_run, int anr_run);
void bp1Set ();
void bpsnbaCheck (int mode, int notch_run);
void bpsnbaSet ();
// NOTCHDB, NBP, SNBA
static void UpdateNBPFiltersLightWeight (RXA& rxa);
static void UpdateNBPFilters(RXA& rxa);
static int NBPAddNotch (RXA& rxa, int notch, double fcenter, double fwidth, int active);
static int NBPGetNotch (RXA& rxa, int notch, double* fcenter, double* fwidth, int* active);
static int NBPDeleteNotch (RXA& rxa, int notch);
static int NBPEditNotch (RXA& rxa, int notch, double fcenter, double fwidth, int active);
static void NBPGetNumNotches (RXA& rxa, int* nnotches);
static void NBPSetTuneFrequency (RXA& rxa, double tunefreq);
static void NBPSetShiftFrequency (RXA& rxa, double shift);
static void NBPSetNotchesRun (RXA& rxa, int run);
static void NBPSetWindow (RXA& rxa, int wintype);
static void NBPSetAutoIncrease (RXA& rxa, int autoincr);
void updateNBPFiltersLightWeight();
void updateNBPFilters();
int nbpAddNotch(int notch, double fcenter, double fwidth, int active);
int nbpGetNotch(int notch, double* fcenter, double* fwidth, int* active);
int nbpDeleteNotch(int notch);
int nbpEditNotch(int notch, double fcenter, double fwidth, int active);
void nbpGetNumNotches(int* nnotches);
void nbpSetTuneFrequency(double tunefreq);
void nbpSetShiftFrequency(double shift);
void nbpSetNotchesRun(int run);
void nbpSetWindow(int wintype);
void nbpSetAutoIncrease(int autoincr);
// AMD
static void SetAMDRun(RXA& rxa, int run);
void setAMDRun(int run);
// SNBA
static void SetSNBARun (RXA& rxa, int run);
void setSNBARun(int run);
// ANF
static void SetANFRun (RXA& rxa, int run);
static void SetANFPosition (RXA& rxa, int position);
void setANFRun(int run);
void setANFPosition(int position);
// ANR
static void SetANRRun (RXA& rxa, int run);
static void SetANRPosition (RXA& rxa, int position);
void setANRRun(int run);
void setANRPosition(int position);
// EMNR
static void SetEMNRRun (RXA& rxa, int run);
static void SetEMNRPosition (RXA& rxa, int position);
void setEMNRRun(int run);
void setEMNRPosition(int position);
// WCPAGC
static void SetAGCThresh(RXA& rxa, double thresh, double size, double rate);
static void GetAGCThresh(RXA& rxa, double *thresh, double size, double rate);
void setAGCThresh(double thresh, double size, double rate);
void getAGCThresh(double *thresh, double size, double rate);
// Collectives
static void SetPassband (RXA& rxa, float f_low, float f_high);
static void SetNC (RXA& rxa, int nc);
static void SetMP (RXA& rxa, int mp);
private:
float* inbuff;
float* midbuff;
float* outbuff;
void setPassband(float f_low, float f_high);
void setNC(int nc);
void setMP(int mp);
};
} // namespace WDSP

File diff suppressed because it is too large Load Diff

View File

@ -157,11 +157,6 @@ public:
GEN *gen0;
GEN *gen1;
USLEW *uslew;
// struct
// {
// CALCC *p;
// CRITICAL_SECTION cs_update;
// } calcc;
struct
{
IQC *p0, *p1;
@ -169,42 +164,42 @@ public:
} iqc;
CFIR *cfir;
static TXA* create_txa (
TXA(
int in_rate, // input samplerate
int out_rate, // output samplerate
int dsp_rate, // sample rate for mainstream dsp processing
int dsp_size // number complex samples processed per buffer in mainstream dsp processing
);
static void destroy_txa (TXA *txa);
static void flush_txa (TXA *txa);
static void xtxa (TXA *txa);
int get_insize() const { return dsp_insize; }
int get_outsize() const { return dsp_outsize; }
float *get_inbuff() { return inbuff; }
float *get_outbuff() { return outbuff; }
static void setInputSamplerate (TXA *txa, int in_rate);
static void setOutputSamplerate (TXA *txa, int out_rate);
static void setDSPSamplerate (TXA *txa, int dsp_rate);
static void setDSPBuffsize (TXA *txa, int dsp_size);
TXA(const TXA&) = delete;
TXA& operator=(const TXA& other) = delete;
~TXA();
void flush();
void execute();
void setInputSamplerate(int _in_rate);
void setOutputSamplerate(int _out_rate);
void setDSPSamplerate(int _dsp_rate);
void setDSPBuffsize(int _dsp_size);
int get_insize() const { return Unit::dsp_insize; }
int get_outsize() const { return Unit::dsp_outsize; }
float *get_inbuff() { return Unit::inbuff; }
float *get_outbuff() { return Unit::outbuff; }
// TXA Properties
static void SetMode (TXA& txa, int mode);
static void SetBandpassFreqs (TXA& txa, float f_low, float f_high);
static void SetBandpassNC (TXA& txa, int nc);
static void SetBandpassMP (TXA& txa, int mp);
void setMode(int mode);
void setBandpassFreqs(float f_low, float f_high);
void setBandpassNC(int nc);
void setBandpassMP(int mp);
// Collectives
static void SetNC (TXA& txa, int nc);
static void SetMP (TXA& txa, int mp);
static void SetFMAFFilter (TXA& txa, float low, float high);
static void SetupBPFilters (TXA& txa);
static int UslewCheck (TXA& txa);
void setNC(int nc);
void setMP(int mp);
void setFMAFFilter(float low, float high);
void setupBPFilters();
int uslewCheck();
private:
static void ResCheck (TXA& txa);
float* inbuff;
float* midbuff;
float* outbuff;
void resCheck();
};
} // namespace WDSP

View File

@ -26,6 +26,7 @@ warren@wpratt.com
*/
#include <cmath>
#include <array>
#include "comm.hpp"
#include "amd.hpp"
@ -116,15 +117,19 @@ void AMD::flush()
void AMD::execute()
{
int i;
double audio;
double vco[2];
double corr[2];
std::array<double, 2> vco;
std::array<double, 2> corr;
double det;
double del_out;
double ai, bi, aq, bq;
double ai_ps, bi_ps, aq_ps, bq_ps;
int j, k;
double ai;
double bi;
double aq;
double bq;
double ai_ps;
double bi_ps;
double aq_ps;
double bq_ps;
if (run)
{
@ -133,7 +138,7 @@ void AMD::execute()
case 0: //AM Demodulator
{
for (i = 0; i < buff_size; i++)
for (int i = 0; i < buff_size; i++)
{
double xr = in_buff[2 * i + 0];
double xi = in_buff[2 * i + 1];
@ -146,8 +151,8 @@ void AMD::execute()
audio += dc_insert - dc;
}
out_buff[2 * i + 0] = audio;
out_buff[2 * i + 1] = audio;
out_buff[2 * i + 0] = (float) audio;
out_buff[2 * i + 1] = (float) audio;
}
break;
@ -155,7 +160,7 @@ void AMD::execute()
case 1: //Synchronous AM Demodulator with Sideband Separation
{
for (i = 0; i < buff_size; i++)
for (int i = 0; i < buff_size; i++)
{
vco[0] = cos(phs);
vco[1] = sin(phs);
@ -174,9 +179,9 @@ void AMD::execute()
dsI = ai;
dsQ = bq;
for (j = 0; j < STAGES; j++)
for (int j = 0; j < STAGES; j++)
{
k = 3 * j;
int k = 3 * j;
a[k + 3] = c0[j] * (a[k] - a[k + 5]) + a[k + 2];
b[k + 3] = c1[j] * (b[k] - b[k + 5]) + b[k + 2];
c[k + 3] = c0[j] * (c[k] - c[k + 5]) + c[k + 2];
@ -188,7 +193,7 @@ void AMD::execute()
bq_ps = c[OUT_IDX];
aq_ps = d[OUT_IDX];
for (j = OUT_IDX + 2; j > 0; j--)
for (int j = OUT_IDX + 2; j > 0; j--)
{
a[j] = a[j - 1];
b[j] = b[j - 1];
@ -217,6 +222,8 @@ void AMD::execute()
audio = (ai_ps + bi_ps) - (aq_ps - bq_ps);
break;
}
default:
break;
}
if (levelfade)
@ -226,8 +233,8 @@ void AMD::execute()
audio += dc_insert - dc;
}
out_buff[2 * i + 0] = audio;
out_buff[2 * i + 1] = audio;
out_buff[2 * i + 0] = (float) audio;
out_buff[2 * i + 1] = (float) audio;
if ((corr[0] == 0.0) && (corr[1] == 0.0))
corr[0] = 1.0;
@ -254,6 +261,8 @@ void AMD::execute()
break;
}
default:
break;
}
}
else if (in_buff != out_buff)

View File

@ -28,15 +28,6 @@ warren@wpratt.com
#ifndef wdsp_amd_hpp
#define wdsp_amd_hpp
// ff defines for sbdemod
#ifndef STAGES
#define STAGES 7
#endif
#ifndef OUT_IDX
#define OUT_IDX (3 * STAGES)
#endif
#include <array>
#include "export.h"
@ -61,13 +52,16 @@ public:
double phs; // pll - phase accumulator
double omega; // pll - locked pll frequency
double fil_out; // pll - filter output
double g1, g2; // pll - filter gain parameters
double g1; // pll - filter gain parameters
double g2; // pll - filter gain parameters
double tauR; // carrier removal time constant
double tauI; // carrier insertion time constant
double mtauR; // carrier removal multiplier
double onem_mtauR; // 1.0 - carrier_removal_multiplier
double mtauI; // carrier insertion multiplier
double onem_mtauI; // 1.0 - carrier_insertion_multiplier
static const int STAGES = 7;
static const int OUT_IDX = 3 * STAGES;
std::array<double, 3*STAGES + 3> a; // Filter a variables
std::array<double, 3*STAGES + 3> b; // Filter b variables
std::array<double, 3*STAGES + 3> c; // Filter c variables

View File

@ -25,6 +25,8 @@ warren@wpratt.com
*/
#include <cstdio>
#include "comm.hpp"
#include "amsq.hpp"
@ -32,12 +34,12 @@ namespace WDSP {
void AMSQ::compute_slews()
{
int i;
double delta, theta;
double delta;
double theta;
delta = PI / (double)ntup;
theta = 0.0;
for (i = 0; i <= ntup; i++)
for (int i = 0; i <= ntup; i++)
{
cup[i] = muted_gain + (1.0 - muted_gain) * 0.5 * (1.0 - cos (theta));
theta += delta;
@ -46,7 +48,7 @@ void AMSQ::compute_slews()
delta = PI / (double)ntdown;
theta = 0.0;
for (i = 0; i <= ntdown; i++)
for (int i = 0; i <= ntdown; i++)
{
cdown[i] = muted_gain + (1.0 - muted_gain) * 0.5 * (1.0 + cos (theta));
theta += delta;
@ -63,11 +65,11 @@ void AMSQ::calc()
// level change
ntup = (int)(tup * rate);
ntdown = (int)(tdown * rate);
cup.resize((ntup + 1) * 2); // (float *)malloc0((ntup + 1) * sizeof(float));
cdown.resize((ntdown + 1) * 2); // (float *)malloc0((ntdown + 1) * sizeof(float));
cup.resize((ntup + 1) * 2);
cdown.resize((ntdown + 1) * 2);
compute_slews();
// control
state = 0;
state = AMSQState::MUTED;
}
AMSQ::AMSQ (
@ -108,26 +110,17 @@ void AMSQ::flush()
{
std::fill(trigsig.begin(), trigsig.end(), 0);
avsig = 0.0;
state = 0;
state = AMSQState::MUTED;
}
enum _amsqstate
{
MUTED,
INCREASE,
UNMUTED,
TAIL,
DECREASE
};
void AMSQ::execute()
{
if (run)
{
int i;
double sig, siglimit;
double sig;
double siglimit;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
double trigr = trigsig[2 * i + 0];
double trigi = trigsig[2 * i + 1];
@ -136,31 +129,31 @@ void AMSQ::execute()
switch (state)
{
case MUTED:
case AMSQState::MUTED:
if (avsig > unmute_thresh)
{
state = INCREASE;
state = AMSQState::INCREASE;
count = ntup;
}
out[2 * i + 0] = muted_gain * in[2 * i + 0];
out[2 * i + 1] = muted_gain * in[2 * i + 1];
out[2 * i + 0] = (float) (muted_gain * in[2 * i + 0]);
out[2 * i + 1] = (float) (muted_gain * in[2 * i + 1]);
break;
case INCREASE:
out[2 * i + 0] = in[2 * i + 0] * cup[ntup - count];
out[2 * i + 1] = in[2 * i + 1] * cup[ntup - count];
case AMSQState::INCREASE:
out[2 * i + 0] = (float) (in[2 * i + 0] * cup[ntup - count]);
out[2 * i + 1] = (float) (in[2 * i + 1] * cup[ntup - count]);
if (count-- == 0)
state = UNMUTED;
state = AMSQState::UNMUTED;
break;
case UNMUTED:
case AMSQState::UNMUTED:
if (avsig < tail_thresh)
{
state = TAIL;
state = AMSQState::TAIL;
if ((siglimit = avsig) > 1.0)
siglimit = 1.0;
@ -173,28 +166,28 @@ void AMSQ::execute()
break;
case TAIL:
case AMSQState::TAIL:
out[2 * i + 0] = in[2 * i + 0];
out[2 * i + 1] = in[2 * i + 1];
if (avsig > unmute_thresh)
{
state = UNMUTED;
state = AMSQState::UNMUTED;
}
else if (count-- == 0)
{
state = DECREASE;
state = AMSQState::DECREASE;
count = ntdown;
}
break;
case DECREASE:
out[2 * i + 0] = in[2 * i + 0] * cdown[ntdown - count];
out[2 * i + 1] = in[2 * i + 1] * cdown[ntdown - count];
case AMSQState::DECREASE:
out[2 * i + 0] = (float) (in[2 * i + 0] * cdown[ntdown - count]);
out[2 * i + 1] = (float) (in[2 * i + 1] * cdown[ntdown - count]);
if (count-- == 0)
state = MUTED;
state = AMSQState::MUTED;
break;
}

View File

@ -36,6 +36,15 @@ namespace WDSP {
class WDSP_API AMSQ
{
public:
enum class AMSQState
{
MUTED,
INCREASE,
UNMUTED,
TAIL,
DECREASE
};
int run; // 0 if squelch system is OFF; 1 if it's ON
int size; // size of input/output buffers
float* in; // squelch input signal buffer
@ -47,7 +56,7 @@ public:
double avm;
double onem_avm;
double avsig;
int state; // state machine control
AMSQState state; // state machine control
int count;
double tup;
double tdown;

View File

@ -36,7 +36,6 @@ namespace WDSP {
void ANB::initBlanker()
{
int i;
trans_count = (int)(tau * samplerate);
if (trans_count < 2)
@ -54,7 +53,7 @@ void ANB::initBlanker()
backmult = exp(-1.0 / (samplerate * backtau));
ombackmult = 1.0 - backmult;
for (i = 0; i <= trans_count; i++)
for (int i = 0; i <= trans_count; i++)
wave[i] = 0.5 * cos(i * coef);
std::fill(dline.begin(), dline.end(), 0);
@ -76,23 +75,44 @@ ANB::ANB (
buffsize(_buffsize),
in(_in),
out(_out),
dline_size((int)((MAX_TAU + MAX_ADVTIME) * MAX_SAMPLERATE) + 1),
samplerate(_samplerate),
tau(_tau),
hangtime(_hangtime),
advtime(_advtime),
backtau(_backtau),
threshold(_threshold),
dtime(0),
htime(0),
itime(0),
atime(0)
threshold(_threshold)
{
tau = tau < 0.0 ? 0.0 : (tau > MAX_TAU ? MAX_TAU : tau);
hangtime = hangtime < 0.0 ? 0.0 : (hangtime > MAX_ADVTIME ? MAX_ADVTIME : hangtime);
advtime = advtime < 0.0 ? 0.0 : (advtime > MAX_ADVTIME ? MAX_ADVTIME : advtime);
samplerate = samplerate < 0.0 ? 0.0 : (samplerate > MAX_SAMPLERATE ? MAX_SAMPLERATE : samplerate);
dtime = 0;
htime = 0;
itime = 0;
atime = 0;
if (tau < 0.0) {
tau = 0.0;
} else if (tau > MAX_TAU) {
tau = MAX_TAU;
}
if (hangtime < 0.0) {
hangtime = 0.0;
} else if (hangtime > MAX_ADVTIME) {
hangtime = MAX_ADVTIME;
}
if (advtime < 0.0) {
advtime = 0.0;
} else if (advtime > MAX_ADVTIME) {
advtime = MAX_ADVTIME;
}
if (samplerate < 0.0) {
samplerate = 0.0;
} else if (samplerate > MAX_SAMPLERATE) {
samplerate = MAX_SAMPLERATE;
}
wave.resize((int)(MAX_SAMPLERATE * MAX_TAU) + 1);
dline_size = (int)((MAX_TAU + MAX_ADVTIME) * MAX_SAMPLERATE) + 1;
dline.resize(dline_size * 2);
initBlanker();
}
@ -106,11 +126,10 @@ void ANB::execute()
{
double scale;
double mag;
int i;
if (run)
{
for (i = 0; i < buffsize; i++)
for (int i = 0; i < buffsize; i++)
{
double xr = in[2 * i + 0];
double xi = in[2 * i + 1];
@ -139,8 +158,8 @@ void ANB::execute()
case 1:
scale = power * (0.5 + wave[dtime]);
out[2 * i + 0] = dline[2 * out_idx + 0] * scale;
out[2 * i + 1] = dline[2 * out_idx + 1] * scale;
out[2 * i + 0] = (float) (dline[2 * out_idx + 0] * scale);
out[2 * i + 1] = (float) (dline[2 * out_idx + 1] * scale);
if (++dtime > trans_count)
{
@ -177,8 +196,8 @@ void ANB::execute()
case 4:
scale = 0.5 - wave[itime];
out[2 * i + 0] = dline[2 * out_idx + 0] * scale;
out[2 * i + 1] = dline[2 * out_idx + 1] * scale;
out[2 * i + 0] = (float) (dline[2 * out_idx + 0] * scale);
out[2 * i + 1] = (float) (dline[2 * out_idx + 1] * scale);
if (count > 0)
{

View File

@ -53,49 +53,53 @@ ANF::ANF(
double _den_mult,
double _lincr,
double _ldecr
)
) :
run(_run),
position(_position),
buff_size(_buff_size),
in_buff(_in_buff),
out_buff(_out_buff),
dline_size(_dline_size),
mask(_dline_size - 1),
n_taps(_n_taps),
delay(_delay),
two_mu(_two_mu),
gamma(_gamma),
lidx(_lidx),
lidx_min(_lidx_min),
lidx_max(_lidx_max),
ngamma(_ngamma),
den_mult(_den_mult),
lincr(_lincr),
ldecr(_ldecr)
{
run = _run;
position = _position;
buff_size = _buff_size;
in_buff = _in_buff;
out_buff = _out_buff;
dline_size = _dline_size;
mask = _dline_size - 1;
n_taps = _n_taps;
delay = _delay;
two_mu = _two_mu;
gamma = _gamma;
in_idx = 0;
lidx = _lidx;
lidx_min = _lidx_min;
lidx_max = _lidx_max;
ngamma = _ngamma;
den_mult = _den_mult;
lincr = _lincr;
ldecr = _ldecr;
std::fill(d.begin(), d.end(), 0);
std::fill(w.begin(), w.end(), 0);
}
void ANF::execute(int _position)
{
int i, j, idx;
double c0, c1;
double y, error, sigma, inv_sigp;
double nel, nev;
int idx;
double c0;
double c1;
double y;
double error;
double sigma;
double inv_sigp;
double nel;
double nev;
if (run && (position == _position))
{
for (i = 0; i < buff_size; i++)
for (int i = 0; i < buff_size; i++)
{
d[in_idx] = in_buff[2 * i + 0];
y = 0;
sigma = 0;
for (j = 0; j < n_taps; j++)
for (int j = 0; j < n_taps; j++)
{
idx = (in_idx + j + delay) & mask;
y += w[j] * d[idx];
@ -105,7 +109,7 @@ void ANF::execute(int _position)
inv_sigp = 1.0 / (sigma + 1e-10);
error = d[in_idx] - y;
out_buff[2 * i + 0] = error;
out_buff[2 * i + 0] = (float) error;
out_buff[2 * i + 1] = 0.0;
if ((nel = error * (1.0 - two_mu * sigma * inv_sigp)) < 0.0)
@ -128,7 +132,7 @@ void ANF::execute(int _position)
c0 = 1.0 - two_mu * ngamma;
c1 = two_mu * error * inv_sigp;
for (j = 0; j < n_taps; j++)
for (int j = 0; j < n_taps; j++)
{
idx = (in_idx + j + delay) & mask;
w[j] = c0 * w[j] + c1 * d[idx];

View File

@ -32,8 +32,6 @@ warren@wpratt.com
#include "export.h"
#define ANF_DLINE_SIZE 2048
namespace WDSP {
class WDSP_API ANF
@ -50,6 +48,7 @@ public:
int delay;
double two_mu;
double gamma;
static const int ANF_DLINE_SIZE = 2048;
std::array<double, ANF_DLINE_SIZE> d;
std::array<double, ANF_DLINE_SIZE> w;
int in_idx;

View File

@ -80,21 +80,26 @@ ANR::ANR(
void ANR::execute(int _position)
{
int i, j, idx;
double c0, c1;
double y, error, sigma, inv_sigp;
double nel, nev;
int idx;
double c0;
double c1;
double y;
double error;
double sigma;
double inv_sigp;
double nel;
double nev;
if (run && (position == _position))
{
for (i = 0; i < buff_size; i++)
for (int i = 0; i < buff_size; i++)
{
d[in_idx] = in_buff[2 * i + 0];
y = 0;
sigma = 0;
for (j = 0; j < n_taps; j++)
for (int j = 0; j < n_taps; j++)
{
idx = (in_idx + j + delay) & mask;
y += w[j] * d[idx];
@ -104,7 +109,7 @@ void ANR::execute(int _position)
inv_sigp = 1.0 / (sigma + 1e-10);
error = d[in_idx] - y;
out_buff[2 * i + 0] = y;
out_buff[2 * i + 0] = (float) y;
out_buff[2 * i + 1] = 0.0;
if ((nel = error * (1.0 - two_mu * sigma * inv_sigp)) < 0.0)
@ -127,7 +132,7 @@ void ANR::execute(int _position)
c0 = 1.0 - two_mu * ngamma;
c1 = two_mu * error * inv_sigp;
for (j = 0; j < n_taps; j++)
for (int j = 0; j < n_taps; j++)
{
idx = (in_idx + j + delay) & mask;
w[j] = c0 * w[j] + c1 * d[idx];

View File

@ -32,8 +32,6 @@ warren@wpratt.com
#include "export.h"
#define ANR_DLINE_SIZE 2048
namespace WDSP {
class WDSP_API ANR
@ -50,6 +48,7 @@ public:
int delay;
double two_mu;
double gamma;
static const int ANR_DLINE_SIZE = 2048;
std::array<double, ANR_DLINE_SIZE> d;
std::array<double, ANR_DLINE_SIZE> w;
int in_idx;

View File

@ -51,21 +51,21 @@ BANDPASS::BANDPASS(
int _samplerate,
int _wintype,
double _gain
)
{
) :
// NOTE: 'nc' must be >= 'size'
run = _run;
position = _position;
size = _size;
nc = _nc;
mp = _mp;
in = _in;
out = _out;
f_low = _f_low;
f_high = _f_high;
samplerate = _samplerate;
wintype = _wintype;
gain = _gain;
run(_run),
position(_position),
size(_size),
nc(_nc),
mp(_mp),
in(_in),
out(_out),
f_low(_f_low),
f_high(_f_high),
samplerate(_samplerate),
wintype(_wintype),
gain(_gain)
{
float* impulse = FIR::fir_bandpass (
nc,
f_low,
@ -136,7 +136,7 @@ void BANDPASS::setSize(int _size)
gain / (double) (2 * size)
);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void BANDPASS::setGain(double _gain, int _update)
@ -152,7 +152,7 @@ void BANDPASS::setGain(double _gain, int _update)
gain / (double) (2 * size)
);
FIRCORE::setImpulse_fircore (fircore, impulse, _update);
delete[] (impulse);
delete[] impulse;
}
void BANDPASS::calcBandpassFilter(double _f_low, double _f_high, double _gain)
@ -172,7 +172,7 @@ void BANDPASS::calcBandpassFilter(double _f_low, double _f_high, double _gain)
gain / (double)(2 * size)
);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
}
@ -197,7 +197,7 @@ void BANDPASS::setBandpassFreqs(double _f_low, double _f_high)
);
FIRCORE::setImpulse_fircore (fircore, impulse, 0);
delete[] (impulse);
delete[] impulse;
f_low = _f_low;
f_high = _f_high;
FIRCORE::setUpdate_fircore (fircore);
@ -220,7 +220,7 @@ void BANDPASS::SetBandpassNC(int _nc)
gain / (double)( 2 * size)
);
FIRCORE::setNc_fircore (fircore, nc, impulse);
delete[] (impulse);
delete[] impulse;
}
}
@ -239,38 +239,4 @@ void BANDPASS::SetBandpassMP(int _mp)
* *
********************************************************************************************************/
//PORT
//void SetTXABandpassFreqs (int channel, float f_low, float f_high)
//{
// float* impulse;
// BANDPASS a;
// a = txa.bp0;
// if ((f_low != a->f_low) || (f_high != a->f_high))
// {
// a->f_low = f_low;
// a->f_high = f_high;
// impulse = fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
// setImpulse_fircore (a->fircore, impulse, 1);
// delete[] (impulse);
// }
// a = txa.bp1;
// if ((f_low != a->f_low) || (f_high != a->f_high))
// {
// a->f_low = f_low;
// a->f_high = f_high;
// impulse = fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
// setImpulse_fircore (a->fircore, impulse, 1);
// delete[] (impulse);
// }
// a = txa.bp2;
// if ((f_low != a->f_low) || (f_high != a->f_high))
// {
// a->f_low = f_low;
// a->f_high = f_high;
// impulse = fir_bandpass (a->nc, a->f_low, a->f_high, a->samplerate, a->wintype, 1, a->gain / (float)(2 * a->size));
// setImpulse_fircore (a->fircore, impulse, 1);
// delete[] (impulse);
// }
//}
} // namespace WDSP

View File

@ -52,7 +52,7 @@ namespace WDSP {
void BPSNBA::calc()
{
buff = new float[size * 2]; // (double *) malloc0 (size * sizeof (complex));
buff.resize(size * 2);
bpsnba = new NBP (
1, // run, always runs (use bpsnba 'run')
run_notches, // run the notches
@ -60,7 +60,7 @@ void BPSNBA::calc()
size, // buffer size
nc, // number of filter coefficients
mp, // minimum phase flag
buff, // pointer to input buffer
buff.data(), // pointer to input buffer
out, // pointer to output buffer
f_low, // lower filter frequency
f_high, // upper filter frequency
@ -116,8 +116,7 @@ BPSNBA::BPSNBA(
void BPSNBA::decalc()
{
delete (bpsnba);
delete[] (buff);
delete bpsnba;
}
BPSNBA::~BPSNBA()
@ -127,7 +126,7 @@ BPSNBA::~BPSNBA()
void BPSNBA::flush()
{
std::fill(buff, buff + size * 2, 0);
std::fill(buff.begin(), buff.end(), 0);
bpsnba->flush();
}
@ -156,7 +155,7 @@ void BPSNBA::setSize(int _size)
void BPSNBA::exec_in(int _position)
{
if (run && position == _position)
std::copy(in, in + size * 2, buff);
std::copy(in, in + size * 2, buff.begin());
}
void BPSNBA::exec_out(int _position)

View File

@ -28,12 +28,15 @@ warren@wpratt.com
#ifndef wdsp_bpsnba_h
#define wdsp_bpsnba_h
#include <vector>
#include "export.h"
namespace WDSP{
class NOTCHDB;
class NBP;
class BPSNBA
class WDSP_API BPSNBA
{
public:
int run; // run the filter
@ -49,7 +52,7 @@ public:
double abs_high_freq; // highest positive freq supported by SNG
double f_low; // low cutoff frequency
double f_high; // high cutoff frequency
float* buff; // internal buffer
std::vector<float> buff; // internal buffer
int wintype; // filter window type
double gain; // filter gain
int autoincr; // use auto increment for notch width

View File

@ -71,22 +71,24 @@ void CBL::execute()
{
if (run)
{
int i;
double tempI, tempQ;
double tempI;
double tempQ;
for (i = 0; i < buff_size; i++)
for (int i = 0; i < buff_size; i++)
{
tempI = in_buff[2 * i + 0];
tempQ = in_buff[2 * i + 1];
out_buff[2 * i + 0] = in_buff[2 * i + 0] - prevIin + mtau * prevIout;
out_buff[2 * i + 1] = in_buff[2 * i + 1] - prevQin + mtau * prevQout;
out_buff[2 * i + 0] = (float) (in_buff[2 * i + 0] - prevIin + mtau * prevIout);
out_buff[2 * i + 1] = (float) (in_buff[2 * i + 1] - prevQin + mtau * prevQout);
prevIin = tempI;
prevQin = tempQ;
prevIout = out_buff[2 * i + 0];
prevQout = out_buff[2 * i + 1];
if (fabs(prevIout = out_buff[2 * i + 0]) < 1.0e-20)
if (fabs(prevIout) < 1.0e-20)
prevIout = 0.0;
if (fabs(prevQout = out_buff[2 * i + 1]) < 1.0e-20)
if (fabs(prevQout) < 1.0e-20)
prevQout = 0.0;
}
}

View File

@ -103,7 +103,7 @@ void COMPRESSOR::SetCompressorRun (TXA& txa, int run)
if (txa.compressor->run != run)
{
txa.compressor->run = run;
TXA::SetupBPFilters (txa);
txa.setupBPFilters();
}
}

View File

@ -97,7 +97,7 @@ void DSPHP::execute()
x1[n] = x0[n];
}
out[i] = y0[nstages - 1];
out[i] = (float) y0[nstages - 1];
}
}
else if (out != in)

View File

@ -50,8 +50,13 @@ public:
double rate;
double fc;
int nstages;
double a1, b0, b1;
std::vector<double> x0, x1, y0, y1;
double a1;
double b0;
double b1;
std::vector<double> x0;
std::vector<double> x1;
std::vector<double> y0;
std::vector<double> y1;
DSPHP(
int run,

View File

@ -72,10 +72,10 @@ EMNR::NPS::NPS(
epsH1(_epsH1)
{
epsH1r = epsH1 / (1.0 + epsH1);
sigma2N.resize(msize); // (float *)malloc0(nps.msize * sizeof(float));
PH1y.resize(msize); // (float *)malloc0(nps.msize * sizeof(float));
Pbar.resize(msize); // (float *)malloc0(nps.msize * sizeof(float));
EN2y.resize(msize); // (float *)malloc0(nps.msize * sizeof(float));
sigma2N.resize(msize);
PH1y.resize(msize);
Pbar.resize(msize);
EN2y.resize(msize);
for (int i = 0; i < msize; i++)
{
@ -86,9 +86,8 @@ EMNR::NPS::NPS(
void EMNR::NPS::LambdaDs()
{
int k;
for (k = 0; k < msize; k++)
for (int k = 0; k < msize; k++)
{
PH1y[k] = 1.0 / (1.0 + (1.0 + epsH1) * exp (- epsH1r * lambda_y[k] / sigma2N[k]));
Pbar[k] = alpha_Pbar * Pbar[k] + (1.0 - alpha_Pbar) * PH1y[k];
@ -122,31 +121,17 @@ EMNR::NP::NP(
lambda_d(_lambda_d)
{
{
double tau = -128.0 / 8000.0 / log(0.7);
alphaCsmooth = exp(-incr / rate / tau);
}
{
double tau = -128.0 / 8000.0 / log(0.96);
alphaMax = exp(-incr / rate / tau);
}
{
double tau = -128.0 / 8000.0 / log(0.7);
alphaCmin = exp(-incr / rate / tau);
}
{
double tau = -128.0 / 8000.0 / log(0.3);
alphaMin_max_value = exp(-incr / rate / tau);
}
double tau0 = -128.0 / 8000.0 / log(0.7);
alphaCsmooth = exp(-incr / rate / tau0);
double tau1 = -128.0 / 8000.0 / log(0.96);
alphaMax = exp(-incr / rate / tau1);
double tau2 = -128.0 / 8000.0 / log(0.7);
alphaCmin = exp(-incr / rate / tau2);
double tau3 = -128.0 / 8000.0 / log(0.3);
alphaMin_max_value = exp(-incr / rate / tau3);
snrq = -incr / (0.064 * rate);
double tau4 = -128.0 / 8000.0 / log(0.8);
betamax = exp(-incr / rate / tau4);
invQeqMax = 0.5;
av = 2.12;
Dtime = 8.0 * 12.0 * 128.0 / 8000.0;
@ -166,68 +151,63 @@ EMNR::NP::NP(
invQbar_points[0] = 0.03;
invQbar_points[1] = 0.05;
invQbar_points[2] = 0.06;
invQbar_points[3] = std::numeric_limits<double>::max();;
invQbar_points[3] = std::numeric_limits<double>::max();
{
double db;
db = 10.0 * log10(8.0) / (12.0 * 128 / 8000);
nsmax[0] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(4.0) / (12.0 * 128 / 8000);
nsmax[1] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(2.0) / (12.0 * 128 / 8000);
nsmax[2] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(1.2) / (12.0 * 128 / 8000);
nsmax[3] = pow(10.0, db / 10.0 * V * incr / rate);
}
double db;
db = 10.0 * log10(8.0) / (12.0 * 128 / 8000);
nsmax[0] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(4.0) / (12.0 * 128 / 8000);
nsmax[1] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(2.0) / (12.0 * 128 / 8000);
nsmax[2] = pow(10.0, db / 10.0 * V * incr / rate);
db = 10.0 * log10(1.2) / (12.0 * 128 / 8000);
nsmax[3] = pow(10.0, db / 10.0 * V * incr / rate);
p.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
alphaOptHat.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
alphaHat.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
sigma2N.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
pbar.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
p2bar.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
Qeq.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
bmin.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
bmin_sub.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
k_mod.resize(msize); // (int *)malloc0(np.msize * sizeof(int));
actmin.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
actmin_sub.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
lmin_flag.resize(msize); // (int *)malloc0(np.msize * sizeof(int));
pmin_u.resize(msize); // (float *)malloc0(np.msize * sizeof(float));
actminbuff.resize(U); // (float**)malloc0(np.U * sizeof(float*));
p.resize(msize);
alphaOptHat.resize(msize);
alphaHat.resize(msize);
sigma2N.resize(msize);
pbar.resize(msize);
p2bar.resize(msize);
Qeq.resize(msize);
bmin.resize(msize);
bmin_sub.resize(msize);
k_mod.resize(msize);
actmin.resize(msize);
actmin_sub.resize(msize);
lmin_flag.resize(msize);
pmin_u.resize(msize);
actminbuff.resize(U);
for (int i = 0; i < U; i++) {
actminbuff[i].resize(msize); // (float *)malloc0(np.msize * sizeof(float));
actminbuff[i].resize(msize);
}
alphaC = 1.0;
subwc = V;
amb_idx = 0;
for (int k = 0; k < msize; k++) {
lambda_y[k] = 0.5;
}
std::copy(lambda_y.begin(), lambda_y.end(), p.begin());
std::copy(lambda_y.begin(), lambda_y.end(), sigma2N.begin());
std::copy(lambda_y.begin(), lambda_y.end(), pbar.begin());
std::copy(lambda_y.begin(), lambda_y.end(), pmin_u.begin());
for (int k = 0; k < msize; k++)
{
int k, ku;
alphaC = 1.0;
subwc = V;
amb_idx = 0;
p2bar[k] = lambda_y[k] * lambda_y[k];
actmin[k] = std::numeric_limits<double>::max();
actmin_sub[k] = std::numeric_limits<double>::max();
for (k = 0; k < msize; k++) {
lambda_y[k] = 0.5;
for (int ku = 0; ku < U; ku++) {
actminbuff[ku][k] = std::numeric_limits<double>::max();
}
std::copy(lambda_y.begin(), lambda_y.end(), p.begin());
std::copy(lambda_y.begin(), lambda_y.end(), sigma2N.begin());
std::copy(lambda_y.begin(), lambda_y.end(), pbar.begin());
std::copy(lambda_y.begin(), lambda_y.end(), pmin_u.begin());
for (k = 0; k < msize; k++)
{
p2bar[k] = lambda_y[k] * lambda_y[k];
actmin[k] = std::numeric_limits<double>::max();
actmin_sub[k] = std::numeric_limits<double>::max();;
for (ku = 0; ku < U; ku++) {
actminbuff[ku][k] = std::numeric_limits<double>::max();;
}
}
std::fill(lmin_flag.begin(), lmin_flag.end(), 0);
}
std::fill(lmin_flag.begin(), lmin_flag.end(), 0);
}
void EMNR::NP::interpM (
@ -249,7 +229,9 @@ void EMNR::NP::interpM (
else
{
int idx = 1;
double xllow, xlhigh, frac;
double xllow;
double xlhigh;
double frac;
while ((x >= xvals[idx]) && (idx < nvals - 1))
idx++;
@ -264,16 +246,23 @@ void EMNR::NP::interpM (
void EMNR::NP::LambdaD()
{
int k;
double f0, f1, f2, f3;
double f0;
double f1;
double f2;
double f3;
double sum_prev_p;
double sum_lambda_y;
double alphaCtilda;
double sum_prev_sigma2N;
double alphaMin, SNR;
double beta, varHat, invQeq;
double alphaMin;
double SNR;
double beta;
double varHat;
double invQeq;
double invQbar;
double bc;
double QeqTilda, QeqTildaSub;
double QeqTilda;
double QeqTildaSub;
double noise_slope_max;
sum_prev_p = 0.0;
@ -369,7 +358,7 @@ void EMNR::NP::LambdaD()
lmin_flag[k] = 0;
actminbuff[amb_idx][k] = actmin[k];
min = std::numeric_limits<double>::max();;
min = std::numeric_limits<double>::max();
for (ku = 0; ku < U; ku++)
{
@ -389,8 +378,8 @@ void EMNR::NP::LambdaD()
}
lmin_flag[k] = 0;
actmin[k] = std::numeric_limits<double>::max();;
actmin_sub[k] = std::numeric_limits<double>::max();;
actmin[k] = std::numeric_limits<double>::max();
actmin_sub[k] = std::numeric_limits<double>::max();
}
if (++amb_idx == U)
@ -433,17 +422,15 @@ EMNR::G::G(
y(_y)
{
lambda_y.resize(msize); // (float *)malloc0(msize * sizeof(float));
lambda_d.resize(msize); // (float *)malloc0(msize * sizeof(float));
prev_gamma.resize(msize); // (float *)malloc0(msize * sizeof(float));
prev_mask.resize(msize); // (float *)malloc0(msize * sizeof(float));
lambda_y.resize(msize);
lambda_d.resize(msize);
prev_gamma.resize(msize);
prev_mask.resize(msize);
gf1p5 = sqrt(PI) / 2.0;
{
double tau = -128.0 / 8000.0 / log(0.98);
alpha = exp(-incr / rate / tau);
}
double tau = -128.0 / 8000.0 / log(0.98);
alpha = exp(-incr / rate / tau);
eps_floor = std::numeric_limits<double>::min();
gamma_max = 1000.0;
@ -480,7 +467,9 @@ EMNR::G::G(
void EMNR::G::calc_gamma0()
{
double gamma, eps_hat, v;
double gamma;
double eps_hat;
double v;
for (int k = 0; k < msize; k++)
{
@ -490,20 +479,15 @@ void EMNR::G::calc_gamma0()
v = (eps_hat / (1.0 + eps_hat)) * gamma;
mask[k] = gf1p5 * sqrt (v) / gamma * exp (- 0.5 * v)
* ((1.0 + v) * bessI0 (0.5 * v) + v * bessI1 (0.5 * v));
{
double v2 = std::min (v, 700.0);
double eta = mask[k] * mask[k] * lambda_y[k] / lambda_d[k];
double eps = eta / (1.0 - q);
double witchHat = (1.0 - q) / q * exp (v2) / (1.0 + eps);
mask[k] *= witchHat / (1.0 + witchHat);
}
double v2 = std::min (v, 700.0);
double eta = mask[k] * mask[k] * lambda_y[k] / lambda_d[k];
double eps = eta / (1.0 - q);
double witchHat = (1.0 - q) / q * exp (v2) / (1.0 + eps);
mask[k] *= witchHat / (1.0 + witchHat);
if (mask[k] > gmax)
mask[k] = gmax;
if (mask[k] != mask[k])
mask[k] = 0.01;
prev_gamma[k] = gamma;
prev_mask[k] = mask[k];
}
@ -511,7 +495,10 @@ void EMNR::G::calc_gamma0()
void EMNR::G::calc_gamma1()
{
double gamma, eps_hat, v, ehr;
double gamma;
double eps_hat;
double v;
double ehr;
for (int k = 0; k < msize; k++)
{
@ -524,9 +511,6 @@ void EMNR::G::calc_gamma1()
if ((mask[k] = ehr * exp (std::min (700.0, 0.5 * e1xb(v)))) > gmax)
mask[k] = gmax;
if (mask[k] != mask[k])
mask[k] = 0.01;
prev_gamma[k] = gamma;
prev_mask[k] = mask[k];
}
@ -534,7 +518,9 @@ void EMNR::G::calc_gamma1()
void EMNR::G::calc_gamma2()
{
double gamma, eps_hat, eps_p;
double gamma;
double eps_hat;
double eps_p;
for (int k = 0; k < msize; k++)
{
@ -572,7 +558,8 @@ void EMNR::G::calc_lambda_y()
double EMNR::G::bessI0 (double x)
{
double res, p;
double res;
double p;
if (x == 0.0)
{
@ -616,7 +603,8 @@ double EMNR::G::bessI0 (double x)
double EMNR::G::bessI1 (double x)
{
double res, p;
double res;
double p;
if (x == 0.0)
{
@ -667,12 +655,17 @@ double EMNR::G::bessI1 (double x)
double EMNR::G::e1xb (double x)
{
double e1, ga, r, t, t0;
int k, m;
double e1;
double ga;
double r;
double t;
double t0;
int k;
int m;
if (x == 0.0)
{
e1 = std::numeric_limits<double>::max();;
e1 = std::numeric_limits<double>::max();
}
else if (x <= 1.0)
{
@ -715,26 +708,25 @@ double EMNR::G::e1xb (double x)
void EMNR::calc_window()
{
int i;
double arg, sum, inv_coherent_gain;
double arg;
double sum;
double inv_coherent_gain;
switch (wintype)
if (wintype == 0)
{
case 0:
arg = 2.0 * PI / (double) fsize;
sum = 0.0;
for (i = 0; i < fsize; i++)
{
window[i] = sqrt (0.54 - 0.46 * cos((float)i * arg));
window[i] = (float) (sqrt (0.54 - 0.46 * cos((float)i * arg)));
sum += window[i];
}
inv_coherent_gain = (double) fsize / sum;
for (i = 0; i < fsize; i++)
window[i] *= inv_coherent_gain;
break;
window[i] *= (float) inv_coherent_gain;
}
}
@ -771,20 +763,20 @@ void EMNR::calc()
init_oainidx = oainidx;
oaoutidx = 0;
msize = fsize / 2 + 1;
window.resize(fsize); // (float *)malloc0(fsize * sizeof(float));
inaccum.resize(iasize); // (float *)malloc0(iasize * sizeof(float));
forfftin.resize(fsize); // (float *)malloc0(fsize * sizeof(float));
forfftout.resize(msize * 2); // (float *)malloc0(msize * sizeof(complex));
mask.resize(msize); // (float *)malloc0(msize * sizeof(float));
window.resize(fsize);
inaccum.resize(iasize);
forfftin.resize(fsize);
forfftout.resize(msize * 2);
mask.resize(msize);
std::fill(mask.begin(), mask.end(), 1.0);
revfftin.resize(msize * 2); // (float *)malloc0(msize * sizeof(complex));
revfftout.resize(fsize); // (float *)malloc0(fsize * sizeof(float));
save.resize(ovrlp); // (float **)malloc0(ovrlp * sizeof(float *));
revfftin.resize(msize * 2);
revfftout.resize(fsize);
save.resize(ovrlp);
for (int i = 0; i < ovrlp; i++)
save[i].resize(fsize); // (float *)malloc0(fsize * sizeof(float));
save[i].resize(fsize);
outaccum.resize(oasize); // (float *)malloc0(oasize * sizeof(float));
outaccum.resize(oasize);
nsamps = 0;
saveidx = 0;
Rfor = fftwf_plan_dft_r2c_1d(
@ -853,10 +845,10 @@ void EMNR::calc()
void EMNR::decalc()
{
delete (ae);
delete (nps);
delete (np);
delete (g);
delete ae;
delete nps;
delete np;
delete g;
fftwf_destroy_plan(Rrev);
fftwf_destroy_plan(Rfor);
@ -896,10 +888,9 @@ EMNR::EMNR(
void EMNR::flush()
{
int i;
std::fill(inaccum.begin(), inaccum.end(), 0);
for (i = 0; i < ovrlp; i++)
for (int i = 0; i < ovrlp; i++)
std::fill(save[i].begin(), save[i].end(), 0);
std::fill(outaccum.begin(), outaccum.end(), 0);
@ -918,9 +909,13 @@ EMNR::~EMNR()
void EMNR::aepf()
{
int k, m;
int N, n;
double sumPre, sumPost, zeta, zetaT;
int k;
int N;
int n;
double sumPre;
double sumPost;
double zeta;
double zetaT;
sumPre = 0.0;
sumPost = 0.0;
@ -947,7 +942,7 @@ void EMNR::aepf()
for (k = n; k < (ae->msize - n); k++)
{
ae->nmask[k] = 0.0;
for (m = k - n; m <= (k + n); m++)
for (int m = k - n; m <= (k + n); m++)
ae->nmask[k] += mask[m];
ae->nmask[k] /= (float)N;
}
@ -957,8 +952,14 @@ void EMNR::aepf()
double EMNR::G::getKey(const std::array<double, 241*241>& type, double gamma, double xi)
{
int ngamma1, ngamma2, nxi1 = 0, nxi2 = 0;
double tg, tx, dg, dx;
int ngamma1;
int ngamma2;
int nxi1 = 0;
int nxi2 = 0;
double tg;
double tx;
double dg;
double dx;
const double dmin = 0.001;
const double dmax = 1000.0;
@ -1005,8 +1006,13 @@ double EMNR::G::getKey(const std::array<double, 241*241>& type, double gamma, do
ix[2] = 241 * nxi1 + ngamma2;
ix[3] = 241 * nxi2 + ngamma2;
for (auto& ixi : ix) {
ixi = ixi < 0 ? 0 : (ixi >= 241*241 ? 241*241 - 1 : ixi);
for (auto& ixi : ix)
{
if (ixi < 0) {
ixi = 0;
} else if (ixi >= 241*241) {
ixi = 241*241 - 1;
}
}
return (1.0 - dg) * (1.0 - dx) * type[ix[0]]
@ -1027,6 +1033,8 @@ void EMNR::calc_gain()
case 1:
nps->LambdaDs();
break;
default:
break;
}
switch (g->gain_method)
@ -1040,6 +1048,8 @@ void EMNR::calc_gain()
case 2:
g->calc_gamma2();
break;
default:
break;
}
if (g->ae_run)
@ -1050,7 +1060,11 @@ void EMNR::execute(int _pos)
{
if (run && _pos == position)
{
int i, j, k, sbuff, sbegin;
int i;
int j;
int k;
int sbuff;
int sbegin;
double g1;
for (i = 0; i < 2 * bsize; i += 2)
@ -1074,8 +1088,8 @@ void EMNR::execute(int _pos)
for (i = 0; i < msize; i++)
{
g1 = gain * mask[i];
revfftin[2 * i + 0] = g1 * forfftout[2 * i + 0];
revfftin[2 * i + 1] = g1 * forfftout[2 * i + 1];
revfftin[2 * i + 0] = (float) (g1 * forfftout[2 * i + 0]);
revfftin[2 * i + 1] = (float) (g1 * forfftout[2 * i + 1]);
}
fftwf_execute (Rrev);

View File

@ -119,7 +119,8 @@ public:
static double e1xb (double x);
static double bessI0 (double x);
static double bessI1 (double x);
} *g;
};
G *g;
struct NP
{
@ -186,7 +187,8 @@ public:
const std::array<double, 18>& xvals,
const std::array<double, 18>& yvals
);
} *np;
};
NP *np;
struct NPS
{
@ -221,8 +223,8 @@ public:
~NPS() = default;
void LambdaDs();
} *nps;
};
NPS *nps;
struct AE
{
int msize;
@ -240,7 +242,8 @@ public:
AE(const AE&) = delete;
AE& operator=(const AE& other) = delete;
~AE() = default;
} *ae;
};
AE *ae;
EMNR(
int run,

View File

@ -42,332 +42,124 @@ namespace WDSP {
********************************************************************************************************/
float* EQ::eq_mults (int size, int nfreqs, float* F, float* G, float samplerate, float scale, int ctfmode, int wintype)
void EQ::eq_mults (std::vector<float>& mults, int size, int nfreqs, float* F, float* G, float samplerate, float scale, int ctfmode, int wintype)
{
float* impulse = EQP::eq_impulse (size + 1, nfreqs, F, G, samplerate, scale, ctfmode, wintype);
float* mults = FIR::fftcv_mults(2 * size, impulse);
delete[] (impulse);
return mults;
float* _mults = FIR::fftcv_mults(2 * size, impulse);
mults.resize(2 * size * 2);
std::copy(_mults, _mults + 2*size*2, mults.begin());
delete[] _mults;
delete[] impulse;
}
void EQ::calc_eq (EQ *a)
void EQ::calc()
{
a->scale = 1.0 / (float)(2 * a->size);
a->infilt = new float[2 * a->size * 2]; // (float *)malloc0(2 * a->size * sizeof(complex));
a->product = new float[2 * a->size * 2]; // (float *)malloc0(2 * a->size * sizeof(complex));
a->CFor = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->infilt, (fftwf_complex *)a->product, FFTW_FORWARD, FFTW_PATIENT);
a->CRev = fftwf_plan_dft_1d(2 * a->size, (fftwf_complex *)a->product, (fftwf_complex *)a->out, FFTW_BACKWARD, FFTW_PATIENT);
a->mults = EQ::eq_mults(a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
scale = (float) (1.0 / (float)(2 * size));
infilt.resize(2 * size * 2);
product.resize(2 * size * 2);
CFor = fftwf_plan_dft_1d(
2 * size,
(fftwf_complex *) infilt.data(),
(fftwf_complex *) product.data(),
FFTW_FORWARD,
FFTW_PATIENT
);
CRev = fftwf_plan_dft_1d(
2 * size,
(fftwf_complex *) product.data(),
(fftwf_complex *) out,
FFTW_BACKWARD,
FFTW_PATIENT
);
EQ::eq_mults(mults, size, nfreqs, F.data(), G.data(), samplerate, scale, ctfmode, wintype);
}
void EQ::decalc_eq (EQ *a)
void EQ::decalc()
{
fftwf_destroy_plan(a->CRev);
fftwf_destroy_plan(a->CFor);
delete[] (a->mults);
delete[] (a->product);
delete[] (a->infilt);
fftwf_destroy_plan(CRev);
fftwf_destroy_plan(CFor);
}
EQ* EQ::create_eq (int run, int size, float *in, float *out, int nfreqs, float* F, float* G, int ctfmode, int wintype, int samplerate)
EQ::EQ(
int _run,
int _size,
float *_in,
float *_out,
int _nfreqs,
const float* _F,
const float* _G,
int _ctfmode,
int _wintype,
int _samplerate
) :
run(_run),
size(_size),
in(_in),
out(_out),
nfreqs(_nfreqs),
ctfmode(_ctfmode),
wintype(_wintype),
samplerate((float) _samplerate)
{
EQ *a = new EQ;
a->run = run;
a->size = size;
a->in = in;
a->out = out;
a->nfreqs = nfreqs;
a->F = new float[a->nfreqs + 1]; // (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = new float[a->nfreqs + 1]; // (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
memcpy (a->F, F, (nfreqs + 1) * sizeof (float));
memcpy (a->G, G, (nfreqs + 1) * sizeof (float));
a->ctfmode = ctfmode;
a->wintype = wintype;
a->samplerate = (float)samplerate;
calc_eq (a);
return a;
F.resize(nfreqs + 1);
G.resize(nfreqs + 1);
std::copy(_F, _F + nfreqs + 1, F.begin());
std::copy(_G, _G + nfreqs + 1, G.begin());
calc();
}
void EQ::destroy_eq (EQ *a)
EQ::~EQ()
{
decalc_eq (a);
delete[] (a->G);
delete[] (a->F);
delete[] (a);
decalc();
}
void EQ::flush_eq (EQ *a)
void EQ::flush()
{
std::fill(a->infilt, a->infilt + 2 * a->size * 2, 0);
std::fill(infilt.begin(), infilt.end(), 0);
}
void EQ::xeq (EQ *a)
void EQ::execute()
{
int i;
float I, Q;
if (a->run)
float I;
float Q;
if (run)
{
std::copy(a->in, a->in + a->size * 2, &(a->infilt[2 * a->size]));
fftwf_execute (a->CFor);
for (i = 0; i < 2 * a->size; i++)
std::copy(in, in + size * 2, &(infilt[2 * size]));
fftwf_execute (CFor);
for (int i = 0; i < 2 * size; i++)
{
I = a->product[2 * i + 0];
Q = a->product[2 * i + 1];
a->product[2 * i + 0] = I * a->mults[2 * i + 0] - Q * a->mults[2 * i + 1];
a->product[2 * i + 1] = I * a->mults[2 * i + 1] + Q * a->mults[2 * i + 0];
I = product[2 * i + 0];
Q = product[2 * i + 1];
product[2 * i + 0] = I * mults[2 * i + 0] - Q * mults[2 * i + 1];
product[2 * i + 1] = I * mults[2 * i + 1] + Q * mults[2 * i + 0];
}
fftwf_execute (a->CRev);
std::copy(&(a->infilt[2 * a->size]), &(a->infilt[2 * a->size]) + a->size * 2, a->infilt);
fftwf_execute (CRev);
std::copy(&(infilt[2 * size]), &(infilt[2 * size]) + size * 2, infilt);
}
else if (a->in != a->out)
std::copy( a->in, a->in + a->size * 2, a->out);
else if (in != out)
std::copy( in, in + size * 2, out);
}
void EQ::setBuffers_eq (EQ *a, float* in, float* out)
void EQ::setBuffers(float* _in, float* _out)
{
decalc_eq (a);
a->in = in;
a->out = out;
calc_eq (a);
decalc();
in = _in;
out = _out;
calc();
}
void EQ::setSamplerate_eq (EQ *a, int rate)
void EQ::setSamplerate(int _rate)
{
decalc_eq (a);
a->samplerate = rate;
calc_eq (a);
decalc();
samplerate = (float) _rate;
calc();
}
void EQ::setSize_eq (EQ *a, int size)
void EQ::setSize(int _size)
{
decalc_eq (a);
a->size = size;
calc_eq (a);
decalc();
size = _size;
calc();
}
/********************************************************************************************************
* *
* Overlap-Save Equalizer: RXA Properties *
* *
********************************************************************************************************/
/* // UNCOMMENT properties when a pointer is in place in rxa
PORT
void SetRXAEQRun (int channel, int run)
{
ch.csDSP.lock();
rxa.eq->run = run;
ch.csDSP.unlock();
}
PORT
void SetRXAEQProfile (int channel, int nfreqs, float* F, float* G)
{
EQ a;
ch.csDSP.lock();
a = rxa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = nfreqs;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
memcpy (a->F, F, (nfreqs + 1) * sizeof (float));
memcpy (a->G, G, (nfreqs + 1) * sizeof (float));
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetRXAEQCtfmode (int channel, int mode)
{
EQ a;
ch.csDSP.lock();
a = rxa.eq;
a->ctfmode = mode;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetRXAEQWintype (int channel, int wintype)
{
EQ a;
ch.csDSP.lock();
a = rxa.eq;
a->wintype = wintype;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetRXAGrphEQ (int channel, int *rxeq)
{ // three band equalizer (legacy compatibility)
EQ a;
ch.csDSP.lock();
a = rxa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = 4;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->F[1] = 150.0;
a->F[2] = 400.0;
a->F[3] = 1500.0;
a->F[4] = 6000.0;
a->G[0] = (float)rxeq[0];
a->G[1] = (float)rxeq[1];
a->G[2] = (float)rxeq[1];
a->G[3] = (float)rxeq[2];
a->G[4] = (float)rxeq[3];
a->ctfmode = 0;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetRXAGrphEQ10 (int channel, int *rxeq)
{ // ten band equalizer (legacy compatibility)
EQ a;
int i;
ch.csDSP.lock();
a = rxa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = 10;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->F[1] = 32.0;
a->F[2] = 63.0;
a->F[3] = 125.0;
a->F[4] = 250.0;
a->F[5] = 500.0;
a->F[6] = 1000.0;
a->F[7] = 2000.0;
a->F[8] = 4000.0;
a->F[9] = 8000.0;
a->F[10] = 16000.0;
for (i = 0; i <= a->nfreqs; i++)
a->G[i] = (float)rxeq[i];
a->ctfmode = 0;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
*/
/********************************************************************************************************
* *
* Overlap-Save Equalizer: TXA Properties *
* *
********************************************************************************************************/
/* // UNCOMMENT properties when a pointer is in place in rxa
PORT
void SetTXAEQRun (int channel, int run)
{
ch.csDSP.lock();
txa.eq->run = run;
ch.csDSP.unlock();
}
PORT
void SetTXAEQProfile (int channel, int nfreqs, float* F, float* G)
{
EQ a;
ch.csDSP.lock();
a = txa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = nfreqs;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
memcpy (a->F, F, (nfreqs + 1) * sizeof (float));
memcpy (a->G, G, (nfreqs + 1) * sizeof (float));
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetTXAEQCtfmode (int channel, int mode)
{
EQ a;
ch.csDSP.lock();
a = txa.eq;
a->ctfmode = mode;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetTXAEQMethod (int channel, int wintype)
{
EQ a;
ch.csDSP.lock();
a = txa.eq;
a->wintype = wintype;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetTXAGrphEQ (int channel, int *txeq)
{ // three band equalizer (legacy compatibility)
EQ a;
ch.csDSP.lock();
a = txa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = 4;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->F[1] = 150.0;
a->F[2] = 400.0;
a->F[3] = 1500.0;
a->F[4] = 6000.0;
a->G[0] = (float)txeq[0];
a->G[1] = (float)txeq[1];
a->G[2] = (float)txeq[1];
a->G[3] = (float)txeq[2];
a->G[4] = (float)txeq[3];
a->ctfmode = 0;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
PORT
void SetTXAGrphEQ10 (int channel, int *txeq)
{ // ten band equalizer (legacy compatibility)
EQ a;
int i;
ch.csDSP.lock();
a = txa.eq;
delete[] (a->G);
delete[] (a->F);
a->nfreqs = 10;
a->F = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->G = (float *) malloc0 ((a->nfreqs + 1) * sizeof (float));
a->F[1] = 32.0;
a->F[2] = 63.0;
a->F[3] = 125.0;
a->F[4] = 250.0;
a->F[5] = 500.0;
a->F[6] = 1000.0;
a->F[7] = 2000.0;
a->F[8] = 4000.0;
a->F[9] = 8000.0;
a->F[10] = 16000.0;
for (i = 0; i <= a->nfreqs; i++)
a->G[i] = (float)txeq[i];
a->ctfmode = 0;
delete[] (a->mults);
a->mults = eq_mults (a->size, a->nfreqs, a->F, a->G, a->samplerate, a->scale, a->ctfmode, a->wintype);
ch.csDSP.unlock();
}
*/
} // namespace WDSP

View File

@ -34,6 +34,8 @@ warren@wpratt.com
#ifndef wdsp_eq_h
#define wdsp_eq_h
#include <vector>
#include "fftw3.h"
#include "export.h"
@ -47,11 +49,11 @@ public:
float* in;
float* out;
int nfreqs;
float* F;
float* G;
float* infilt;
float* product;
float* mults;
std::vector<float> F;
std::vector<float> G;
std::vector<float> infilt;
std::vector<float> product;
std::vector<float> mults;
float scale;
int ctfmode;
int wintype;
@ -59,19 +61,32 @@ public:
fftwf_plan CFor;
fftwf_plan CRev;
static EQ* create_eq (int run, int size, float *in, float *out, int nfreqs, float* F, float* G, int ctfmode, int wintype, int samplerate);
// static float* eq_mults (int size, int nfreqs, float* F, float* G, float samplerate, float scale, int ctfmode, int wintype);
static void destroy_eq (EQ *a);
static void flush_eq (EQ *a);
static void xeq (EQ *a);
static void setBuffers_eq (EQ *a, float* in, float* out);
static void setSamplerate_eq (EQ *a, int rate);
static void setSize_eq (EQ *a, int size);
EQ(
int run,
int size,
float *in,
float *out,
int nfreqs,
const float* F,
const float* G,
int ctfmode,
int wintype,
int samplerate
);
EQ(const EQ&) = delete;
EQ& operator=(EQ& other) = delete;
~EQ() = default;
void flush();
void execute();
void setBuffers(float* in, float* out);
void setSamplerate(int rate);
void setSize(int size);
private:
static float* eq_mults (int size, int nfreqs, float* F, float* G, float samplerate, float scale, int ctfmode, int wintype);
static void calc_eq (EQ *a);
static void decalc_eq (EQ *a);
static void eq_mults (std::vector<float>& mults, int size, int nfreqs, float* F, float* G, float samplerate, float scale, int ctfmode, int wintype);
void calc();
void decalc();
};
} // namespace WDSP

View File

@ -42,22 +42,27 @@ int EQP::fEQcompare (const void * a, const void * b)
return 1;
}
float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate, double scale, int ctfmode, int wintype)
float* EQP::eq_impulse (int N, int nfreqs, const float* F, const float* G, double samplerate, double scale, int ctfmode, int wintype)
{
float* fp = new float[nfreqs + 2]; // (float *) malloc0 ((nfreqs + 2) * sizeof (float));
float* gp = new float[nfreqs + 2]; // (float *) malloc0 ((nfreqs + 2) * sizeof (float));
float* A = new float[N / 2 + 1]; // (float *) malloc0 ((N / 2 + 1) * sizeof (float));
float* sary = new float[2 * nfreqs]; // (float *) malloc0 (2 * nfreqs * sizeof (float));
double gpreamp, f, frac;
std::vector<float> fp(nfreqs + 2);
std::vector<float> gp(nfreqs + 2);
std::vector<float> A(N / 2 + 1);
float* sary = new float[2 * nfreqs];
double gpreamp;
double f;
double frac;
float* impulse;
int i, j, mid;
int i;
int j;
int mid;
fp[0] = 0.0;
fp[nfreqs + 1] = 1.0;
gpreamp = G[0];
for (i = 1; i <= nfreqs; i++)
{
fp[i] = 2.0 * F[i] / samplerate;
fp[i] = (float) (2.0 * F[i] / samplerate);
if (fp[i] < 0.0)
fp[i] = 0.0;
@ -97,7 +102,7 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
j++;
frac = (f - fp[j]) / (fp[j + 1] - fp[j]);
A[i] = pow (10.0, 0.05 * (frac * gp[j + 1] + (1.0 - frac) * gp[j] + gpreamp)) * scale;
A[i] = (float) (pow (10.0, 0.05 * (frac * gp[j + 1] + (1.0 - frac) * gp[j] + gpreamp)) * scale);
}
}
else
@ -110,14 +115,19 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
j++;
frac = (f - fp[j]) / (fp[j + 1] - fp[j]);
A[i] = pow (10.0, 0.05 * (frac * gp[j + 1] + (1.0 - frac) * gp[j] + gpreamp)) * scale;
A[i] = (float) (pow (10.0, 0.05 * (frac * gp[j + 1] + (1.0 - frac) * gp[j] + gpreamp)) * scale);
}
}
if (ctfmode == 0)
{
int k, low, high;
double lowmag, highmag, flow4, fhigh4;
int k;
int low;
int high;
double lowmag;
double highmag;
double flow4;
double fhigh4;
if (N & 1)
{
@ -134,7 +144,7 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
f = (double)k / (double)mid;
lowmag *= (f * f * f * f) / flow4;
if (lowmag < 1.0e-20) lowmag = 1.0e-20;
A[k] = lowmag;
A[k] = (float) lowmag;
}
k = high;
@ -144,7 +154,7 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
f = (double)k / (double)mid;
highmag *= fhigh4 / (f * f * f * f);
if (highmag < 1.0e-20) highmag = 1.0e-20;
A[k] = highmag;
A[k] = (float) highmag;
}
}
else
@ -162,7 +172,7 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
f = (double)k / (double)mid;
lowmag *= (f * f * f * f) / flow4;
if (lowmag < 1.0e-20) lowmag = 1.0e-20;
A[k] = lowmag;
A[k] = (float) lowmag;
}
k = high;
@ -172,21 +182,17 @@ float* EQP::eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate
f = (double)k / (double)mid;
highmag *= fhigh4 / (f * f * f * f);
if (highmag < 1.0e-20) highmag = 1.0e-20;
A[k] = highmag;
A[k] = (float) highmag;
}
}
}
if (N & 1)
impulse = FIR::fir_fsamp_odd(N, A, 1, 1.0, wintype);
impulse = FIR::fir_fsamp_odd(N, A.data(), 1, 1.0, wintype);
else
impulse = FIR::fir_fsamp(N, A, 1, 1.0, wintype);
impulse = FIR::fir_fsamp(N, A.data(), 1, 1.0, wintype);
// print_impulse("eq.txt", N, impulse, 1, 0);
delete[] (sary);
delete[] (A);
delete[] (gp);
delete[] (fp);
delete[] sary;
return impulse;
}
@ -220,8 +226,8 @@ EQP::EQP(
in = _in;
out = _out;
nfreqs = _nfreqs;
F.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
G.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
F.resize(nfreqs + 1);
G.resize(nfreqs + 1);
std::copy(_F, _F + (_nfreqs + 1), F.begin());
std::copy(_G, _G + (_nfreqs + 1), G.begin());
ctfmode = _ctfmode;
@ -229,7 +235,7 @@ EQP::EQP(
samplerate = (double) _samplerate;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
fircore = FIRCORE::create_fircore (size, in, out, nc, mp, impulse);
delete[] (impulse);
delete[] impulse;
}
EQP::~EQP()
@ -263,7 +269,7 @@ void EQP::setSamplerate(int rate)
samplerate = rate;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void EQP::setSize(int _size)
@ -273,7 +279,7 @@ void EQP::setSize(int _size)
FIRCORE::setSize_fircore (fircore, size);
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
/********************************************************************************************************
@ -296,7 +302,7 @@ void EQP::setNC(int _nc)
nc = _nc;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setNc_fircore (fircore, nc, impulse);
delete[] (impulse);
delete[] impulse;
}
}
@ -313,13 +319,13 @@ void EQP::setProfile(int _nfreqs, const float* _F, const float* _G)
{
float* impulse;
nfreqs = _nfreqs;
F.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
G.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
F.resize(nfreqs + 1);
G.resize(nfreqs + 1);
std::copy(_F, _F + (_nfreqs + 1), F.begin());
std::copy(_G, _G + (_nfreqs + 1), G.begin());
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void EQP::setCtfmode(int _mode)
@ -328,7 +334,7 @@ void EQP::setCtfmode(int _mode)
ctfmode = _mode;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void EQP::setWintype(int _wintype)
@ -337,15 +343,15 @@ void EQP::setWintype(int _wintype)
wintype = _wintype;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void EQP::setGrphEQ(int *rxeq)
void EQP::setGrphEQ(const int *rxeq)
{ // three band equalizer (legacy compatibility)
float* impulse;
nfreqs = 4;
F.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
G.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
F.resize(nfreqs + 1);
G.resize(nfreqs + 1);
F[1] = 150.0;
F[2] = 400.0;
F[3] = 1500.0;
@ -358,16 +364,15 @@ void EQP::setGrphEQ(int *rxeq)
ctfmode = 0;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void EQP::setGrphEQ10(int *rxeq)
void EQP::setGrphEQ10(const int *rxeq)
{ // ten band equalizer (legacy compatibility)
float* impulse;
int i;
nfreqs = 10;
F.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
G.resize(nfreqs + 1); // (float *) malloc0 ((nfreqs + 1) * sizeof (float));
F.resize(nfreqs + 1);
G.resize(nfreqs + 1);
F[1] = 32.0;
F[2] = 63.0;
F[3] = 125.0;
@ -378,13 +383,12 @@ void EQP::setGrphEQ10(int *rxeq)
F[8] = 4000.0;
F[9] = 8000.0;
F[10] = 16000.0;
for (i = 0; i <= nfreqs; i++)
for (int i = 0; i <= nfreqs; i++)
G[i] = (float)rxeq[i];
ctfmode = 0;
impulse = eq_impulse (nc, nfreqs, F.data(), G.data(), samplerate, 1.0 / (2.0 * size), ctfmode, wintype);
// print_impulse ("rxeq.txt", nc, impulse, 1, 0);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
} // namespace WDSP

View File

@ -89,10 +89,10 @@ public:
void setProfile(int nfreqs, const float* F, const float* G);
void setCtfmode(int mode);
void setWintype(int wintype);
void setGrphEQ(int *rxeq);
void setGrphEQ10(int *rxeq);
void setGrphEQ(const int *rxeq);
void setGrphEQ10(const int *rxeq);
static float* eq_impulse (int N, int nfreqs, float* F, float* G, double samplerate, double scale, int ctfmode, int wintype);
static float* eq_impulse (int N, int nfreqs, const float* F, const float* G, double samplerate, double scale, int ctfmode, int wintype);
private:
static int fEQcompare (const void * a, const void * b);

View File

@ -27,6 +27,7 @@ warren@pratt.one
#define _CRT_SECURE_NO_WARNINGS
#include <limits>
#include <vector>
#include "fftw3.h"
#include "comm.hpp"
@ -36,132 +37,131 @@ namespace WDSP {
float* FIR::fftcv_mults (int NM, float* c_impulse)
{
float* mults = new float[NM * 2];
float* cfft_impulse = new float[NM * 2];
auto mults = new float[NM * 2];
std::vector<float> cfft_impulse(NM * 2);
fftwf_plan ptmp = fftwf_plan_dft_1d(
NM,
(fftwf_complex *) cfft_impulse,
(fftwf_complex *) cfft_impulse.data(),
(fftwf_complex *) mults,
FFTW_FORWARD,
FFTW_PATIENT
);
std::fill(cfft_impulse, cfft_impulse + NM * 2, 0);
std::fill(cfft_impulse.begin(), cfft_impulse.end(), 0);
// store complex coefs right-justified in the buffer
std::copy(c_impulse, c_impulse + (NM / 2 + 1) * 2, &(cfft_impulse[NM - 2]));
fftwf_execute (ptmp);
fftwf_destroy_plan (ptmp);
delete[] cfft_impulse;
return mults;
}
float* FIR::get_fsamp_window(int N, int wintype)
{
int i;
double arg0, arg1;
float* window = new float[N]; // (float *) malloc0 (N * sizeof(float));
double arg0;
double arg1;
auto window = new float[N];
switch (wintype)
{
case 0:
arg0 = 2.0 * PI / ((double)N - 1.0);
for (i = 0; i < N; i++)
for (int i = 0; i < N; i++)
{
arg1 = cos(arg0 * (double)i);
window[i] = +0.21747
double val = +0.21747
+ arg1 * (-0.45325
+ arg1 * (+0.28256
+ arg1 * (-0.04672)));
window[i] = (float) val;
}
break;
case 1:
arg0 = 2.0 * PI / ((double)N - 1.0);
for (i = 0; i < N; ++i)
for (int i = 0; i < N; ++i)
{
arg1 = cos(arg0 * (double)i);
window[i] = +6.3964424114390378e-02
double val = +6.3964424114390378e-02
+ arg1 * (-2.3993864599352804e-01
+ arg1 * (+3.5015956323820469e-01
+ arg1 * (-2.4774111897080783e-01
+ arg1 * (+8.5438256055858031e-02
+ arg1 * (-1.2320203369293225e-02
+ arg1 * (+4.3778825791773474e-04))))));
window[i] = (float) val;
}
break;
default:
for (i = 0; i < N; i++)
for (int i = 0; i < N; i++)
window[i] = 1.0;
}
return window;
}
float* FIR::fir_fsamp_odd (int N, float* A, int rtype, double scale, int wintype)
float* FIR::fir_fsamp_odd (int N, const float* A, int rtype, double scale, int wintype)
{
int i, j;
int mid = (N - 1) / 2;
double mag, phs;
float* window;
float *fcoef = new float[N * 2];
float *c_impulse = new float[N * 2];
double mag;
double phs;
std::vector<float> fcoef(N * 2);
auto *c_impulse = new float[N * 2];
fftwf_plan ptmp = fftwf_plan_dft_1d(
N,
(fftwf_complex *)fcoef,
(fftwf_complex *)fcoef.data(),
(fftwf_complex *)c_impulse,
FFTW_BACKWARD,
FFTW_PATIENT
);
double local_scale = 1.0 / (double) N;
for (i = 0; i <= mid; i++)
for (int i = 0; i <= mid; i++)
{
mag = A[i] * local_scale;
phs = - (double)mid * TWOPI * (double)i / (double)N;
fcoef[2 * i + 0] = mag * cos (phs);
fcoef[2 * i + 1] = mag * sin (phs);
fcoef[2 * i + 0] = (float) (mag * cos (phs));
fcoef[2 * i + 1] = (float) (mag * sin (phs));
}
for (i = mid + 1, j = 0; i < N; i++, j++)
for (int i = mid + 1, j = 0; i < N; i++, j++)
{
fcoef[2 * i + 0] = + fcoef[2 * (mid - j) + 0];
fcoef[2 * i + 1] = - fcoef[2 * (mid - j) + 1];
}
fftwf_execute (ptmp);
fftwf_destroy_plan (ptmp);
delete[] fcoef;
window = get_fsamp_window(N, wintype);
float* window = get_fsamp_window(N, wintype);
switch (rtype)
{
case 0:
for (i = 0; i < N; i++)
c_impulse[i] = scale * c_impulse[2 * i] * window[i];
for (int i = 0; i < N; i++)
c_impulse[i] = (float) (scale * c_impulse[2 * i] * window[i]);
break;
case 1:
for (i = 0; i < N; i++)
for (int i = 0; i < N; i++)
{
c_impulse[2 * i + 0] *= scale * window[i];
c_impulse[2 * i + 0] *= (float) (scale * window[i]);
c_impulse[2 * i + 1] = 0.0;
}
break;
default:
break;
}
delete[] window;
return c_impulse;
}
float* FIR::fir_fsamp (int N, float* A, int rtype, double scale, int wintype)
float* FIR::fir_fsamp (int N, const float* A, int rtype, double scale, int wintype)
{
int n, i, j, k;
double sum;
float* window;
float *c_impulse = new float[N * 2]; // (float *) malloc0 (N * sizeof (complex));
auto c_impulse = new float[N * 2];
if (N & 1)
{
int M = (N - 1) / 2;
for (n = 0; n < M + 1; n++)
for (int n = 0; n < M + 1; n++)
{
sum = 0.0;
for (k = 1; k < M + 1; k++)
for (int k = 1; k < M + 1; k++)
sum += 2.0 * A[k] * cos(TWOPI * (n - M) * k / N);
c_impulse[2 * n + 0] = (1.0 / N) * (A[0] + sum);
c_impulse[2 * n + 0] = (float) ((1.0 / N) * (A[0] + sum));
c_impulse[2 * n + 1] = 0.0;
}
for (n = M + 1, j = 1; n < N; n++, j++)
for (int n = M + 1, j = 1; n < N; n++, j++)
{
c_impulse[2 * n + 0] = c_impulse[2 * (M - j) + 0];
c_impulse[2 * n + 1] = 0.0;
@ -170,34 +170,36 @@ float* FIR::fir_fsamp (int N, float* A, int rtype, double scale, int wintype)
else
{
double M = (double)(N - 1) / 2.0;
for (n = 0; n < N / 2; n++)
for (int n = 0; n < N / 2; n++)
{
sum = 0.0;
for (k = 1; k < N / 2; k++)
for (int k = 1; k < N / 2; k++)
sum += 2.0 * A[k] * cos(TWOPI * (n - M) * k / N);
c_impulse[2 * n + 0] = (1.0 / N) * (A[0] + sum);
c_impulse[2 * n + 0] = (float) ((1.0 / N) * (A[0] + sum));
c_impulse[2 * n + 1] = 0.0;
}
for (n = N / 2, j = 1; n < N; n++, j++)
for (int n = N / 2, j = 1; n < N; n++, j++)
{
c_impulse[2 * n + 0] = c_impulse[2 * (N / 2 - j) + 0];
c_impulse[2 * n + 1] = 0.0;
}
}
window = get_fsamp_window (N, wintype);
float* window = get_fsamp_window (N, wintype);
switch (rtype)
{
case 0:
for (i = 0; i < N; i++)
c_impulse[i] = scale * c_impulse[2 * i] * window[i];
for (int i = 0; i < N; i++)
c_impulse[i] = (float) (scale * c_impulse[2 * i] * window[i]);
break;
case 1:
for (i = 0; i < N; i++)
for (int i = 0; i < N; i++)
{
c_impulse[2 * i + 0] *= scale * window[i];
c_impulse[2 * i + 0] *= (float) (scale * window[i]);
c_impulse[2 * i + 1] = 0.0;
}
break;
default:
break;
}
delete[] window;
return c_impulse;
@ -205,45 +207,42 @@ float* FIR::fir_fsamp (int N, float* A, int rtype, double scale, int wintype)
float* FIR::fir_bandpass (int N, double f_low, double f_high, double samplerate, int wintype, int rtype, double scale)
{
float *c_impulse = new float[N * 2]; // (float *) malloc0 (N * sizeof (complex));
auto *c_impulse = new float[N * 2];
double ft = (f_high - f_low) / (2.0 * samplerate);
double ft_rad = TWOPI * ft;
double w_osc = PI * (f_high + f_low) / samplerate;
int i, j;
double m = 0.5 * (double)(N - 1);
double delta = PI / m;
double cosphi;
double posi, posj;
double sinc, window, coef;
double posi;
double posj;
double sinc;
double window;
double coef;
if (N & 1)
{
switch (rtype)
{
case 0:
c_impulse[N >> 1] = scale * 2.0 * ft;
c_impulse[N >> 1] = (float) (scale * 2.0 * ft);
break;
case 1:
c_impulse[N - 1] = scale * 2.0 * ft;
c_impulse[N - 1] = (float) (scale * 2.0 * ft);
c_impulse[ N ] = 0.0;
break;
default:
break;
}
}
for (i = (N + 1) / 2, j = N / 2 - 1; i < N; i++, j--)
for (int i = (N + 1) / 2, j = N / 2 - 1; i < N; i++, j--)
{
posi = (double)i - m;
posj = (double)j - m;
sinc = sin (ft_rad * posi) / (PI * posi);
switch (wintype)
if (wintype == 1) // Blackman-Harris 7-term
{
case 0: // Blackman-Harris 4-term
cosphi = cos (delta * i);
window = + 0.21747
+ cosphi * ( - 0.45325
+ cosphi * ( + 0.28256
+ cosphi * ( - 0.04672 )));
break;
case 1: // Blackman-Harris 7-term
cosphi = cos (delta * i);
window = + 6.3964424114390378e-02
+ cosphi * ( - 2.3993864599352804e-01
@ -252,20 +251,31 @@ float* FIR::fir_bandpass (int N, double f_low, double f_high, double samplerate,
+ cosphi * ( + 8.5438256055858031e-02
+ cosphi * ( - 1.2320203369293225e-02
+ cosphi * ( + 4.3778825791773474e-04 ))))));
break;
}
else // Blackman-Harris 4-term
{
cosphi = cos (delta * i);
window = + 0.21747
+ cosphi * ( - 0.45325
+ cosphi * ( + 0.28256
+ cosphi * ( - 0.04672 )));
}
coef = scale * sinc * window;
switch (rtype)
{
case 0:
c_impulse[i] = + coef * cos (posi * w_osc);
c_impulse[j] = + coef * cos (posj * w_osc);
c_impulse[i] = (float) (+ coef * cos (posi * w_osc));
c_impulse[j] = (float) (+ coef * cos (posj * w_osc));
break;
case 1:
c_impulse[2 * i + 0] = + coef * cos (posi * w_osc);
c_impulse[2 * i + 1] = - coef * sin (posi * w_osc);
c_impulse[2 * j + 0] = + coef * cos (posj * w_osc);
c_impulse[2 * j + 1] = - coef * sin (posj * w_osc);
c_impulse[2 * i + 0] = (float) (+ coef * cos (posi * w_osc));
c_impulse[2 * i + 1] = (float) (- coef * sin (posi * w_osc));
c_impulse[2 * j + 0] = (float) (+ coef * cos (posj * w_osc));
c_impulse[2 * j + 1] = (float) (- coef * sin (posj * w_osc));
break;
default:
break;
}
}
@ -282,11 +292,17 @@ float *FIR::fir_read (int N, const char *filename, int rtype, float scale)
// NOTE: The number of values in the file must NOT exceed those implied by N and rtype
{
FILE *file;
int i;
float I, Q;
float *c_impulse = new float[N * 2]; // (float *) malloc0 (N * sizeof (complex));
float I;
float Q;
auto c_impulse = new float[N * 2];
std::fill(c_impulse, c_impulse + N*2, 0);
file = fopen (filename, "r");
for (i = 0; i < N; i++)
if (!file) {
return c_impulse;
}
for (int i = 0; i < N; i++)
{
// read in the complex impulse response
// NOTE: IF the freq response is symmetrical about 0, the imag coeffs will all be zero.
@ -309,6 +325,8 @@ float *FIR::fir_read (int N, const char *filename, int rtype, float scale)
c_impulse[2 * i + 1] = - scale * Q;
break;
}
default:
break;
}
}
fclose (file);
@ -317,77 +335,74 @@ float *FIR::fir_read (int N, const char *filename, int rtype, float scale)
void FIR::analytic (int N, float* in, float* out)
{
if (N < 1) {
if (N < 2) {
return;
}
int i;
double inv_N = 1.0 / (double) N;
double two_inv_N = 2.0 * inv_N;
float* x = new float[N * 2]; // (float *) malloc0 (N * sizeof (complex));
std::vector<float> x(N * 2);
fftwf_plan pfor = fftwf_plan_dft_1d (
N,
(fftwf_complex *) in,
(fftwf_complex *) x,
(fftwf_complex *) x.data(),
FFTW_FORWARD,
FFTW_PATIENT
);
fftwf_plan prev = fftwf_plan_dft_1d (
N,
(fftwf_complex *) x,
(fftwf_complex *) x.data(),
(fftwf_complex *) out,
FFTW_BACKWARD,
FFTW_PATIENT
);
fftwf_execute (pfor);
x[0] *= inv_N;
x[1] *= inv_N;
x[0] *= (float) inv_N;
x[1] *= (float) inv_N;
for (i = 1; i < N / 2; i++)
for (int i = 1; i < N / 2; i++)
{
x[2 * i + 0] *= two_inv_N;
x[2 * i + 1] *= two_inv_N;
x[2 * i + 0] *= (float) two_inv_N;
x[2 * i + 1] *= (float) two_inv_N;
}
x[N + 0] *= inv_N;
x[N + 1] *= inv_N;
x[N + 0] *= (float) inv_N;
x[N + 1] *= (float) inv_N;
memset (&x[N + 2], 0, (N - 2) * sizeof (float));
fftwf_execute (prev);
fftwf_destroy_plan (prev);
fftwf_destroy_plan (pfor);
delete[] x;
}
void FIR::mp_imp (int N, float* fir, float* mpfir, int pfactor, int polarity)
{
int i;
int size = N * pfactor;
double inv_PN = 1.0 / (float)size;
float* firpad = new float[size * 2]; // (float *) malloc0 (size * sizeof (complex));
float* firfreq = new float[size * 2]; // (float *) malloc0 (size * sizeof (complex));
double* mag = new double[size]; // (float *) malloc0 (size * sizeof (float));
float* ana = new float[size * 2]; // (float *) malloc0 (size * sizeof (complex));
float* impulse = new float[size * 2]; // (float *) malloc0 (size * sizeof (complex));
float* newfreq = new float[size * 2]; // (float *) malloc0 (size * sizeof (complex));
std::copy(fir, fir + N * 2, firpad);
double inv_PN = 1.0 / (double)size;
std::vector<float> firpad(size * 2);
std::vector<float> firfreq(size * 2);
std::vector<double> mag(size);
std::vector<float> ana(size * 2);
std::vector<float> impulse(size * 2);
std::vector<float> newfreq(size * 2);
std::copy(fir, fir + N * 2, firpad.begin());
fftwf_plan pfor = fftwf_plan_dft_1d (
size,
(fftwf_complex *) firpad,
(fftwf_complex *) firfreq,
(fftwf_complex *) firpad.data(),
(fftwf_complex *) firfreq.data(),
FFTW_FORWARD,
FFTW_PATIENT);
fftwf_plan prev = fftwf_plan_dft_1d (
size,
(fftwf_complex *) newfreq,
(fftwf_complex *) impulse,
(fftwf_complex *) newfreq.data(),
(fftwf_complex *) impulse.data(),
FFTW_BACKWARD,
FFTW_PATIENT
);
// print_impulse("orig_imp.txt", N, fir, 1, 0);
fftwf_execute (pfor);
for (i = 0; i < size; i++)
{
@ -395,33 +410,27 @@ void FIR::mp_imp (int N, float* fir, float* mpfir, int pfactor, int polarity)
double xi = firfreq[2 * i + 1];
mag[i] = sqrt (xr*xr + xi*xi) * inv_PN;
if (mag[i] > 0.0)
ana[2 * i + 0] = log (mag[i]);
ana[2 * i + 0] = (float) log (mag[i]);
else
ana[2 * i + 0] = log (std::numeric_limits<float>::min());
}
analytic (size, ana, ana);
analytic (size, ana.data(), ana.data());
for (i = 0; i < size; i++)
{
newfreq[2 * i + 0] = + mag[i] * cos (ana[2 * i + 1]);
newfreq[2 * i + 0] = (float) (+ mag[i] * cos (ana[2 * i + 1]));
if (polarity)
newfreq[2 * i + 1] = + mag[i] * sin (ana[2 * i + 1]);
newfreq[2 * i + 1] = (float) (+ mag[i] * sin (ana[2 * i + 1]));
else
newfreq[2 * i + 1] = - mag[i] * sin (ana[2 * i + 1]);
newfreq[2 * i + 1] = (float) (- mag[i] * sin (ana[2 * i + 1]));
}
fftwf_execute (prev);
if (polarity)
std::copy(&impulse[2 * (pfactor - 1) * N], &impulse[2 * (pfactor - 1) * N] + N * 2, mpfir);
else
std::copy(impulse, impulse + N * 2, mpfir);
// print_impulse("min_imp.txt", N, mpfir, 1, 0);
std::copy(impulse.begin(), impulse.end(), mpfir);
fftwf_destroy_plan (prev);
fftwf_destroy_plan (pfor);
delete[] (newfreq);
delete[] (impulse);
delete[] (ana);
delete[] (mag);
delete[] (firfreq);
delete[] (firpad);
}
// impulse response of a zero frequency filter comprising a cascade of two resonators,
@ -432,15 +441,15 @@ float* FIR::zff_impulse(int nc, float scale)
int n_resdet = nc / 2 - 1; // size of single zero-frequency resonator with detrender
int n_dresdet = 2 * n_resdet - 1; // size of two cascaded units; when we convolve these we get 2 * n - 1 length
// allocate the single and make the values
float* resdet = new float[n_resdet]; // (float*)malloc0 (n_resdet * sizeof(float));
std::vector<float> resdet(n_resdet); // (float*)malloc0 (n_resdet * sizeof(float));
for (int i = 1, j = 0, k = n_resdet - 1; i < nc / 4; i++, j++, k--)
resdet[j] = resdet[k] = (float)(i * (i + 1) / 2);
resdet[nc / 4 - 1] = (float)(nc / 4 * (nc / 4 + 1) / 2);
// print_impulse ("resdet", n_resdet, resdet, 0, 0);
// allocate the float and complex versions and make the values
float* dresdet = new float[n_dresdet]; // (float*)malloc0 (n_dresdet * sizeof(float));
float div = (float)((nc / 2 + 1) * (nc / 2 + 1)); // calculate divisor
float* c_dresdet = new float[nc * 2]; // (float*)malloc0 (nc * sizeof(complex));
std::vector<float> dresdet(n_dresdet);
auto div = (float) ((nc / 2 + 1) * (nc / 2 + 1)); // calculate divisor
auto c_dresdet = new float[nc * 2];
for (int n = 0; n < n_dresdet; n++) // convolve to make the cascade
{
for (int k = 0; k < n_resdet; k++)
@ -450,10 +459,7 @@ float* FIR::zff_impulse(int nc, float scale)
c_dresdet[2 * n + 0] = dresdet[n] * scale;
c_dresdet[2 * n + 1] = 0.0;
}
// print_impulse("dresdet", n_dresdet, dresdet, 0, 0);
// print_impulse("c_dresdet", nc, c_dresdet, 1, 0);
delete[] (dresdet);
delete[] (resdet);
return c_dresdet;
}

View File

@ -35,8 +35,8 @@ class WDSP_API FIR
{
public:
static float* fftcv_mults (int NM, float* c_impulse);
static float* fir_fsamp_odd (int N, float* A, int rtype, double scale, int wintype);
static float* fir_fsamp (int N, float* A, int rtype, double scale, int wintype);
static float* fir_fsamp_odd (int N, const float* A, int rtype, double scale, int wintype);
static float* fir_fsamp (int N, const float* A, int rtype, double scale, int wintype);
static float* fir_bandpass (int N, double f_low, double f_high, double samplerate, int wintype, int rtype, double scale);
static void mp_imp (int N, float* fir, float* mpfir, int pfactor, int polarity);

View File

@ -25,6 +25,8 @@ warren@wpratt.com
*/
#include <array>
#include "comm.hpp"
#include "fircore.hpp"
#include "fcurve.hpp"
@ -91,8 +93,8 @@ void FMD::calc()
void FMD::decalc()
{
delete (plim);
delete (sntch);
delete plim;
delete sntch;
}
FMD::FMD(
@ -144,14 +146,24 @@ FMD::FMD(
float* impulse;
calc();
// de-emphasis filter
audio.resize(size * 2); // (float *) malloc0 (size * sizeof (complex));
impulse = FCurve::fc_impulse (nc_de, f_low, f_high, +20.0 * log10(f_high / f_low), 0.0, 1, rate, 1.0 / (2.0 * size), 0, 0);
audio.resize(size * 2);
impulse = FCurve::fc_impulse (
nc_de,
(float) f_low,
(float) f_high,
(float) (+20.0 * log10(f_high / f_low)),
0.0, 1,
(float) rate,
(float) (1.0 / (2.0 * size)),
0,
0
);
pde = FIRCORE::create_fircore (size, audio.data(), out, nc_de, mp_de, impulse);
delete[] (impulse);
delete[] impulse;
// audio filter
impulse = FIR::fir_bandpass(nc_aud, 0.8 * f_low, 1.1 * f_high, rate, 0, 1, afgain / (2.0 * size));
paud = FIRCORE::create_fircore (size, out, out, nc_aud, mp_aud, impulse);
delete[] (impulse);
delete[] impulse;
}
FMD::~FMD()
@ -179,8 +191,11 @@ void FMD::execute()
if (run)
{
int i;
double det, del_out;
double vco[2], corr[2];
double det;
double del_out;
std::array<double, 2> vco;
std::array<double, 2> corr;
for (i = 0; i < size; i++)
{
// pll
@ -200,7 +215,7 @@ void FMD::execute()
while (phs < 0.0) phs += TWOPI;
// dc removal, gain, & demod output
fmdc = mtau * fmdc + onem_mtau * fil_out;
audio[2 * i + 0] = again * (fil_out - fmdc);
audio[2 * i + 0] = (float) (again * (fil_out - fmdc));
audio[2 * i + 1] = audio[2 * i + 0];
}
// de-emphasis
@ -212,7 +227,7 @@ void FMD::execute()
if (lim_run)
{
for (i = 0; i < 2 * size; i++)
out[i] *= lim_pre_gain;
out[i] *= (float) lim_pre_gain;
plim->execute();
}
}
@ -238,13 +253,24 @@ void FMD::setSamplerate(int _rate)
rate = _rate;
calc();
// de-emphasis filter
impulse = FCurve::fc_impulse (nc_de, f_low, f_high, +20.0 * log10(f_high / f_low), 0.0, 1, rate, 1.0 / (2.0 * size), 0, 0);
impulse = FCurve::fc_impulse (
nc_de,
(float) f_low,
(float) f_high,
(float) (+20.0 * log10(f_high / f_low)),
0.0,
1,
(float) rate,
(float) (1.0 / (2.0 * size)),
0,
0
);
FIRCORE::setImpulse_fircore (pde, impulse, 1);
delete[] (impulse);
delete[] impulse;
// audio filter
impulse = FIR::fir_bandpass(nc_aud, 0.8 * f_low, 1.1 * f_high, rate, 0, 1, afgain / (2.0 * size));
FIRCORE::setImpulse_fircore (paud, impulse, 1);
delete[] (impulse);
delete[] impulse;
plim->setSamplerate((int) rate);
}
@ -254,17 +280,28 @@ void FMD::setSize(int _size)
decalc();
size = _size;
calc();
audio.resize(size * 2); // (float *) malloc0 (size * sizeof (complex));
audio.resize(size * 2);
// de-emphasis filter
FIRCORE::destroy_fircore (pde);
impulse = FCurve::fc_impulse (nc_de, f_low, f_high, +20.0 * log10(f_high / f_low), 0.0, 1, rate, 1.0 / (2.0 * size), 0, 0);
impulse = FCurve::fc_impulse (
nc_de,
(float) f_low,
(float) f_high,
(float) (+20.0 * log10(f_high / f_low)),
0.0,
1,
(float) rate,
(float) (1.0 / (2.0 * size)),
0,
0
);
pde = FIRCORE::create_fircore (size, audio.data(), out, nc_de, mp_de, impulse);
delete[] (impulse);
delete[] impulse;
// audio filter
FIRCORE::destroy_fircore (paud);
impulse = FIR::fir_bandpass(nc_aud, 0.8 * f_low, 1.1 * f_high, rate, 0, 1, afgain / (2.0 * size));
paud = FIRCORE::create_fircore (size, out, out, nc_aud, mp_aud, impulse);
delete[] (impulse);
delete[] impulse;
plim->setSize(size);
}
@ -286,9 +323,9 @@ void FMD::setCTCSSFreq(double freq)
sntch->setFreq(ctcss_freq);
}
void FMD::setCTCSSRun(int run)
void FMD::setCTCSSRun(int _run)
{
sntch_run = run;
sntch_run = _run;
sntch->setRun(sntch_run);
}
@ -299,9 +336,20 @@ void FMD::setNCde(int nc)
if (nc_de != nc)
{
nc_de = nc;
impulse = FCurve::fc_impulse (nc_de, f_low, f_high, +20.0 * log10(f_high / f_low), 0.0, 1, rate, 1.0 / (2.0 * size), 0, 0);
impulse = FCurve::fc_impulse (
nc_de,
(float) f_low,
(float) f_high,
(float) (+20.0 * log10(f_high / f_low)),
0.0,
1,
(float) rate,
(float) (1.0 / (2.0 * size)),
0,
0
);
FIRCORE::setNc_fircore (pde, nc_de, impulse);
delete[] (impulse);
delete[] impulse;
}
}
@ -323,7 +371,7 @@ void FMD::setNCaud(int nc)
nc_aud = nc;
impulse = FIR::fir_bandpass(nc_aud, 0.8 * f_low, 1.1 * f_high, rate, 0, 1, afgain / (2.0 * size));
FIRCORE::setNc_fircore (paud, nc_aud, impulse);
delete[] (impulse);
delete[] impulse;
}
}
@ -336,10 +384,10 @@ void FMD::setMPaud(int mp)
}
}
void FMD::setLimRun(int run)
void FMD::setLimRun(int _run)
{
if (lim_run != run) {
lim_run = run;
if (lim_run != _run) {
lim_run = _run;
}
}
@ -364,13 +412,24 @@ void FMD::setAFFilter(double low, double high)
f_low = low;
f_high = high;
// de-emphasis filter
impulse = FCurve::fc_impulse (nc_de, f_low, f_high, +20.0 * log10(f_high / f_low), 0.0, 1, rate, 1.0 / (2.0 * size), 0, 0);
impulse = FCurve::fc_impulse (
nc_de,
(float) f_low,
(float) f_high,
(float) (+20.0 * log10(f_high / f_low)),
0.0,
1,
(float) rate,
(float) (1.0 / (2.0 * size)),
0,
0
);
FIRCORE::setImpulse_fircore (pde, impulse, 1);
delete[] (impulse);
delete[] impulse;
// audio filter
impulse = FIR::fir_bandpass (nc_aud, 0.8 * f_low, 1.1 * f_high, rate, 0, 1, afgain / (2.0 * size));
FIRCORE::setImpulse_fircore (paud, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
}

View File

@ -34,22 +34,23 @@ namespace WDSP {
void FMSQ::calc()
{
double delta, theta;
double delta;
double theta;
float* impulse;
int i;
// noise filter
noise.resize(2 * size * 2); // (float *)malloc0(2 * size * sizeof(complex));
noise.resize(2 * size * 2);
F[0] = 0.0;
F[1] = fc;
F[2] = *pllpole;
F[1] = (float) fc;
F[2] = (float) *pllpole;
F[3] = 20000.0;
G[0] = 0.0;
G[1] = 0.0;
G[2] = 3.0;
G[3] = +20.0 * log10(20000.0 / *pllpole);
G[3] = (float) (+20.0 * log10(20000.0 / *pllpole));
impulse = EQP::eq_impulse (nc, 3, F.data(), G.data(), rate, 1.0 / (2.0 * size), 0, 0);
p = FIRCORE::create_fircore (size, trigger, noise.data(), nc, mp, impulse);
delete[] (impulse);
delete[] impulse;
// noise averaging
avm = exp(-1.0 / (rate * avtau));
onem_avm = 1.0 - avm;
@ -60,8 +61,8 @@ void FMSQ::calc()
// level change
ntup = (int)(tup * rate);
ntdown = (int)(tdown * rate);
cup.resize(ntup + 1); // (float *)malloc0 ((ntup + 1) * sizeof(float));
cdown.resize(ntdown + 1); //(float *)malloc0 ((ntdown + 1) * sizeof(float));
cup.resize(ntup + 1);
cdown.resize(ntdown + 1);
delta = PI / (double) ntup;
theta = 0.0;
@ -80,7 +81,7 @@ void FMSQ::calc()
theta += delta;
}
// control
state = 0;
state = FMSQState::MUTED;
ready = 0;
ramp = 0.0;
rstep = 1.0 / rate;
@ -145,29 +146,20 @@ void FMSQ::flush()
FIRCORE::flush_fircore (p);
avnoise = 100.0;
longnoise = 1.0;
state = 0;
state = FMSQState::MUTED;
ready = 0;
ramp = 0.0;
}
enum _fmsqstate
{
MUTED,
INCREASE,
UNMUTED,
TAIL,
DECREASE
};
void FMSQ::execute()
{
if (run)
{
int i;
double _noise, lnlimit;
double _noise;
double lnlimit;
FIRCORE::xfircore (p);
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
double noise0 = noise[2 * i + 0];
double noise1 = noise[2 * i + 1];
@ -183,10 +175,10 @@ void FMSQ::execute()
switch (state)
{
case MUTED:
case FMSQState::MUTED:
if (avnoise < unmute_thresh && ready)
{
state = INCREASE;
state = FMSQState::INCREASE;
count = ntup;
}
@ -195,19 +187,19 @@ void FMSQ::execute()
break;
case INCREASE:
outsig[2 * i + 0] = insig[2 * i + 0] * cup[ntup - count];
outsig[2 * i + 1] = insig[2 * i + 1] * cup[ntup - count];
case FMSQState::INCREASE:
outsig[2 * i + 0] = (float) (insig[2 * i + 0] * cup[ntup - count]);
outsig[2 * i + 1] = (float) (insig[2 * i + 1] * cup[ntup - count]);
if (count-- == 0)
state = UNMUTED;
state = FMSQState::UNMUTED;
break;
case UNMUTED:
case FMSQState::UNMUTED:
if (avnoise > tail_thresh)
{
state = TAIL;
state = FMSQState::TAIL;
if ((lnlimit = longnoise) > 1.0)
lnlimit = 1.0;
@ -220,28 +212,28 @@ void FMSQ::execute()
break;
case TAIL:
case FMSQState::TAIL:
outsig[2 * i + 0] = insig[2 * i + 0];
outsig[2 * i + 1] = insig[2 * i + 1];
if (avnoise < unmute_thresh)
{
state = UNMUTED;
state = FMSQState::UNMUTED;
}
else if (count-- == 0)
{
state = DECREASE;
state = FMSQState::DECREASE;
count = ntdown;
}
break;
case DECREASE:
outsig[2 * i + 0] = insig[2 * i + 0] * cdown[ntdown - count];
outsig[2 * i + 1] = insig[2 * i + 1] * cdown[ntdown - count];
case FMSQState::DECREASE:
outsig[2 * i + 0] = (float) (insig[2 * i + 0] * cdown[ntdown - count]);
outsig[2 * i + 1] = (float) (insig[2 * i + 1] * cdown[ntdown - count]);
if (count-- == 0)
state = MUTED;
state = FMSQState::MUTED;
break;
}
@ -301,7 +293,7 @@ void FMSQ::setNC(int _nc)
nc = _nc;
impulse = EQP::eq_impulse (nc, 3, F.data(), G.data(), rate, 1.0 / (2.0 * size), 0, 0);
FIRCORE::setNc_fircore (p, nc, impulse);
delete[] (impulse);
delete[] impulse;
}
}

View File

@ -40,6 +40,15 @@ class FIRCORE;
class WDSP_API FMSQ
{
public:
enum class FMSQState
{
MUTED,
INCREASE,
UNMUTED,
TAIL,
DECREASE
};
int run; // 0 if squelch system is OFF; 1 if it's ON
int size; // size of input/output buffers
float* insig; // squelch input signal buffer
@ -59,7 +68,7 @@ public:
double longavm;
double onem_longavm;
double longnoise;
int state; // state machine control
FMSQState state; // state machine control
int count;
double tup;
double tdown;

View File

@ -78,8 +78,8 @@ void GEN::calc_triangle()
void GEN::calc_pulse ()
{
int i;
float delta, theta;
double delta;
double theta;
pulse.pperiod = 1.0 / pulse.pf;
pulse.tphs = 0.0;
pulse.tdelta = TWOPI * pulse.tf / rate;
@ -93,11 +93,11 @@ void GEN::calc_pulse ()
pulse.pnoff = 0;
pulse.pcount = pulse.pnoff;
pulse.state = 0;
pulse.ctrans = new float[pulse.pntrans + 1]; // (float *) malloc0 ((pulse.pntrans + 1) * sizeof (float));
pulse.state = PState::OFF;
pulse.ctrans = new double[pulse.pntrans + 1];
delta = PI / (float)pulse.pntrans;
theta = 0.0;
for (i = 0; i <= pulse.pntrans; i++)
for (int i = 0; i <= pulse.pntrans; i++)
{
pulse.ctrans[i] = 0.5 * (1.0 - cos (theta));
theta += delta;
@ -143,7 +143,7 @@ GEN::GEN(
tt.f1 = + 900.0;
tt.f2 = + 1700.0;
// noise
srand ((unsigned int) time (0));
srand ((unsigned int) time (nullptr));
noise.mag = 1.0;
// sweep
sweep.mag = 1.0;
@ -172,17 +172,9 @@ GEN::~GEN()
void GEN::flush()
{
pulse.state = 0;
pulse.state = PState::OFF;
}
enum pstate
{
OFF,
UP,
ON,
DOWN
};
void GEN::execute()
{
if (run)
@ -191,14 +183,14 @@ void GEN::execute()
{
case 0: // tone
{
int i;
float t1, t2;
float cosphase = cos (tone.phs);
float sinphase = sin (tone.phs);
for (i = 0; i < size; i++)
double t1;
double t2;
double cosphase = cos (tone.phs);
double sinphase = sin (tone.phs);
for (int i = 0; i < size; i++)
{
out[2 * i + 0] = + tone.mag * cosphase;
out[2 * i + 1] = - tone.mag * sinphase;
out[2 * i + 0] = (float) (+ tone.mag * cosphase);
out[2 * i + 1] = (float) (- tone.mag * sinphase);
t1 = cosphase;
t2 = sinphase;
cosphase = t1 * tone.cosdelta - t2 * tone.sindelta;
@ -211,16 +203,16 @@ void GEN::execute()
}
case 1: // two-tone
{
int i;
float tcos, tsin;
float cosphs1 = cos (tt.phs1);
float sinphs1 = sin (tt.phs1);
float cosphs2 = cos (tt.phs2);
float sinphs2 = sin (tt.phs2);
for (i = 0; i < size; i++)
double tcos;
double tsin;
double cosphs1 = cos (tt.phs1);
double sinphs1 = sin (tt.phs1);
double cosphs2 = cos (tt.phs2);
double sinphs2 = sin (tt.phs2);
for (int i = 0; i < size; i++)
{
out[2 * i + 0] = + tt.mag1 * cosphs1 + tt.mag2 * cosphs2;
out[2 * i + 1] = - tt.mag1 * sinphs1 - tt.mag2 * sinphs2;
out[2 * i + 0] = (float) (+ tt.mag1 * cosphs1 + tt.mag2 * cosphs2);
out[2 * i + 1] = (float) (- tt.mag1 * sinphs1 - tt.mag2 * sinphs2);
tcos = cosphs1;
tsin = sinphs1;
cosphs1 = tcos * tt.cosdelta1 - tsin * tt.sindelta1;
@ -240,29 +232,30 @@ void GEN::execute()
}
case 2: // noise
{
int i;
float r1, r2, c, rad;
for (i = 0; i < size; i++)
double r1;
double r2;
double c;
double rad;
for (int i = 0; i < size; i++)
{
do
{
r1 = 2.0 * (float)rand() / (float)RAND_MAX - 1.0;
r2 = 2.0 * (float)rand() / (float)RAND_MAX - 1.0;
r1 = 2.0 * (double)rand() / (double)RAND_MAX - 1.0;
r2 = 2.0 * (double)rand() / (double)RAND_MAX - 1.0;
c = r1 * r1 + r2 * r2;
} while (c >= 1.0);
rad = sqrt (-2.0 * log (c) / c);
out[2 * i + 0] = noise.mag * rad * r1;
out[2 * i + 1] = noise.mag * rad * r2;
out[2 * i + 0] = (float) (noise.mag * rad * r1);
out[2 * i + 1] = (float) (noise.mag * rad * r2);
}
break;
}
case 3: // sweep
{
int i;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
out[2 * i + 0] = + sweep.mag * cos(sweep.phs);
out[2 * i + 1] = - sweep.mag * sin(sweep.phs);
out[2 * i + 0] = (float) (+ sweep.mag * cos(sweep.phs));
out[2 * i + 1] = (float) (- sweep.mag * sin(sweep.phs));
sweep.phs += sweep.dphs;
sweep.dphs += sweep.d2phs;
if (sweep.phs >= TWOPI) sweep.phs -= TWOPI;
@ -274,11 +267,10 @@ void GEN::execute()
}
case 4: // sawtooth (audio only)
{
int i;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
if (saw.t > saw.period) saw.t -= saw.period;
out[2 * i + 0] = saw.mag * (saw.t * saw.f - 1.0);
out[2 * i + 0] = (float) (saw.mag * (saw.t * saw.f - 1.0));
out[2 * i + 1] = 0.0;
saw.t += saw.delta;
}
@ -286,13 +278,12 @@ void GEN::execute()
break;
case 5: // triangle (audio only)
{
int i;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
if (tri.t > tri.period) tri.t1 = tri.t -= tri.period;
if (tri.t > tri.half) tri.t1 -= tri.delta;
else tri.t1 += tri.delta;
out[2 * i + 0] = tri.mag * (4.0 * tri.t1 * tri.f - 1.0);
out[2 * i + 0] = (float) (tri.mag * (4.0 * tri.t1 * tri.f - 1.0));
out[2 * i + 1] = 0.0;
tri.t += tri.delta;
}
@ -300,44 +291,44 @@ void GEN::execute()
break;
case 6: // pulse (audio only)
{
int i;
float t1, t2;
float cosphase = cos (pulse.tphs);
float sinphase = sin (pulse.tphs);
for (i = 0; i < size; i++)
double t1;
double t2;
double cosphase = cos (pulse.tphs);
double sinphase = sin (pulse.tphs);
for (int i = 0; i < size; i++)
{
if (pulse.pnoff != 0)
switch (pulse.state)
{
case OFF:
case PState::OFF:
out[2 * i + 0] = 0.0;
if (--pulse.pcount == 0)
{
pulse.state = UP;
pulse.state = PState::UP;
pulse.pcount = pulse.pntrans;
}
break;
case UP:
out[2 * i + 0] = pulse.mag * cosphase * pulse.ctrans[pulse.pntrans - pulse.pcount];
case PState::UP:
out[2 * i + 0] = (float) (pulse.mag * cosphase * pulse.ctrans[pulse.pntrans - pulse.pcount]);
if (--pulse.pcount == 0)
{
pulse.state = ON;
pulse.state = PState::ON;
pulse.pcount = pulse.pnon;
}
break;
case ON:
out[2 * i + 0] = pulse.mag * cosphase;
case PState::ON:
out[2 * i + 0] = (float) (pulse.mag * cosphase);
if (--pulse.pcount == 0)
{
pulse.state = DOWN;
pulse.state = PState::DOWN;
pulse.pcount = pulse.pntrans;
}
break;
case DOWN:
out[2 * i + 0] = pulse.mag * cosphase * pulse.ctrans[pulse.pcount];
case PState::DOWN:
out[2 * i + 0] = (float) (pulse.mag * cosphase * pulse.ctrans[pulse.pcount]);
if (--pulse.pcount == 0)
{
pulse.state = OFF;
pulse.state = PState::OFF;
pulse.pcount = pulse.pnoff;
}
break;
@ -432,9 +423,9 @@ void GEN::SetPreSweepFreq(float freq1, float freq2)
calc_sweep();
}
void GEN::SetPreSweepRate(float rate)
void GEN::SetPreSweepRate(float _rate)
{
sweep.sweeprate = rate;
sweep.sweeprate = _rate;
calc_sweep();
}

View File

@ -37,88 +37,103 @@ class TXA;
class WDSP_API GEN
{
public:
enum class PState
{
OFF,
UP,
ON,
DOWN
};
int run; // run
int size; // number of samples per buffer
float* in; // input buffer (retained in case I want to mix in a generated signal)
float* out; // output buffer
float rate; // sample rate
double rate; // sample rate
int mode;
struct _tone
{
float mag;
float freq;
float phs;
float delta;
float cosdelta;
float sindelta;
} tone;
double mag;
double freq;
double phs;
double delta;
double cosdelta;
double sindelta;
};
_tone tone;
struct _tt
{
float mag1;
float mag2;
float f1;
float f2;
float phs1;
float phs2;
float delta1;
float delta2;
float cosdelta1;
float cosdelta2;
float sindelta1;
float sindelta2;
} tt;
double mag1;
double mag2;
double f1;
double f2;
double phs1;
double phs2;
double delta1;
double delta2;
double cosdelta1;
double cosdelta2;
double sindelta1;
double sindelta2;
};
_tt tt;
struct _noise
{
float mag;
} noise;
double mag;
};
_noise noise;
struct _sweep
{
float mag;
float f1;
float f2;
float sweeprate;
float phs;
float dphs;
float d2phs;
float dphsmax;
} sweep;
double mag;
double f1;
double f2;
double sweeprate;
double phs;
double dphs;
double d2phs;
double dphsmax;
};
_sweep sweep;
struct _saw
{
float mag;
float f;
float period;
float delta;
float t;
} saw;
double mag;
double f;
double period;
double delta;
double t;
};
_saw saw;
struct _tri
{
float mag;
float f;
float period;
float half;
float delta;
float t;
float t1;
} tri;
double mag;
double f;
double period;
double half;
double delta;
double t;
double t1;
};
_tri tri;
struct _pulse
{
float mag;
float pf;
float pdutycycle;
float ptranstime;
float* ctrans;
double mag;
double pf;
double pdutycycle;
double ptranstime;
double* ctrans;
int pcount;
int pnon;
int pntrans;
int pnoff;
float pperiod;
float tf;
float tphs;
float tdelta;
float tcosdelta;
float tsindelta;
int state;
} pulse;
double pperiod;
double tf;
double tphs;
double tdelta;
double tcosdelta;
double tsindelta;
PState state;
};
_pulse pulse;
GEN(
int run,

View File

@ -51,20 +51,20 @@ METER::METER(
int _enum_pk,
int _enum_gain,
double* _pgain
)
) :
run(_run),
prun(_prun),
size(_size),
buff(_buff),
rate((double) _rate),
tau_average(_tau_av),
tau_peak_decay(_tau_decay),
result(_result),
enum_av(_enum_av),
enum_pk(_enum_pk),
enum_gain(_enum_gain),
pgain(_pgain)
{
run = _run;
prun = _prun;
size = _size;
buff = _buff;
rate = (double) _rate;
tau_average = _tau_av;
tau_peak_decay = _tau_decay;
result = _result;
enum_av = _enum_av;
enum_pk = _enum_pk;
enum_gain = _enum_gain;
pgain = _pgain;
calc();
}
@ -75,7 +75,7 @@ void METER::flush()
result[enum_av] = -100.0;
result[enum_pk] = -100.0;
if ((pgain != 0) && (enum_gain >= 0))
if ((pgain != nullptr) && (enum_gain >= 0))
result[enum_gain] = 0.0;
}
@ -83,18 +83,17 @@ void METER::execute()
{
int srun;
if (prun != 0)
srun = *(prun);
if (prun != nullptr)
srun = *prun;
else
srun = 1;
if (run && srun)
{
int i;
double smag;
double np = 0.0;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
double xr = buff[2 * i + 0];
double xi = buff[2 * i + 1];
@ -112,7 +111,7 @@ void METER::execute()
result[enum_av] = 10.0 * MemLog::mlog10 (avg <= 0 ? 1.0e-20 : avg);
result[enum_pk] = 10.0 * MemLog::mlog10 (peak <= 0 ? 1.0e-20 : peak);
if ((pgain != 0) && (enum_gain >= 0))
if ((pgain != nullptr) && (enum_gain >= 0))
result[enum_gain] = 20.0 * MemLog::mlog10 (*pgain <= 0 ? 1.0e-20 : *pgain);
}
else
@ -150,7 +149,7 @@ void METER::setSize(int _size)
* *
********************************************************************************************************/
double METER::getMeter(int mt)
double METER::getMeter(int mt) const
{
return result[mt];
}

View File

@ -75,7 +75,7 @@ public:
void setBuffers(float* in);
void setSamplerate(int rate);
void setSize(int size);
double getMeter(int mt);
double getMeter(int mt) const;
private:
void calc();

View File

@ -39,8 +39,8 @@ namespace WDSP {
void MPEAK::calc()
{
tmp.resize(size * 2); // (float *) malloc0 (size * sizeof (complex));
mix.resize(size * 2); // (float *) malloc0 (size * sizeof (complex));
tmp.resize(size * 2);
mix.resize(size * 2);
for (int i = 0; i < npeaks; i++)
{
pfil[i] = new SPEAK(
@ -85,15 +85,15 @@ MPEAK::MPEAK(
rate = _rate;
npeaks = _npeaks;
nstages = _nstages;
enable.resize(npeaks); // (int *) malloc0 (npeaks * sizeof (int));
f.resize(npeaks); // (float *) malloc0 (npeaks * sizeof (float));
bw.resize(npeaks); // (float *) malloc0 (npeaks * sizeof (float));
gain.resize(npeaks); // (float *) malloc0 (npeaks * sizeof (float));
enable.resize(npeaks);
f.resize(npeaks);
bw.resize(npeaks);
gain.resize(npeaks);
std::copy(_enable, _enable + npeaks, enable.begin());
std::copy(_f, _f + npeaks, f.begin());
std::copy(_bw, _bw + npeaks, bw.begin());
std::copy(_gain, _gain + npeaks, gain.begin());
pfil.resize(npeaks); // (SPEAK *) malloc0 (npeaks * sizeof (SPEAK));
pfil.resize(npeaks);
calc();
}

View File

@ -25,6 +25,8 @@ warren@wpratt.com
*/
#include <array>
#include "comm.hpp"
#include "fir.hpp"
#include "fircore.hpp"
@ -44,16 +46,17 @@ NOTCHDB::NOTCHDB(int _master_run, int _maxnotches)
master_run = _master_run;
maxnotches = _maxnotches;
nn = 0;
fcenter.resize(maxnotches); // (float *) malloc0 (maxnotches * sizeof (float));
fwidth.resize(maxnotches); // (float *) malloc0 (maxnotches * sizeof (float));
nlow.resize(maxnotches); // (float *) malloc0 (maxnotches * sizeof (float));
nhigh.resize(maxnotches); // (float *) malloc0 (maxnotches * sizeof (float));
active.resize(maxnotches); // (int *) malloc0 (maxnotches * sizeof (int ));
fcenter.resize(maxnotches);
fwidth.resize(maxnotches);
nlow.resize(maxnotches);
nhigh.resize(maxnotches);
active.resize(maxnotches);
}
int NOTCHDB::addNotch(int notch, double _fcenter, double _fwidth, int _active)
{
int i, j;
int i;
int j;
int rval;
if (notch <= nn && nn < maxnotches)
@ -104,7 +107,8 @@ int NOTCHDB::getNotch(int _notch, double* _fcenter, double* _fwidth, int* _activ
int NOTCHDB::deleteNotch(int _notch)
{
int i, j;
int i;
int j;
int rval;
if (_notch < nn)
@ -143,7 +147,7 @@ int NOTCHDB::editNotch(int _notch, double _fcenter, double _fwidth, int _active)
return rval;
}
void NOTCHDB::getNumNotches(int* _nnotches)
void NOTCHDB::getNumNotches(int* _nnotches) const
{
*_nnotches = nn;
}
@ -154,25 +158,24 @@ void NOTCHDB::getNumNotches(int* _nnotches)
* *
********************************************************************************************************/
float* NBP::fir_mbandpass (int N, int nbp, double* flow, double* fhigh, double rate, double scale, int wintype)
float* NBP::fir_mbandpass (int N, int nbp, const double* flow, const double* fhigh, double rate, double scale, int wintype)
{
int i, k;
float* impulse = new float[N * 2]; // (float *) malloc0 (N * sizeof (complex));
float* imp;
for (k = 0; k < nbp; k++)
auto* impulse = new float[N * 2];
std::fill(impulse, impulse + N*2, 0);
for (int k = 0; k < nbp; k++)
{
imp = FIR::fir_bandpass (N, flow[k], fhigh[k], rate, wintype, 1, scale);
for (i = 0; i < N; i++)
float* imp = FIR::fir_bandpass (N, flow[k], fhigh[k], rate, wintype, 1, scale);
for (int i = 0; i < N; i++)
{
impulse[2 * i + 0] += imp[2 * i + 0];
impulse[2 * i + 1] += imp[2 * i + 1];
}
delete[] (imp);
delete[] imp;
}
return impulse;
}
double NBP::min_notch_width()
double NBP::min_notch_width() const
{
double min_width;
switch (wintype)
@ -183,6 +186,8 @@ double NBP::min_notch_width()
case 1:
min_width = 2200.0 / (nc / 256) * (rate / 48000);
break;
default:
min_width = 0;
}
return min_width;
}
@ -204,10 +209,12 @@ int NBP::make_nbp (
)
{
int nbp;
int nnbp, adds;
int i, j, k;
double nl, nh;
int* del = new int[1024]; // (int *) malloc0 (1024 * sizeof (int));
int nnbp;
int adds;
int i;
double nl;
double nh;
std::array<int, 1024> del;
if (fhigh > flow)
{
bplow[0] = flow;
@ -217,11 +224,10 @@ int NBP::make_nbp (
else
{
nbp = 0;
delete[] del;
return nbp;
}
*havnotch = 0;
for (k = 0; k < nn; k++)
for (int k = 0; k < nn; k++)
{
if (autoincr && width[k] < minwidth)
{
@ -270,7 +276,7 @@ int NBP::make_nbp (
if (del[i] == 1)
{
nnbp--;
for (j = i; j < nnbp; j++)
for (int j = i; j < nnbp; j++)
{
bplow[j] = bplow[j + 1];
bphigh[j] = bphigh[j + 1];
@ -281,14 +287,13 @@ int NBP::make_nbp (
nbp = nnbp;
}
}
delete[] (del);
return nbp;
}
void NBP::calc_lightweight()
{ // calculate and set new impulse response; used when changing tune freq or shift freq
int i;
double fl, fh;
double fl;
double fh;
double offset;
NOTCHDB *b = notchdb;
if (fnfrun)
@ -314,7 +319,7 @@ void NBP::calc_lightweight()
// when tuning, no need to recalc filter if there were not and are not any notches in passband
if (hadnotch || havnotch)
{
for (i = 0; i < numpb; i++)
for (int i = 0; i < numpb; i++)
{
bplow[i] -= offset;
bphigh[i] -= offset;
@ -330,7 +335,7 @@ void NBP::calc_lightweight()
);
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
// print_impulse ("nbp.txt", size + 1, impulse, 1, 0);
delete[](impulse);
delete[] impulse;
}
hadnotch = havnotch;
}
@ -340,8 +345,8 @@ void NBP::calc_lightweight()
void NBP::calc_impulse ()
{ // calculates impulse response; for create_fircore() and parameter changes
int i;
float fl, fh;
double fl;
double fh;
double offset;
NOTCHDB *b = notchdb;
@ -365,7 +370,7 @@ void NBP::calc_impulse ()
bphigh,
&havnotch
);
for (i = 0; i < numpb; i++)
for (int i = 0; i < numpb; i++)
{
bplow[i] -= offset;
bphigh[i] -= offset;
@ -429,12 +434,12 @@ NBP::NBP(
maxpb(_maxpb),
notchdb(_notchdb)
{
bplow.resize(maxpb); // (float *) malloc0 (maxpb * sizeof (float));
bphigh.resize(maxpb); // (float *) malloc0 (maxpb * sizeof (float));
bplow.resize(maxpb);
bphigh.resize(maxpb);
calc_impulse ();
fircore = FIRCORE::create_fircore (size, in, out, nc, mp, impulse);
// print_impulse ("nbp.txt", size + 1, impulse, 1, 0);
delete[](impulse);
delete[]impulse;
}
NBP::~NBP()
@ -467,7 +472,7 @@ void NBP::setSamplerate(int _rate)
rate = _rate;
calc_impulse ();
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void NBP::setSize(int _size)
@ -477,14 +482,14 @@ void NBP::setSize(int _size)
FIRCORE::setSize_fircore (fircore, size);
calc_impulse ();
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
void NBP::setNc()
{
calc_impulse();
FIRCORE::setNc_fircore (fircore, nc, impulse);
delete[] (impulse);
delete[] impulse;
}
void NBP::setMp()
@ -513,7 +518,7 @@ void NBP::SetFreqs(double _flow, double _fhigh)
fhigh = _fhigh;
calc_impulse();
FIRCORE::setImpulse_fircore (fircore, impulse, 1);
delete[] (impulse);
delete[] impulse;
}
}
@ -536,7 +541,7 @@ void NBP::SetMP(int _mp)
}
}
void NBP::GetMinNotchWidth(double* minwidth)
void NBP::GetMinNotchWidth(double* minwidth) const
{
*minwidth = min_notch_width();
}

View File

@ -59,7 +59,7 @@ public:
int getNotch (int notch, double* fcenter, double* fwidth, int* active);
int deleteNotch (int notch);
int editNotch (int notch, double fcenter, double fwidth, int active);
void getNumNotches (int* nnotches);
void getNumNotches (int* nnotches) const;
};
@ -125,12 +125,12 @@ public:
void SetFreqs(double flow, double fhigh);
void SetNC(int nc);
void SetMP(int mp);
void GetMinNotchWidth(double* minwidth);
void GetMinNotchWidth(double* minwidth) const;
void calc_lightweight();
private:
static float* fir_mbandpass (int N, int nbp, double* flow, double* fhigh, double rate, double scale, int wintype);
double min_notch_width ();
static float* fir_mbandpass (int N, int nbp, const double* flow, const double* fhigh, double rate, double scale, int wintype);
double min_notch_width () const;
static int make_nbp (
int nn,
std::vector<int>& active,

View File

@ -150,8 +150,10 @@ void NOB::execute()
double mag;
int bf_idx;
int ff_idx;
int lidx, tidx;
int i, j, k;
int lidx;
int tidx;
int j;
int k;
int bfboutidx;
int ffboutidx;
int hcount;
@ -161,7 +163,7 @@ void NOB::execute()
if (run)
{
for (i = 0; i < buffsize; i++)
for (int i = 0; i < buffsize; i++)
{
dline[2 * in_idx + 0] = in[2 * i + 0];
dline[2 * in_idx + 1] = in[2 * i + 1];
@ -189,8 +191,8 @@ void NOB::execute()
{
case 0: // normal output & impulse setup
{
out[2 * i + 0] = dline[2 * out_idx + 0];
out[2 * i + 1] = dline[2 * out_idx + 1];
out[2 * i + 0] = (float) (dline[2 * out_idx + 0]);
out[2 * i + 1] = (float) (dline[2 * out_idx + 1]);
Ilast = dline[2 * out_idx + 0];
Qlast = dline[2 * out_idx + 1];
@ -210,7 +212,6 @@ void NOB::execute()
do
{
len = 0;
hcount = 0;
while ((imp[tidx] > 0 || hcount > 0) && blank_count < max_imp_seq)
@ -314,7 +315,7 @@ void NOB::execute()
switch (mode)
{
case 0: // zero
default: // zero
deltaI = 0.0;
deltaQ = 0.0;
I = 0.0;
@ -367,8 +368,8 @@ void NOB::execute()
case 1: // slew output in advance of blanking period
{
scale = 0.5 + awave[time];
out[2 * i + 0] = Ilast * scale + (1.0 - scale) * I;
out[2 * i + 1] = Qlast * scale + (1.0 - scale) * Q;
out[2 * i + 0] = (float) (Ilast * scale + (1.0 - scale) * I);
out[2 * i + 1] = (float) (Qlast * scale + (1.0 - scale) * Q);
if (++time == adv_slew_count)
{
@ -385,8 +386,8 @@ void NOB::execute()
case 2: // initial advance period
{
out[2 * i + 0] = I;
out[2 * i + 1] = Q;
out[2 * i + 0] = (float) I;
out[2 * i + 1] = (float) Q;
I += deltaI;
Q += deltaQ;
@ -401,8 +402,8 @@ void NOB::execute()
case 3: // impulse & hang period
{
out[2 * i + 0] = I;
out[2 * i + 1] = Q;
out[2 * i + 0] = (float) I;
out[2 * i + 1] = (float) Q;
I += deltaI;
Q += deltaQ;
@ -425,8 +426,8 @@ void NOB::execute()
case 4: // slew output after blanking period
{
scale = 0.5 - hwave[time];
out[2 * i + 0] = Inext * scale + (1.0 - scale) * I;
out[2 * i + 1] = Qnext * scale + (1.0 - scale) * Q;
out[2 * i + 0] = (float) (Inext * scale + (1.0 - scale) * I);
out[2 * i + 1] = (float) (Qnext * scale + (1.0 - scale) * Q);
if (++time == hang_slew_count)
state = 0;
@ -437,8 +438,8 @@ void NOB::execute()
case 5:
{
scale = 0.5 + awave[time];
out[2 * i + 0] = Ilast * scale;
out[2 * i + 1] = Qlast * scale;
out[2 * i + 0] = (float) (Ilast * scale);
out[2 * i + 1] = (float) (Qlast * scale);
if (++time == adv_slew_count)
{
@ -542,8 +543,8 @@ void NOB::execute()
case 9:
{
scale = 0.5 - hwave[time];
out[2 * i + 0] = Inext * scale;
out[2 * i + 1] = Qnext * scale;
out[2 * i + 0] = (float) (Inext * scale);
out[2 * i + 1] = (float) (Qnext * scale);
if (++time >= hang_slew_count)
{

View File

@ -149,7 +149,7 @@ void OSCTRL::SetosctrlRun (TXA& txa, int run)
if (txa.osctrl->run != run)
{
txa.osctrl->run = run;
TXA::SetupBPFilters (txa);
txa.setupBPFilters();
}
}

View File

@ -55,24 +55,26 @@ PANEL::PANEL(
void PANEL::flush()
{
// There is no data to be reset internally
}
void PANEL::execute()
{
int i;
double I, Q;
double I;
double Q;
double gainI = gain1 * gain2I;
double gainQ = gain1 * gain2Q;
// inselect is either 0(neither), 1(Q), 2(I), or 3(both)
switch (copy)
{
case 0: // no copy
default: // 0 (default) no copy
for (i = 0; i < size; i++)
{
I = in[2 * i + 0] * (inselect >> 1);
Q = in[2 * i + 1] * (inselect & 1);
out[2 * i + 0] = gainI * I;
out[2 * i + 1] = gainQ * Q;
out[2 * i + 0] = (float) (gainI * I);
out[2 * i + 1] = (float) (gainQ * Q);
}
break;
case 1: // copy I to Q (then Q == I)
@ -80,17 +82,17 @@ void PANEL::execute()
{
I = in[2 * i + 0] * (inselect >> 1);
Q = I;
out[2 * i + 0] = gainI * I;
out[2 * i + 1] = gainQ * Q;
out[2 * i + 0] = (float) (gainI * I);
out[2 * i + 1] = (float) (gainQ * Q);
}
break;
case 2: // copy Q to I (then I == Q)
for (i = 0; i < size; i++)
{
Q = in[2 * i + 1] * (inselect & 1);
Q = in[2 * i + 1] * (inselect & 1);
I = Q;
out[2 * i + 0] = gainI * I;
out[2 * i + 1] = gainQ * Q;
out[2 * i + 0] = (float) (gainI * I);
out[2 * i + 1] = (float) (gainQ * Q);
}
break;
case 3: // reverse (I=>Q and Q=>I)
@ -98,8 +100,8 @@ void PANEL::execute()
{
Q = in[2 * i + 0] * (inselect >> 1);
I = in[2 * i + 1] * (inselect & 1);
out[2 * i + 0] = gainI * I;
out[2 * i + 1] = gainQ * Q;
out[2 * i + 0] = (float) (gainI * I);
out[2 * i + 1] = (float) (gainQ * Q);
}
break;
}
@ -113,6 +115,7 @@ void PANEL::setBuffers(float* _in, float* _out)
void PANEL::setSamplerate(int)
{
// There is no sample rate to be set for this component
}
void PANEL::setSize(int _size)
@ -149,21 +152,22 @@ void PANEL::setGain2(double _gainI, double _gainQ)
void PANEL::setPan(double _pan)
{
double gain1, gain2;
double _gain1;
double _gain2;
if (_pan <= 0.5)
{
gain1 = 1.0;
gain2 = sin (_pan * PI);
_gain1 = 1.0;
_gain2 = sin (_pan * PI);
}
else
{
gain1 = sin (_pan * PI);
gain2 = 1.0;
_gain1 = sin (_pan * PI);
_gain2 = 1.0;
}
gain2I = gain1;
gain2Q = gain2;
gain2I = _gain1;
gain2Q = _gain2;
}
void PANEL::setCopy(int _copy)

View File

@ -39,10 +39,10 @@ namespace WDSP {
void PHROT::calc()
{
double g;
x0.resize(nstages); // (float *) malloc0 (nstages * sizeof (float));
x1.resize(nstages); // (float *) malloc0 (nstages * sizeof (float));
y0.resize(nstages); // (float *) malloc0 (nstages * sizeof (float));
y1.resize(nstages); // (float *) malloc0 (nstages * sizeof (float));
x0.resize(nstages);
x1.resize(nstages);
y0.resize(nstages);
y1.resize(nstages);
g = tan (PI * fc / (float)rate);
b0 = (g - 1.0) / (g + 1.0);
b1 = 1.0;
@ -58,7 +58,6 @@ PHROT::PHROT(
double _fc,
int _nstages
) :
reverse(0),
run(_run),
size(_size),
in(_in),
@ -67,6 +66,7 @@ PHROT::PHROT(
fc(_fc),
nstages(_nstages)
{
reverse = 0;
calc();
}
@ -102,7 +102,7 @@ void PHROT::execute()
x1[n] = x0[n];
}
out[2 * i + 0] = y0[nstages - 1];
out[2 * i + 0] = (float) y0[nstages - 1];
}
}
else if (out != in)
@ -135,9 +135,9 @@ void PHROT::setSize(int _size)
* *
********************************************************************************************************/
void PHROT::setRun(int run)
void PHROT::setRun(int _run)
{
run = run;
run = _run;
if (run)
flush();

View File

@ -54,8 +54,13 @@ public:
double fc;
int nstages;
// normalized such that a0 = 1
double a1, b0, b1;
std::vector<double> x0, x1, y0, y1;
double a1;
double b0;
double b1;
std::vector<double> x0;
std::vector<double> x1;
std::vector<double> y0;
std::vector<double> y1;
PHROT(
int run,

View File

@ -39,11 +39,14 @@ namespace WDSP {
void RESAMPLE::calc()
{
int x, y, z;
int i, j, k;
int x;
int y;
int z;
int i;
int min_rate;
double full_rate;
double fc_norm_high, fc_norm_low;
double fc_norm_high;
double fc_norm_low;
float* impulse;
fc = fcin;
ncoef = ncoefin;
@ -84,22 +87,22 @@ void RESAMPLE::calc()
ncoef = (ncoef / L + 1) * L;
cpp = ncoef / L;
h.resize(ncoef); // (float *)malloc0(ncoef * sizeof(float));
h.resize(ncoef);
impulse = FIR::fir_bandpass(ncoef, fc_norm_low, fc_norm_high, 1.0, 1, 0, gain * (double)L);
i = 0;
for (j = 0; j < L; j++)
for (int j = 0; j < L; j++)
{
for (k = 0; k < ncoef; k += L)
for (int k = 0; k < ncoef; k += L)
h[i++] = impulse[j + k];
}
ringsize = cpp;
ring.resize(ringsize); // (float *)malloc0(ringsize * sizeof(complex));
ring.resize(ringsize);
idx_in = ringsize - 1;
phnum = 0;
delete[] (impulse);
delete[] impulse;
}
RESAMPLE::RESAMPLE (
@ -141,11 +144,12 @@ int RESAMPLE::execute()
if (run)
{
int i, j, n;
int n;
int idx_out;
double I, Q;
double I;
double Q;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
ring[2 * idx_in + 0] = in[2 * i + 0];
ring[2 * idx_in + 1] = in[2 * i + 1];
@ -156,7 +160,7 @@ int RESAMPLE::execute()
Q = 0.0;
n = cpp * phnum;
for (j = 0; j < cpp; j++)
for (int j = 0; j < cpp; j++)
{
if ((idx_out = idx_in + j) >= ringsize)
idx_out -= ringsize;
@ -165,8 +169,8 @@ int RESAMPLE::execute()
Q += h[n + j] * ring[2 * idx_out + 1];
}
out[2 * outsamps + 0] = I;
out[2 * outsamps + 1] = Q;
out[2 * outsamps + 0] = (float) I;
out[2 * outsamps + 1] = (float) Q;
outsamps++;
phnum += M;
}
@ -231,25 +235,24 @@ void RESAMPLE::setBandwidth(double _fc_low, double _fc_high)
// exported calls
void* RESAMPLE::createV (int in_rate, int out_rate)
RESAMPLE* RESAMPLE::Create(int in_rate, int out_rate)
{
return (void *) new RESAMPLE(1, 0, 0, 0, in_rate, out_rate, 0.0, 0, 1.0);
return new RESAMPLE(1, 0, nullptr, nullptr, in_rate, out_rate, 0.0, 0, 1.0);
}
void RESAMPLE::executeV (float* input, float* output, int numsamps, int* outsamps, void* ptr)
void RESAMPLE::Execute(float* input, float* output, int numsamps, int* outsamps, RESAMPLE* ptr)
{
RESAMPLE *a = (RESAMPLE*) ptr;
a->in = input;
a->out = output;
a->size = numsamps;
*outsamps = a->execute();
ptr->in = input;
ptr->out = output;
ptr->size = numsamps;
*outsamps = ptr->execute();
}
void RESAMPLE::destroyV (void* ptr)
void RESAMPLE::Destroy(RESAMPLE* ptr)
{
delete ( (RESAMPLE*) ptr );
delete ptr;
}
} // namespace WDSP

View File

@ -87,10 +87,10 @@ public:
void setOutRate(int rate);
void setFCLow(double fc_low);
void setBandwidth(double fc_low, double fc_high);
// Exported calls
static void* createV (int in_rate, int out_rate);
static void executeV (float* input, float* output, int numsamps, int* outsamps, void* ptr);
static void destroyV (void* ptr);
// Static methods
static RESAMPLE* Create (int in_rate, int out_rate);
static void Execute (float* input, float* output, int numsamps, int* outsamps, RESAMPLE* ptr);
static void Destroy (RESAMPLE* ptr);
private:
void calc();

View File

@ -36,29 +36,20 @@ SENDER::SENDER(int _run, int _flag, int _mode, int _size, float* _in) :
flag(_flag),
mode(_mode),
size(_size),
in(_in)
in(_in),
spectrumProbe(nullptr)
{
spectrumProbe = nullptr;
}
void SENDER::flush()
{
// There is no internal data to be reset
}
void SENDER::execute()
{
if (run && flag)
{
switch (mode)
{
case 0:
{
if (spectrumProbe) {
spectrumProbe->proceed(in, size);
}
break;
}
}
if (run && flag && (mode == 0) && spectrumProbe) {
spectrumProbe->proceed(in, size);
}
}

View File

@ -50,9 +50,9 @@ SHIFT::SHIFT (
in(_in),
out(_out),
rate((double) _rate),
shift(_fshift)
shift(_fshift),
phase(0.0)
{
phase = 0.0;
calc();
}
@ -65,17 +65,19 @@ void SHIFT::execute()
{
if (run)
{
int i;
double I1, Q1, t1, t2;
double I1;
double Q1;
double t1;
double t2;
double cos_phase = cos (phase);
double sin_phase = sin (phase);
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
I1 = in[2 * i + 0];
Q1 = in[2 * i + 1];
out[2 * i + 0] = I1 * cos_phase - Q1 * sin_phase;
out[2 * i + 1] = I1 * sin_phase + Q1 * cos_phase;
out[2 * i + 0] = (float) (I1 * cos_phase - Q1 * sin_phase);
out[2 * i + 1] = (float) (I1 * sin_phase + Q1 * cos_phase);
t1 = cos_phase;
t2 = sin_phase;
cos_phase = t1 * cos_delta - t2 * sin_delta;

View File

@ -34,23 +34,25 @@ namespace WDSP {
void SIPHON::build_window()
{
int i;
double arg0, cosphi;
double sum, scale;
double arg0;
double cosphi;
double sum;
float scale;
arg0 = 2.0 * PI / ((double) fftsize - 1.0);
sum = 0.0;
for (i = 0; i < fftsize; i++)
{
cosphi = cos (arg0 * (float)i);
window[i] = + 6.3964424114390378e-02
window[i] = (float) (+ 6.3964424114390378e-02
+ cosphi * ( - 2.3993864599352804e-01
+ cosphi * ( + 3.5015956323820469e-01
+ cosphi * ( - 2.4774111897080783e-01
+ cosphi * ( + 8.5438256055858031e-02
+ cosphi * ( - 1.2320203369293225e-02
+ cosphi * ( + 4.3778825791773474e-04 ))))));
+ cosphi * ( + 4.3778825791773474e-04 )))))));
sum += window[i];
}
scale = 1.0 / sum;
scale = 1.0f / (float) sum;
for (i = 0; i < fftsize; i++)
window[i] *= scale;
}
@ -65,23 +67,23 @@ SIPHON::SIPHON(
int _sipsize,
int _fftsize,
int _specmode
)
) :
run(_run),
position(_position),
mode(_mode),
disp(_disp),
insize(_insize),
in(_in),
sipsize(_sipsize), // NOTE: sipsize MUST BE A POWER OF TWO!!
fftsize(_fftsize),
specmode(_specmode)
{
run = _run;
position = _position;
mode = _mode;
disp = _disp;
insize = _insize;
in = _in;
sipsize = _sipsize; // NOTE: sipsize MUST BE A POWER OF TWO!!
fftsize = _fftsize;
specmode = _specmode;
sipbuff.resize(sipsize * 2); // (float *) malloc0 (sipsize * sizeof (complex));
sipbuff.resize(sipsize * 2);
idx = 0;
sipout.resize(sipsize * 2); // (float *) malloc0 (sipsize * sizeof (complex));
specout.resize(fftsize * 2); // (float *) malloc0 (fftsize * sizeof (complex));
sipout.resize(sipsize * 2);
specout.resize(fftsize * 2);
sipplan = fftwf_plan_dft_1d (fftsize, (fftwf_complex *) sipout.data(), (fftwf_complex *) specout.data(), FFTW_FORWARD, FFTW_PATIENT);
window.resize(fftsize * 2); // (float *) malloc0 (fftsize * sizeof (complex));
window.resize(fftsize * 2);
build_window();
}
@ -100,35 +102,28 @@ void SIPHON::flush()
void SIPHON::execute(int pos)
{
int first, second;
int first;
int second;
if (run && position == pos)
if (run && (position == pos) && (mode == 0))
{
switch (mode)
if (insize >= sipsize)
std::copy(&(in[2 * (insize - sipsize)]), &(in[2 * (insize - sipsize)]) + sipsize * 2, sipbuff.begin());
else
{
case 0:
if (insize >= sipsize)
std::copy(&(in[2 * (insize - sipsize)]), &(in[2 * (insize - sipsize)]) + sipsize * 2, sipbuff.begin());
if (insize > (sipsize - idx))
{
first = sipsize - idx;
second = insize - first;
}
else
{
if (insize > (sipsize - idx))
{
first = sipsize - idx;
second = insize - first;
}
else
{
first = insize;
second = 0;
}
std::copy(in, in + first * 2, sipbuff.begin() + 2 * idx);
std::copy(in + 2 * first, in + 2 * first + second * 2, sipbuff.begin());
if ((idx += insize) >= sipsize) idx -= sipsize;
first = insize;
second = 0;
}
break;
case 1:
// Spectrum0 (1, disp, 0, 0, in);
break;
std::copy(in, in + first * 2, sipbuff.begin() + 2 * idx);
std::copy(in + 2 * first, in + 2 * first + second * 2, sipbuff.begin());
if ((idx += insize) >= sipsize) idx -= sipsize;
}
}
}
@ -168,8 +163,7 @@ void SIPHON::suck()
void SIPHON::sip_spectrum()
{
int i;
for (i = 0; i < fftsize; i++)
for (int i = 0; i < fftsize; i++)
{
sipout[2 * i + 0] *= window[i];
sipout[2 * i + 1] *= window[i];
@ -189,7 +183,7 @@ void SIPHON::getaSipF(float* _out, int _size)
suck ();
for (int i = 0; i < _size; i++) {
_out[i] = (float) sipout[2 * i + 0];
_out[i] = sipout[2 * i + 0];
}
}
@ -200,8 +194,8 @@ void SIPHON::getaSipF1(float* _out, int _size)
for (int i = 0; i < _size; i++)
{
_out[2 * i + 0] = (float) sipout[2 * i + 0];
_out[2 * i + 1] = (float) sipout[2 * i + 1];
_out[2 * i + 0] = sipout[2 * i + 0];
_out[2 * i + 1] = sipout[2 * i + 1];
}
}
@ -236,7 +230,11 @@ void SIPHON::setSipSpecmode(int _mode)
void SIPHON::getSpecF1(float* _out)
{ // return spectrum magnitudes in dB
int i, j, mid, m, n;
int i;
int j;
int mid;
int m;
int n;
outsize = fftsize;
suck();
sip_spectrum();
@ -262,68 +260,5 @@ void SIPHON::getSpecF1(float* _out)
}
}
/********************************************************************************************************
* *
* CALLS FOR EXTERNAL USE *
* *
********************************************************************************************************/
/*
#define MAX_EXT_SIPHONS (2) // maximum number of Siphons called from outside wdsp
__declspec (align (16)) SIPHON psiphon[MAX_EXT_SIPHONS]; // array of pointers for Siphons used EXTERNAL to wdsp
PORT
void create_siphonEXT (int id, int run, int insize, int sipsize, int fftsize, int specmode)
{
psiphon[id] = create_siphon (run, 0, 0, 0, insize, 0, sipsize, fftsize, specmode);
}
PORT
void destroy_siphonEXT (int id)
{
destroy_siphon (psiphon[id]);
}
PORT
void flush_siphonEXT (int id)
{
flush_siphon (psiphon[id]);
}
PORT
void xsiphonEXT (int id, float* buff)
{
SIPHON a = psiphon[id];
a->in = buff;
xsiphon (a, 0);
}
PORT
void GetaSipF1EXT (int id, float* out, int size)
{ // return raw samples as floats
SIPHON a = psiphon[id];
int i;
a->update.lock();
a->outsize = size;
suck (a);
a->update.unlock();
for (i = 0; i < size; i++)
{
out[2 * i + 0] = (float)a->sipout[2 * i + 0];
out[2 * i + 1] = (float)a->sipout[2 * i + 1];
}
}
PORT
void SetSiphonInsize (int id, int size)
{
SIPHON a = psiphon[id];
a->update.lock();
a->insize = size;
a->update.unlock();
}
*/
} // namespace WDSP

View File

@ -87,11 +87,6 @@ public:
void setSipDisplay(int disp);
void getSpecF1(float* out);
void setSipSpecmode(int mode);
// Calls for External Use
// static void create_siphonEXT (int id, int run, int insize, int sipsize, int fftsize, int specmode);
// static void destroy_siphonEXT (int id);
// static void xsiphonEXT (int id, float* buff);
// static void SetSiphonInsize (int id, int size);
private:
void build_window();

View File

@ -103,7 +103,7 @@ void USLEW::flush_uslew (USLEW *a)
void USLEW::xuslew (USLEW *a)
{
if (!a->runmode && TXA::UslewCheck (*a->txa))
if (!a->runmode && a->txa->uslewCheck())
a->runmode = 1;
long upslew = *a->ch_upslew;

View File

@ -25,6 +25,8 @@ warren@wpratt.com
*/
#include <array>
#include "comm.hpp"
#include "resample.hpp"
#include "lmath.hpp"
@ -36,10 +38,75 @@ warren@wpratt.com
#include "emnr.hpp"
#include "snba.hpp"
#define MAXIMP 256
namespace WDSP {
SNBA::Exec::Exec(int xsize, int _asize, int _npasses) :
asize(_asize),
npasses(_npasses)
{
a.resize(xsize);
v.resize(xsize);
detout.resize(xsize);
savex.resize(xsize);
xHout.resize(xsize);
unfixed.resize(xsize);
}
void SNBA::Exec::fluxh()
{
std::fill (a.begin(), a.end(), 0);
std::fill (v.begin(), v.end(), 0);
std::fill (detout.begin(), detout.end(), 0);
std::fill (savex.begin(), savex.end(), 0);
std::fill (xHout.begin(), xHout.end(), 0);
std::fill (unfixed.begin(), unfixed.end(), 0);
}
SNBA::Det::Det(
int _xsize,
double _k1,
double _k2,
int _b,
int _pre,
int _post
) :
k1(_k1),
k2(_k2),
b(_b),
pre(_pre),
post(_post)
{
vp.resize(_xsize);
vpwr.resize(_xsize);
}
void SNBA::Det::flush()
{
std::fill(vp.begin(), vp.end(), 0);
std::fill(vpwr.begin(), vpwr.end(), 0);
}
SNBA::Wrk::Wrk(
int xsize,
int asize
) :
xHat_a1rows_max(xsize + asize),
xHat_a2cols_max(xsize + 2 * asize)
{
xHat_r.resize(xsize);
xHat_ATAI.resize(xsize * xsize);
xHat_A1.resize(xHat_a1rows_max * xsize);
xHat_A2.resize(xHat_a1rows_max * xHat_a2cols_max);
xHat_P1.resize(xsize * xHat_a2cols_max);
xHat_P2.resize(xsize);
trI_y.resize(xsize - 1);
trI_v.resize(xsize - 1);
dR_z.resize(xsize - 2);
asolve_r.resize(asize + 1);
asolve_z.resize(asize + 1);
}
void SNBA::calc()
{
if (inrate >= internalrate)
@ -47,8 +114,8 @@ void SNBA::calc()
else
isize = bsize * (internalrate / inrate);
inbuff = new float[isize * 2]; // (double *) malloc0 (isize * sizeof (complex));
outbuff = new float[isize * 2]; // (double *) malloc0 (isize * sizeof (complex));
inbuff = new float[isize * 2];
outbuff = new float[isize * 2];
if (inrate != internalrate)
resamprun = 1;
@ -88,7 +155,7 @@ void SNBA::calc()
iainidx = 0;
iaoutidx = 0;
inaccum = new double[iasize * 2]; // (double *) malloc0 (iasize * sizeof (double));
inaccum.resize(iasize * 2);
nsamps = 0;
if (incr > isize)
@ -105,7 +172,7 @@ void SNBA::calc()
}
init_oaoutidx = oaoutidx;
outaccum = new double[oasize * 2]; // (double *) malloc0 (oasize * sizeof (double));
outaccum.resize(oasize * 2);
}
SNBA::SNBA(
@ -140,15 +207,12 @@ SNBA::SNBA(
iasize(0),
iainidx(0),
iaoutidx(0),
inaccum(nullptr),
xbase(nullptr),
xaux(nullptr),
nsamps(0),
oasize(0),
oainidx(0),
oaoutidx(0),
init_oaoutidx(0),
outaccum(nullptr),
resamprun(0),
isize(0),
inresamp(nullptr),
@ -156,80 +220,29 @@ SNBA::SNBA(
inbuff(nullptr),
outbuff(nullptr),
out_low_cut(_out_low_cut),
out_high_cut(_out_high_cut)
out_high_cut(_out_high_cut),
exec(_xsize, _asize, _npasses),
sdet(_xsize, _k1, _k2, _b, _pre, _post),
wrk(_xsize, _asize)
{
exec.asize = _asize;
exec.npasses = _npasses;
sdet.k1 = _k1;
sdet.k2 = _k2;
sdet.b = _b;
sdet.pre = _pre;
sdet.post = _post;
scan.pmultmin = _pmultmin;
calc();
xbase = new double[2 * xsize]; // (double *) malloc0 (2 * xsize * sizeof (double));
xaux = xbase + xsize;
exec.a = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
exec.v = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
exec.detout = new int[xsize]; //(int *) malloc0 (xsize * sizeof (int));
exec.savex = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
exec.xHout = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
exec.unfixed = new int[xsize]; //(int *) malloc0 (xsize * sizeof (int));
sdet.vp = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
sdet.vpwr = new double[xsize]; //(double *) malloc0 (xsize * sizeof (double));
wrk.xHat_a1rows_max = xsize + exec.asize;
wrk.xHat_a2cols_max = xsize + 2 * exec.asize;
wrk.xHat_r = new double[xsize]; // (double *) malloc0 (xsize * sizeof(double));
wrk.xHat_ATAI = new double[xsize * xsize]; // (double *) malloc0 (xsize * xsize * sizeof(double));
wrk.xHat_A1 = new double[wrk.xHat_a1rows_max * xsize]; // (double *) malloc0 (wrk.xHat_a1rows_max * xsize * sizeof(double));
wrk.xHat_A2 = new double[wrk.xHat_a1rows_max * wrk.xHat_a2cols_max]; // (double *) malloc0 (wrk.xHat_a1rows_max * wrk.xHat_a2cols_max * sizeof(double));
wrk.xHat_P1 = new double[xsize * wrk.xHat_a2cols_max]; // (double *) malloc0 (xsize * wrk.xHat_a2cols_max * sizeof(double));
wrk.xHat_P2 = new double[xsize]; // (double *) malloc0 (xsize * sizeof(double));
wrk.trI_y = new double[xsize - 1]; // (double *) malloc0 ((xsize - 1) * sizeof(double));
wrk.trI_v = new double[xsize - 1]; // (double *) malloc0 ((xsize - 1) * sizeof(double));
wrk.dR_z = new double[xsize - 2]; // (double *) malloc0 ((xsize - 2) * sizeof(double));
wrk.asolve_r = new double[exec.asize + 1]; // (double *) malloc0 ((exec.asize + 1) * sizeof(double));
wrk.asolve_z = new double[exec.asize + 1]; // (double *) malloc0 ((exec.asize + 1) * sizeof(double));
xbase.resize(2 * xsize);
xaux = &xbase[xsize];
}
void SNBA::decalc()
{
delete (outresamp);
delete (inresamp);
delete[] (outbuff);
delete[] (inbuff);
delete[] (outaccum);
delete[] (inaccum);
delete outresamp;
delete inresamp;
delete[] outbuff;
delete[] inbuff;
}
SNBA::~SNBA()
{
delete[] (wrk.xHat_r);
delete[] (wrk.xHat_ATAI);
delete[] (wrk.xHat_A1);
delete[] (wrk.xHat_A2);
delete[] (wrk.xHat_P1);
delete[] (wrk.xHat_P2);
delete[] (wrk.trI_y);
delete[] (wrk.trI_v);
delete[] (wrk.dR_z);
delete[] (wrk.asolve_r);
delete[] (wrk.asolve_z);
delete[] (sdet.vpwr);
delete[] (sdet.vp);
delete[] (exec.unfixed);
delete[] (exec.xHout);
delete[] (exec.savex);
delete[] (exec.detout);
delete[] (exec.v);
delete[] (exec.a);
delete[] (xbase);
decalc();
}
@ -241,20 +254,13 @@ void SNBA::flush()
oainidx = 0;
oaoutidx = init_oaoutidx;
memset (inaccum, 0, iasize * sizeof (double));
memset (outaccum, 0, oasize * sizeof (double));
memset (xaux, 0, xsize * sizeof (double));
memset (exec.a, 0, xsize * sizeof (double));
memset (exec.v, 0, xsize * sizeof (double));
memset (exec.detout, 0, xsize * sizeof (int));
memset (exec.savex, 0, xsize * sizeof (double));
memset (exec.xHout, 0, xsize * sizeof (double));
memset (exec.unfixed, 0, xsize * sizeof (int));
memset (sdet.vp, 0, xsize * sizeof (double));
memset (sdet.vpwr, 0, xsize * sizeof (double));
std::fill(inbuff, inbuff + isize * 2, 0);
std::fill(outbuff, outbuff + isize * 2, 0);
exec.fluxh();
sdet.flush();
std::fill(inaccum.begin(), inaccum.end(), 0);
std::fill(outaccum.begin(), outaccum.end(), 0);
std::fill(xaux, xaux + xsize, 0);
std::fill(inbuff, inbuff + isize * 2, 0);
std::fill(outbuff, outbuff + isize * 2, 0);
inresamp->flush();
outresamp->flush();
@ -282,27 +288,26 @@ void SNBA::setSize(int size)
calc();
}
void SNBA::ATAc0 (int n, int nr, double* A, double* r)
void SNBA::ATAc0 (int n, int nr, std::vector<double>& A, std::vector<double>& r)
{
int i, j;
memset(r, 0, n * sizeof (double));
std::fill(r.begin(), r.begin() + n, 0);
for (i = 0; i < n; i++)
for (int i = 0; i < n; i++)
{
for (j = 0; j < nr; j++)
for (int j = 0; j < nr; j++)
r[i] += A[j * n + i] * A[j * n + 0];
}
}
void SNBA::multA1TA2(double* a1, double* a2, int m, int n, int q, double* c)
void SNBA::multA1TA2(std::vector<double>& a1, std::vector<double>& a2, int m, int n, int q, std::vector<double>& c)
{
int i, j, k;
int k;
int p = q - m;
memset (c, 0, m * n * sizeof (double));
std::fill(c.begin(), c.begin() + m*n, 0);
for (i = 0; i < m; i++)
for (int i = 0; i < m; i++)
{
for (j = 0; j < n; j++)
for (int j = 0; j < n; j++)
{
if (j < p)
{
@ -318,12 +323,12 @@ void SNBA::multA1TA2(double* a1, double* a2, int m, int n, int q, double* c)
}
}
void SNBA::multXKE(double* a, double* xk, int m, int q, int p, double* vout)
void SNBA::multXKE(std::vector<double>& a, const double* xk, int m, int q, int p, std::vector<double>& vout)
{
int i, k;
memset (vout, 0, m * sizeof (double));
int k;
std::fill(vout.begin(), vout.begin() + m, 0);
for (i = 0; i < m; i++)
for (int i = 0; i < m; i++)
{
for (k = i; k < p; k++)
vout[i] += a[i * q + k] * xk[k];
@ -332,14 +337,13 @@ void SNBA::multXKE(double* a, double* xk, int m, int q, int p, double* vout)
}
}
void SNBA::multAv(double* a, double* v, int m, int q, double* vout)
void SNBA::multAv(std::vector<double>& a, std::vector<double>& v, int m, int q, std::vector<double>& vout)
{
int i, k;
memset (vout, 0, m * sizeof (double));
std::fill(vout.begin(), vout.begin() + m, 0);
for (i = 0; i < m; i++)
for (int i = 0; i < m; i++)
{
for (k = 0; k < q; k++)
for (int k = 0; k < q; k++)
vout[i] += a[i * q + k] * v[k];
}
}
@ -347,29 +351,31 @@ void SNBA::multAv(double* a, double* v, int m, int q, double* vout)
void SNBA::xHat(
int xusize,
int asize,
double* xk,
double* a,
double* xout,
double* r,
double* ATAI,
double* A1,
double* A2,
double* P1,
double* P2,
double* trI_y,
double* trI_v,
double* dR_z
const double* xk,
std::vector<double>& a,
std::vector<double>& xout,
std::vector<double>& r,
std::vector<double>& ATAI,
std::vector<double>& A1,
std::vector<double>& A2,
std::vector<double>& P1,
std::vector<double>& P2,
std::vector<double>& trI_y,
std::vector<double>& trI_v,
std::vector<double>& dR_z
)
{
int i, j, k;
int i;
int j;
int k;
int a1rows = xusize + asize;
int a2cols = xusize + 2 * asize;
memset (r, 0, xusize * sizeof(double)); // work space
memset (ATAI, 0, xusize * xusize * sizeof(double)); // work space
memset (A1, 0, a1rows * xusize * sizeof(double)); // work space
memset (A2, 0, a1rows * a2cols * sizeof(double)); // work space
memset (P1, 0, xusize * a2cols * sizeof(double)); // work space
memset (P2, 0, xusize * sizeof(double)); // work space
std::fill (r.begin(), r.begin() + xusize, 0); // work space
std::fill (ATAI.begin(), ATAI.begin() + xusize * xusize, 0); // work space
std::fill (A1.begin(), A1.begin() + a1rows * xusize, 0); // work space
std::fill (A2.begin(), A2.begin() + a1rows * a2cols, 0); // work space
std::fill (P1.begin(), P1.begin() + xusize * a2cols, 0); // work space
std::fill (P2.begin(), P2.begin() + xusize, 0); // work space
for (i = 0; i < xusize; i++)
{
@ -380,28 +386,29 @@ void SNBA::xHat(
}
for (i = 0; i < asize; i++)
{
for (k = asize - i - 1, j = 0; k < asize; k++, j++)
A2[j * a2cols + i] = a[k];
}
{
for (k = asize - i - 1, j = 0; k < asize; k++, j++)
A2[j * a2cols + i] = a[k];
}
for (i = asize + xusize; i < 2 * asize + xusize; i++)
{
A2[(i - asize) * a2cols + i] = - 1.0;
for (j = i - asize + 1, k = 0; j < xusize + asize; j++, k++)
A2[j * a2cols + i] = a[k];
}
{
A2[(i - asize) * a2cols + i] = - 1.0;
for (j = i - asize + 1, k = 0; j < xusize + asize; j++, k++)
A2[j * a2cols + i] = a[k];
}
ATAc0(xusize, xusize + asize, A1, r);
LMathd::trI(xusize, r, ATAI, trI_y, trI_v, dR_z);
LMathd::trI(xusize, r.data(), ATAI.data(), trI_y.data(), trI_v.data(), dR_z.data());
multA1TA2(A1, A2, xusize, 2 * asize + xusize, xusize + asize, P1);
multXKE(P1, xk, xusize, xusize + 2 * asize, asize, P2);
multAv(ATAI, P2, xusize, xusize, xout);
}
void SNBA::invf(int xsize, int asize, double* a, double* x, double* v)
void SNBA::invf(int xsize, int asize, std::vector<double>& a, const double* x, std::vector<double>& v)
{
int i, j;
memset (v, 0, xsize * sizeof (double));
int i;
int j;
std::fill(v.begin(), v.begin() + xsize, 0);
for (i = asize; i < xsize - asize; i++)
{
@ -417,12 +424,16 @@ void SNBA::invf(int xsize, int asize, double* a, double* x, double* v)
}
}
void SNBA::det(int asize, double* v, int* detout)
void SNBA::det(int asize, std::vector<double>& v, std::vector<int>& detout)
{
int i, j;
int i;
int j;
double medpwr;
double t1, t2;
int bstate, bcount, bsamp;
double t1;
double t2;
int bstate;
int bcount;
int bsamp;
for (i = asize, j = 0; i < xsize; i++, j++)
{
@ -430,7 +441,7 @@ void SNBA::det(int asize, double* v, int* detout)
sdet.vp[j] = sdet.vpwr[i];
}
LMathd::median(xsize - asize, sdet.vp, &medpwr);
LMathd::median(xsize - asize, sdet.vp.data(), &medpwr);
t1 = sdet.k1 * medpwr;
t2 = 0.0;
@ -492,6 +503,8 @@ void SNBA::det(int asize, double* v, int* detout)
bstate = 1;
}
break;
default:
break;
}
}
@ -525,24 +538,26 @@ int SNBA::scanFrame(
int xsize,
int pval,
double pmultmin,
int* det,
int* bimp,
int* limp,
int* befimp,
int* aftimp,
int* p_opt,
std::vector<int>& det,
std::array<int, MAXIMP>& bimp,
std::array<int, MAXIMP>& limp,
std::array<int, MAXIMP>& befimp,
std::array<int, MAXIMP>& aftimp,
std::array<int, MAXIMP>& p_opt,
int* next
)
{
int inflag = 0;
int i = 0, j = 0, k = 0;
int i = 0;
int j = 0;
int k = 0;
int nimp = 0;
double td;
int ti;
double merit[MAXIMP] = { 0 };
int nextlist[MAXIMP];
memset (befimp, 0, MAXIMP * sizeof (int));
memset (aftimp, 0, MAXIMP * sizeof (int));
std::array<double, MAXIMP> merit = { 0 };
std::array<int, MAXIMP> nextlist;
std::fill(befimp.begin(), befimp.end(), 0);
std::fill(aftimp.begin(), aftimp.end(), 0);
while (i < xsize && nimp < MAXIMP)
{
@ -555,7 +570,8 @@ int SNBA::scanFrame(
}
else if (det[i] == 1)
{
limp[nimp - 1]++;
if (nimp > 0)
limp[nimp - 1]++;
}
else
{
@ -634,22 +650,20 @@ int SNBA::scanFrame(
void SNBA::execFrame(double* x)
{
int i, k;
int pass;
int nimp;
int bimp[MAXIMP];
int limp[MAXIMP];
int befimp[MAXIMP];
int aftimp[MAXIMP];
int p_opt[MAXIMP];
std::array<int, MAXIMP> bimp;
std::array<int, MAXIMP> limp;
std::array<int, MAXIMP> befimp;
std::array<int, MAXIMP> aftimp;
std::array<int, MAXIMP> p_opt;
int next = 0;
int p;
memcpy (exec.savex, x, xsize * sizeof (double));
LMathd::asolve(xsize, exec.asize, x, exec.a, wrk.asolve_r, wrk.asolve_z);
std::copy(x, x + xsize, exec.savex.begin());
LMathd::asolve(xsize, exec.asize, x, exec.a.data(), wrk.asolve_r.data(), wrk.asolve_z.data());
invf(xsize, exec.asize, exec.a, x, exec.v);
det(exec.asize, exec.v, exec.detout);
for (i = 0; i < xsize; i++)
for (int i = 0; i < xsize; i++)
{
if (exec.detout[i] != 0)
x[i] = 0.0;
@ -657,18 +671,18 @@ void SNBA::execFrame(double* x)
nimp = scanFrame(xsize, exec.asize, scan.pmultmin, exec.detout, bimp, limp, befimp, aftimp, p_opt, &next);
for (pass = 0; pass < exec.npasses; pass++)
for (int pass = 0; pass < exec.npasses; pass++)
{
memcpy (exec.unfixed, exec.detout, xsize * sizeof (int));
std::copy(exec.detout.begin(), exec.detout.end(), exec.unfixed.begin());
for (k = 0; k < nimp; k++)
for (int k = 0; k < nimp; k++)
{
if (k > 0)
scanFrame(xsize, exec.asize, scan.pmultmin, exec.unfixed, bimp, limp, befimp, aftimp, p_opt, &next);
if ((p = p_opt[next]) > 0)
{
LMathd::asolve(xsize, p, x, exec.a, wrk.asolve_r, wrk.asolve_z);
LMathd::asolve(xsize, p, x, exec.a.data(), wrk.asolve_r.data(), wrk.asolve_z.data());
xHat(
limp[next],
p,
@ -685,7 +699,7 @@ void SNBA::execFrame(double* x)
wrk.trI_v,
wrk.dR_z
);
memcpy (&x[bimp[next]], exec.xHout, limp[next] * sizeof (double));
std::copy(exec.xHout.begin(), exec.xHout.begin() + limp[next], &x[bimp[next]]);
memset (&exec.unfixed[bimp[next]], 0, limp[next] * sizeof (int));
}
else
@ -719,12 +733,12 @@ void SNBA::execute()
nsamps -= incr;
memcpy (&outaccum[oainidx], xaux, incr * sizeof (double));
oainidx = (oainidx + incr) % oasize;
memmove (xbase, &xbase[incr], (2 * xsize - incr) * sizeof (double));
std::copy(&xbase[incr], &xbase[incr] + (2 * xsize - incr), xbase.begin());
}
for (i = 0; i < isize; i++)
{
outbuff[2 * i + 0] = outaccum[oaoutidx];
outbuff[2 * i + 0] = (float) outaccum[oaoutidx];
outbuff[2 * i + 1] = 0.0;
oaoutidx = (oaoutidx + 1) % oasize;
}
@ -792,7 +806,8 @@ void SNBA::setPmultmin(double pmultmin)
void SNBA::setOutputBandwidth(double flow, double fhigh)
{
double f_low, f_high;
double f_low = 0;
double f_high = 0;
if (flow >= 0 && fhigh >= 0)
{
@ -802,7 +817,7 @@ void SNBA::setOutputBandwidth(double flow, double fhigh)
if (flow > out_high_cut)
flow = out_high_cut;
f_low = std::max ( out_low_cut, flow);
f_low = std::max (out_low_cut, flow);
f_high = std::min (out_high_cut, fhigh);
}
else if (flow <= 0 && fhigh <= 0)
@ -813,7 +828,7 @@ void SNBA::setOutputBandwidth(double flow, double fhigh)
if (fhigh < -out_high_cut)
fhigh = -out_high_cut;
f_low = std::max ( out_low_cut, -fhigh);
f_low = std::max (out_low_cut, -fhigh);
f_high = std::min (out_high_cut, -flow);
}
else if (flow < 0 && fhigh > 0)

View File

@ -28,11 +28,15 @@ warren@wpratt.com
#ifndef wdsp_snba_h
#define wdsp_snba_h
#include <vector>
#include "export.h"
namespace WDSP{
class RESAMPLE;
class SNBA
class WDSP_API SNBA
{
public:
int run;
@ -47,15 +51,15 @@ public:
int iasize;
int iainidx;
int iaoutidx;
double* inaccum;
double* xbase;
std::vector<double> inaccum;
std::vector<double> xbase;
double* xaux;
int nsamps;
int oasize;
int oainidx;
int oaoutidx;
int init_oaoutidx;
double* outaccum;
std::vector<double> outaccum;
int resamprun;
int isize;
RESAMPLE *inresamp;
@ -64,48 +68,72 @@ public:
float* outbuff;
double out_low_cut;
double out_high_cut;
static const int MAXIMP = 256;
struct _exec
struct Exec
{
int asize;
double* a;
double* v;
int* detout;
double* savex;
double* xHout;
int* unfixed;
std::vector<double> a;
std::vector<double> v;
std::vector<int> detout;
std::vector<double> savex;
std::vector<double> xHout;
std::vector<int> unfixed;
int npasses;
} exec;
struct _det
Exec(int xsize, int _asize, int _npasses);
void fluxh();
};
Exec exec;
struct Det
{
double k1;
double k2;
int b;
int pre;
int post;
double* vp;
double* vpwr;
} sdet;
struct _scan
std::vector<double> vp;
std::vector<double> vpwr;
Det(
int xsize,
double k1,
double k2,
int b,
int pre,
int post
);
void flush();
};
Det sdet;
struct Scan
{
double pmultmin;
} scan;
struct _wrk
};
Scan scan;
struct Wrk
{
int xHat_a1rows_max;
int xHat_a2cols_max;
double* xHat_r;
double* xHat_ATAI;
double* xHat_A1;
double* xHat_A2;
double* xHat_P1;
double* xHat_P2;
double* trI_y;
double* trI_v;
double* dR_z;
double* asolve_r;
double* asolve_z;
} wrk;
std::vector<double> xHat_r;
std::vector<double> xHat_ATAI;
std::vector<double> xHat_A1;
std::vector<double> xHat_A2;
std::vector<double> xHat_P1;
std::vector<double> xHat_P2;
std::vector<double> trI_y;
std::vector<double> trI_v;
std::vector<double> dR_z;
std::vector<double> asolve_r;
std::vector<double> asolve_z;
Wrk(
int xsize,
int asize
);
void flush();
};
Wrk wrk;
SNBA(
int run,
@ -149,40 +177,40 @@ public:
private:
void calc();
void decalc();
static void ATAc0 (int n, int nr, double* A, double* r);
static void multA1TA2(double* a1, double* a2, int m, int n, int q, double* c);
static void multXKE(double* a, double* xk, int m, int q, int p, double* vout);
static void multAv(double* a, double* v, int m, int q, double* vout);
static void ATAc0 (int n, int nr, std::vector<double>& A, std::vector<double>& r);
static void multA1TA2(std::vector<double>& a1, std::vector<double>& a2, int m, int n, int q, std::vector<double>& c);
static void multXKE(std::vector<double>& a, const double* xk, int m, int q, int p, std::vector<double>& vout);
static void multAv(std::vector<double>& a, std::vector<double>& v, int m, int q, std::vector<double>& vout);
static void xHat(
int xusize,
int asize,
double* xk,
double* a,
double* xout,
double* r,
double* ATAI,
double* A1,
double* A2,
double* P1,
double* P2,
double* trI_y,
double* trI_v,
double* dR_z
const double* xk,
std::vector<double>& a,
std::vector<double>& xout,
std::vector<double>& r,
std::vector<double>& ATAI,
std::vector<double>& A1,
std::vector<double>& A2,
std::vector<double>& P1,
std::vector<double>& P2,
std::vector<double>& trI_y,
std::vector<double>& trI_v,
std::vector<double>& dR_z
);
static void invf(int xsize, int asize, double* a, double* x, double* v);
static void invf(int xsize, int asize, std::vector<double>& a, const double* x, std::vector<double>& v);
static int scanFrame(
int xsize,
int pval,
double pmultmin,
int* det,
int* bimp,
int* limp,
int* befimp,
int* aftimp,
int* p_opt,
std::vector<int>& det,
std::array<int, MAXIMP>& bimp,
std::array<int, MAXIMP>& limp,
std::array<int, MAXIMP>& befimp,
std::array<int, MAXIMP>& aftimp,
std::array<int, MAXIMP>& p_opt,
int* next
);
void det(int asize, double* v, int* detout);
void det(int asize, std::vector<double>& v, std::vector<int>& detout);
void execFrame(double* x);
};

View File

@ -38,8 +38,11 @@ namespace WDSP {
void SNOTCH::calc()
{
double fn, qk, qr, csn;
fn = f / (double) rate;
double fn;
double qk;
double qr;
double csn;
fn = f / rate;
csn = cos (TWOPI * fn);
qr = 1.0 - 3.0 * bw;
qk = (1.0 - 2.0 * qr * csn + qr * qr) / (2.0 * (1.0 - csn));
@ -80,11 +83,10 @@ void SNOTCH::execute()
{
if (run)
{
int i;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
x0 = in[2 * i + 0];
out[2 * i + 0] = a0 * x0 + a1 * x1 + a2 * x2 + b1 * y1 + b2 * y2;
out[2 * i + 0] = (float) (a0 * x0 + a1 * x1 + a2 * x2 + b1 * y1 + b2 * y2);
y2 = y1;
y1 = out[2 * i + 0];
x2 = x1;

View File

@ -48,8 +48,16 @@ public:
double rate;
double f;
double bw;
double a0, a1, a2, b1, b2;
double x0, x1, x2, y1, y2;
double a0;
double a1;
double a2;
double b1;
double b2;
double x0;
double x1;
double x2;
double y1;
double y2;
SNOTCH(
int run,
@ -62,7 +70,7 @@ public:
);
SNOTCH(const SNOTCH&) = delete;
SNOTCH& operator=(SNOTCH& other) = delete;
~SNOTCH() {}
~SNOTCH() = default;
void flush();
void execute();

View File

@ -39,31 +39,38 @@ namespace WDSP {
void SPEAK::calc()
{
double ratio;
double f_corr, g_corr, bw_corr, bw_parm, A, f_min;
double f_corr;
double g_corr;
double bw_corr;
double bw_parm;
double A;
double f_min;
switch (design)
{
case 0:
ratio = bw / f;
switch (nstages)
if (nstages == 4)
{
case 4:
bw_parm = 2.4;
f_corr = 1.0 - 0.160 * ratio + 1.440 * ratio * ratio;
g_corr = 1.0 - 1.003 * ratio + 3.990 * ratio * ratio;
break;
default:
}
else
{
bw_parm = 1.0;
f_corr = 1.0;
g_corr = 1.0;
break;
}
{
double fn, qk, qr, csn;
double fn;
double qk;
double qr;
double csn;
fgain = gain / g_corr;
fn = f / (double)rate / f_corr;
fn = f / rate / f_corr;
csn = cos (TWOPI * fn);
qr = 1.0 - 3.0 * bw / (double)rate * bw_parm;
qr = 1.0 - 3.0 * bw / rate * bw_parm;
qk = (1.0 - 2.0 * qr * csn + qr * qr) / (2.0 * (1.0 - csn));
a0 = 1.0 - qk;
a1 = 2.0 * (qk - qr) * csn;
@ -76,26 +83,27 @@ void SPEAK::calc()
case 1:
if (f < 200.0) f = 200.0;
ratio = bw / f;
switch (nstages)
if (nstages == 4)
{
case 4:
bw_parm = 5.0;
bw_corr = 1.13 * ratio - 0.956 * ratio * ratio;
A = 2.5;
f_min = 50.0;
break;
default:
}
else
{
bw_parm = 1.0;
bw_corr = 1.0;
g_corr = 1.0;
A = 2.5;
f_min = 50.0;
break;
}
{
double w0, sn, c, den;
double w0;
double sn;
double c;
double den;
if (f < f_min) f = f_min;
w0 = TWOPI * f / (double)rate;
w0 = TWOPI * f / rate;
sn = sin (w0);
cbw = bw_corr * f;
c = sn * sinh(0.5 * log((f + 0.5 * cbw * bw_parm) / (f - 0.5 * cbw * bw_parm)) * w0 / sn);
@ -108,6 +116,8 @@ void SPEAK::calc()
fgain = gain / pow (A * A, (double)nstages);
}
break;
default:
break;
}
flush();
}
@ -135,12 +145,12 @@ SPEAK::SPEAK(
nstages(_nstages),
design(_design)
{
x0.resize(nstages * 2); // (float *) malloc0 (nstages * sizeof (complex));
x1.resize(nstages * 2); // (float *) malloc0 (nstages * sizeof (complex));
x2.resize(nstages * 2); //(float *) malloc0 (nstages * sizeof (complex));
y0.resize(nstages * 2); // (float *) malloc0 (nstages * sizeof (complex));
y1.resize(nstages * 2); // (float *) malloc0 (nstages * sizeof (complex));
y2.resize(nstages * 2); // (float *) malloc0 (nstages * sizeof (complex));
x0.resize(nstages * 2);
x1.resize(nstages * 2);
x2.resize(nstages * 2);
y0.resize(nstages * 2);
y1.resize(nstages * 2);
y2.resize(nstages * 2);
calc();
}
@ -178,7 +188,7 @@ void SPEAK::execute()
x1[2 * n + j] = x0[2 * n + j];
}
out[2 * i + j] = y0[2 * (nstages - 1) + j];
out[2 * i + j] = (float) y0[2 * (nstages - 1) + j];
}
}
}

View File

@ -55,8 +55,17 @@ public:
double fgain;
int nstages;
int design;
double a0, a1, a2, b1, b2;
std::vector<double> x0, x1, x2, y0, y1, y2;
double a0;
double a1;
double a2;
double b1;
double b2;
std::vector<double> x0;
std::vector<double> x1;
std::vector<double> x2;
std::vector<double> y0;
std::vector<double> y1;
std::vector<double> y2;
SPEAK(
int run,
@ -72,7 +81,7 @@ public:
);
SPEAK(const SPEAK&) = delete;
SPEAK& operator=(const SPEAK& other) = delete;
~SPEAK() {}
~SPEAK() = default;
void flush();
void execute();

View File

@ -81,14 +81,13 @@ void SPHP::execute()
{
if (run)
{
int i, j, n;
for (i = 0; i < size; i++)
for (int i = 0; i < size; i++)
{
for (j = 0; j < 2; j++)
for (int j = 0; j < 2; j++)
{
x0[j] = in[2 * i + j];
for (n = 0; n < nstages; n++)
for (int n = 0; n < nstages; n++)
{
if (n > 0)
x0[2 * n + j] = y0[2 * (n - 1) + j];
@ -100,7 +99,7 @@ void SPHP::execute()
x1[2 * n + j] = x0[2 * n + j];
}
out[2 * i + j] = y0[2 * (nstages - 1) + j];
out[2 * i + j] = (float) y0[2 * (nstages - 1) + j];
}
}
}

View File

@ -56,7 +56,7 @@ FTOV::FTOV(
in = _in;
out = _out;
eps = 0.01;
ring.resize(rsize); // (int*) malloc0 (rsize * sizeof (int));
ring.resize(rsize);
rptr = 0;
inlast = 0.0;
rcount = 0;
@ -91,7 +91,7 @@ void FTOV::execute()
rcount++; // increment the count
}
if (++rptr == rsize) rptr = 0; // increment and wrap the pointer as needed
out[0] = std::min (1.0, (double)rcount / div); // calculate the output sample
out[0] = (float) std::min (1.0, (double)rcount / div); // calculate the output sample
inlast = in[size - 1]; // save the last input sample for next buffer
for (int i = 1; i < size; i++)
{
@ -107,7 +107,7 @@ void FTOV::execute()
rcount++; // increment the count
}
if (++rptr == rsize) rptr = 0; // increment and wrap the pointer as needed
out[i] = std::min(1.0, (double)rcount / div); // calculate the output sample
out[i] = (float) std::min(1.0, (double)rcount / div); // calculate the output sample
}
}
}
@ -118,7 +118,8 @@ void FTOV::execute()
void SSQL::compute_slews()
{
double delta, theta;
double delta;
double theta;
delta = PI / (double) ntup;
theta = 0.0;
for (int i = 0; i <= ntup; i++)
@ -137,15 +138,33 @@ void SSQL::compute_slews()
void SSQL::calc()
{
b1 = new float[size * 2]; // (float*) malloc0 (size * sizeof (complex));
dcbl = new CBL(1, size, in, b1, 0, rate, 0.02);
ibuff = new float[size]; // (float*) malloc0 (size * sizeof (float));
ftovbuff = new float[size]; // (float*) malloc0(size * sizeof (float));
cvtr = new FTOV(1, size, rate, ftov_rsize, ftov_fmax, ibuff, ftovbuff);
lpbuff = new float[size]; // (float*) malloc0 (size * sizeof (float));
filt = new DBQLP(1, size, ftovbuff, lpbuff, rate, 11.3, 1.0, 1.0, 1);
wdbuff = new int[size]; // (int*) malloc0 (size * sizeof (int));
tr_signal = new int[size]; // (int*) malloc0 (size * sizeof (int));
b1.resize(size * 2);
dcbl = new CBL(1, size, in, b1.data(), 0, rate, 0.02);
ibuff.resize(size);
ftovbuff.resize(size);
cvtr = new FTOV(
1,
size,
rate,
ftov_rsize,
ftov_fmax,
ibuff.data(),
ftovbuff.data()
);
lpbuff.resize(size);
filt = new DBQLP(
1,
size,
ftovbuff.data(),
lpbuff.data(),
rate,
11.3,
1.0,
1.0,
1
);
wdbuff.resize(size);
tr_signal.resize(size);
// window detector
wdmult = exp (-1.0 / (rate * wdtau));
wdaverage = 0.0;
@ -156,27 +175,19 @@ void SSQL::calc()
// level change
ntup = (int)(tup * rate);
ntdown = (int)(tdown * rate);
cup = new float[ntup + 1]; // (float*) malloc0 ((ntup + 1) * sizeof (float));
cdown = new float[ntdown + 1]; // (float*) malloc0 ((ntdown + 1) * sizeof (float));
cup.resize(ntup + 1);
cdown.resize(ntdown + 1);
compute_slews();
// control
state = 0;
state = SSQLState::MUTED;
count = 0;
}
void SSQL::decalc()
{
delete[] (tr_signal);
delete[] (wdbuff);
delete (filt);
delete[] (lpbuff);
delete (cvtr);
delete[] (ftovbuff);
delete[] (ibuff);
delete (dcbl);
delete[] (b1);
delete[] (cdown);
delete[] (cup);
delete filt;
delete cvtr;
delete dcbl;
}
SSQL::SSQL(
@ -223,24 +234,17 @@ SSQL::~SSQL()
void SSQL::flush()
{
std::fill(b1, b1 + size * 2, 0);
std::fill(b1.begin(), b1.end(), 0);
dcbl->flush();
memset (ibuff, 0, size * sizeof (float));
memset (ftovbuff, 0, size * sizeof (float));
std::fill(ibuff.begin(), ibuff.end(), 0);
std::fill(ftovbuff.begin(), ftovbuff.end(), 0);
cvtr->flush();
memset (lpbuff, 0, size * sizeof (float));
std::fill(lpbuff.begin(), lpbuff.end(), 0);
filt->flush();
memset (wdbuff, 0, size * sizeof (int));
memset (tr_signal, 0, size * sizeof (int));
std::fill(wdbuff.begin(), wdbuff.end(), 0);
std::fill(tr_signal.begin(), tr_signal.end(), 0);
}
enum _ssqlstate
{
MUTED,
INCREASE,
UNMUTED,
DECREASE
};
void SSQL::execute()
{
@ -277,35 +281,35 @@ void SSQL::execute()
{
switch (state)
{
case MUTED:
case SSQLState::MUTED:
if (tr_signal[i] == 1)
{
state = INCREASE;
state = SSQLState::INCREASE;
count = ntup;
}
out[2 * i + 0] = muted_gain * in[2 * i + 0];
out[2 * i + 1] = muted_gain * in[2 * i + 1];
out[2 * i + 0] = (float) (muted_gain * in[2 * i + 0]);
out[2 * i + 1] = (float) (muted_gain * in[2 * i + 1]);
break;
case INCREASE:
out[2 * i + 0] = in[2 * i + 0] * cup[ntup - count];
out[2 * i + 1] = in[2 * i + 1] * cup[ntup - count];
case SSQLState::INCREASE:
out[2 * i + 0] = (float) (in[2 * i + 0] * cup[ntup - count]);
out[2 * i + 1] = (float) (in[2 * i + 1] * cup[ntup - count]);
if (count-- == 0)
state = UNMUTED;
state = SSQLState::UNMUTED;
break;
case UNMUTED:
case SSQLState::UNMUTED:
if (tr_signal[i] == 0)
{
state = DECREASE;
state = SSQLState::DECREASE;
count = ntdown;
}
out[2 * i + 0] = in[2 * i + 0];
out[2 * i + 1] = in[2 * i + 1];
break;
case DECREASE:
out[2 * i + 0] = in[2 * i + 0] * cdown[ntdown - count];
out[2 * i + 1] = in[2 * i + 1] * cdown[ntdown - count];
case SSQLState::DECREASE:
out[2 * i + 0] = (float) (in[2 * i + 0] * cdown[ntdown - count]);
out[2 * i + 1] = (float) (in[2 * i + 1] * cdown[ntdown - count]);
if (count-- == 0)
state = MUTED;
state = SSQLState::MUTED;
break;
}
}

View File

@ -75,26 +75,33 @@ class DBQLP;
class WDSP_API SSQL // Syllabic Squelch
{
public:
enum class SSQLState
{
MUTED,
INCREASE,
UNMUTED,
DECREASE
};
int run; // 0 if squelch system is OFF; 1 if it's ON
int size; // size of input/output buffers
float* in; // squelch input signal buffer
float* out; // squelch output signal buffer
int rate; // sample rate
int state; // state machine control
SSQLState state; // state machine control
int count; // count variable for raised cosine transitions
double tup; // time for turn-on transition
double tdown; // time for turn-off transition
int ntup; // number of samples for turn-on transition
int ntdown; // number of samples for turn-off transition
float* cup; // coefficients for up-slew
float* cdown; // coefficients for down-slew
std::vector<double> cup; // coefficients for up-slew
std::vector<double> cdown; // coefficients for down-slew
double muted_gain; // audio gain while muted; 0.0 for complete silence
float* b1; // buffer to hold output of dc-block function
float* ibuff; // buffer containing only 'I' component
float* ftovbuff; // buffer containing output of f to v converter
float* lpbuff; // buffer containing output of low-pass filter
int* wdbuff; // buffer containing output of window detector
std::vector<float> b1; // buffer to hold output of dc-block function
std::vector<float> ibuff; // buffer containing only 'I' component
std::vector<float> ftovbuff; // buffer containing output of f to v converter
std::vector<float> lpbuff; // buffer containing output of low-pass filter
std::vector<int> wdbuff; // buffer containing output of window detector
CBL *dcbl; // pointer to DC Blocker data structure
FTOV *cvtr; // pointer to F to V Converter data structure
DBQLP *filt; // pointer to Bi-Quad Low-Pass Filter data structure
@ -114,7 +121,7 @@ public:
double tr_voltage; // trigger voltage
double mute_mult; // multiplier for successive voltage calcs when muted
double unmute_mult; // multiplier for successive voltage calcs when unmuted
int* tr_signal; // trigger signal, 0 or 1
std::vector<int> tr_signal; // trigger signal, 0 or 1
SSQL(
int run,

148
wdsp/unit.cpp Normal file
View File

@ -0,0 +1,148 @@
/* unit.hpp
This file is part of a program that implements a Software-Defined Radio.
Copyright (C) 2013 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 <algorithm>
#include "unit.hpp"
namespace WDSP {
Unit::Unit(
int _in_rate, // input samplerate
int _out_rate, // output samplerate
int _dsp_rate, // sample rate for mainstream dsp processing
int _dsp_size // number complex samples processed per buffer in mainstream dsp processing
) :
in_rate{_in_rate},
out_rate(_out_rate),
dsp_rate(_dsp_rate),
dsp_size(_dsp_size)
{
if (_in_rate >= _dsp_rate)
dsp_insize = _dsp_size * (_in_rate / _dsp_rate);
else
dsp_insize = _dsp_size / (_dsp_rate / _in_rate);
if (_out_rate >= _dsp_rate)
dsp_outsize = _dsp_size * (_out_rate / _dsp_rate);
else
dsp_outsize = _dsp_size / (_dsp_rate / _out_rate);
// buffers
inbuff = new float[1 * dsp_insize * 2];
outbuff = new float[1 * dsp_outsize * 2];
midbuff = new float[2 * dsp_size * 2];
}
Unit::~Unit()
{
delete[] inbuff;
delete[] outbuff;
delete[] midbuff;
}
void Unit::flushBuffers()
{
std::fill(inbuff, inbuff + 1 * dsp_insize * 2, 0);
std::fill(outbuff, outbuff + 1 * dsp_outsize * 2, 0);
std::fill(midbuff, midbuff + 2 * dsp_size * 2, 0);
}
void Unit::setBuffersInputSamplerate(int _in_rate)
{
if (_in_rate >= dsp_rate)
dsp_insize = dsp_size * (_in_rate / dsp_rate);
else
dsp_insize = dsp_size / (dsp_rate / _in_rate);
in_rate = _in_rate;
// buffers
delete[] (inbuff);
inbuff = new float[1 * dsp_insize * 2];
}
void Unit::setBuffersOutputSamplerate(int _out_rate)
{
if (_out_rate >= dsp_rate)
dsp_outsize = dsp_size * (_out_rate / dsp_rate);
else
dsp_outsize = dsp_size / (dsp_rate / _out_rate);
out_rate = _out_rate;
// buffers
delete[] outbuff;
outbuff = new float[1 * dsp_outsize * 2];
}
void Unit::setBuffersDSPSamplerate(int _dsp_rate)
{
if (in_rate >= _dsp_rate)
dsp_insize = dsp_size * (in_rate / _dsp_rate);
else
dsp_insize = dsp_size / (_dsp_rate / in_rate);
if (out_rate >= _dsp_rate)
dsp_outsize = dsp_size * (out_rate / _dsp_rate);
else
dsp_outsize = dsp_size / (_dsp_rate / out_rate);
dsp_rate = _dsp_rate;
// buffers
delete[] inbuff;
inbuff = new float[1 * dsp_insize * 2];
delete[] outbuff;
outbuff = new float[1 * dsp_outsize * 2];
}
void Unit::setBuffersDSPBuffsize(int _dsp_size)
{
if (in_rate >= dsp_rate)
dsp_insize = _dsp_size * (in_rate / dsp_rate);
else
dsp_insize = _dsp_size / (dsp_rate / in_rate);
if (out_rate >= dsp_rate)
dsp_outsize = _dsp_size * (out_rate / dsp_rate);
else
dsp_outsize = _dsp_size / (dsp_rate / out_rate);
dsp_size = _dsp_size;
// buffers
delete[]inbuff;
inbuff = new float[1 * dsp_insize * 2];
delete[] midbuff;
midbuff = new float[2 * dsp_size * 2];
delete[] outbuff;
outbuff = new float[1 * dsp_outsize * 2];
}
} // namespace

View File

@ -42,6 +42,25 @@ public:
int dsp_insize; // size (complex samples) of the input buffer
int dsp_outsize; // size (complex samples) of the output buffer
int state; // 0 for unit OFF; 1 for unit ON
Unit(
int _in_rate, // input samplerate
int _out_rate, // output samplerate
int _dsp_rate, // sample rate for mainstream dsp processing
int _dsp_size // number complex samples processed per buffer in mainstream dsp processing
);
~Unit();
void flushBuffers();
void setBuffersInputSamplerate(int _in_rate);
void setBuffersOutputSamplerate(int _out_rate);
void setBuffersDSPSamplerate(int _dsp_rate);
void setBuffersDSPBuffsize(int _dsp_size);
protected:
float* inbuff;
float* midbuff;
float* outbuff;
};
} // namespace WDSP

View File

@ -146,7 +146,8 @@ void WCPAGC::flush()
void WCPAGC::execute()
{
int i, j, k;
int i;
int k;
double mult;
if (run)
@ -155,8 +156,8 @@ void WCPAGC::execute()
{
for (i = 0; i < io_buffsize; i++)
{
out[2 * i + 0] = fixed_gain * in[2 * i + 0];
out[2 * i + 1] = fixed_gain * in[2 * i + 1];
out[2 * i + 0] = (float) (fixed_gain * in[2 * i + 0]);
out[2 * i + 1] = (float) (fixed_gain * in[2 * i + 1]);
}
return;
@ -173,8 +174,10 @@ void WCPAGC::execute()
out_sample[0] = ring[2 * out_index + 0];
out_sample[1] = ring[2 * out_index + 1];
abs_out_sample = abs_ring[out_index];
double xr = ring[2 * in_index + 0] = in[2 * i + 0];
double xi = ring[2 * in_index + 1] = in[2 * i + 1];
ring[2 * in_index + 0] = in[2 * i + 0];
ring[2 * in_index + 1] = in[2 * i + 1];
double xr = ring[2 * in_index + 0];
double xi = ring[2 * in_index + 1];
if (pmode == 0)
abs_ring[in_index] = std::max(fabs(xr), fabs(xi));
@ -189,7 +192,7 @@ void WCPAGC::execute()
ring_max = 0.0;
k = out_index;
for (j = 0; j < attack_buffsize; j++)
for (int j = 0; j < attack_buffsize; j++)
{
if (++k == ring_buffsize)
k = 0;
@ -323,6 +326,8 @@ void WCPAGC::execute()
}
break;
}
default:
break;
}
if (volts < min_volts)
@ -330,8 +335,8 @@ void WCPAGC::execute()
gain = volts * inv_out_target;
mult = (out_target - slope_constant * std::min (0.0, log10(inv_max_input * volts))) / volts;
out[2 * i + 0] = out_sample[0] * mult;
out[2 * i + 1] = out_sample[1] * mult;
out[2 * i + 0] = (float) (out_sample[0] * mult);
out[2 * i + 1] = (float) (out_sample[1] * mult);
}
}
else if (out != in)
@ -406,7 +411,7 @@ void WCPAGC::setMode(int _mode)
void WCPAGC::setFixed(double _fixed_agc)
{
fixed_gain = pow (10.0, (double) _fixed_agc / 20.0);
fixed_gain = pow (10.0, _fixed_agc / 20.0);
loadWcpAGC();
}
@ -428,7 +433,7 @@ void WCPAGC::setHang(int _hang)
loadWcpAGC();
}
void WCPAGC::getHangLevel(double *hangLevel)
void WCPAGC::getHangLevel(double *hangLevel) const
//for line on bandscope
{
*hangLevel = 20.0 * log10(hang_level / 0.637);
@ -437,7 +442,8 @@ void WCPAGC::getHangLevel(double *hangLevel)
void WCPAGC::setHangLevel(double _hangLevel)
//for line on bandscope
{
double convert, tmp;
double convert;
double tmp;
if (max_input > min_volts)
{
@ -451,7 +457,7 @@ void WCPAGC::setHangLevel(double _hangLevel)
loadWcpAGC();
}
void WCPAGC::getHangThreshold(int *hangthreshold)
void WCPAGC::getHangThreshold(int *hangthreshold) const
//for slider in setup
{
*hangthreshold = (int) (100.0 * hang_thresh);
@ -464,7 +470,7 @@ void WCPAGC::setHangThreshold(int _hangthreshold)
loadWcpAGC();
}
void WCPAGC::getTop(double *max_agc)
void WCPAGC::getTop(double *max_agc) const
//for AGC Max Gain in setup
{
*max_agc = 20 * log10 (max_gain);
@ -473,7 +479,7 @@ void WCPAGC::getTop(double *max_agc)
void WCPAGC::setTop(double _max_agc)
//for AGC Max Gain in setup
{
max_gain = pow (10.0, (double) _max_agc / 20.0);
max_gain = pow (10.0, _max_agc / 20.0);
loadWcpAGC();
}
@ -489,9 +495,9 @@ void WCPAGC::setMaxInputLevel(double _level)
loadWcpAGC();
}
void WCPAGC::setRun(int state)
void WCPAGC::setRun(int _state)
{
run = state;
run = _state;
}
} // namespace WDSP

View File

@ -145,11 +145,11 @@ public:
void setAttack(int attack);
void setDecay(int decay);
void setHang(int hang);
void getHangLevel(double *hangLevel);
void getHangLevel(double *hangLevel) const;
void setHangLevel(double hangLevel);
void getHangThreshold(int *hangthreshold);
void getHangThreshold(int *hangthreshold) const;
void setHangThreshold(int hangthreshold);
void getTop(double *max_agc);
void getTop(double *max_agc) const;
void setTop(double max_agc);
void setSlope(int slope);
void setMaxInputLevel(double level);