New Analog Audio Filters WIP;

This commit is contained in:
lwvmobile 2024-03-19 11:03:55 -04:00
parent df28ebe573
commit 2b7de90474
7 changed files with 289 additions and 166 deletions

View File

@ -86,6 +86,35 @@
extern volatile uint8_t exitflag; //fix for issue #136
//new audio filter stuff from: https://github.com/NedSimao/FilteringLibrary
typedef struct {
float coef[2];
float v_out[2];
}LPFilter;
typedef struct {
float coef;
float v_out[2];
float v_in[2];
}HPFilter;
typedef struct {
LPFilter lpf;
HPFilter hpf;
float out_in;
}PBFilter;
typedef struct {
float alpha;
float beta;
float vin[3];
float vout[3];
}NOTCHFilter;
//end new filters
//group csv import struct
typedef struct
{
@ -285,6 +314,7 @@ typedef struct
int rtl_started;
long int rtl_rms;
int monitor_input_audio;
int analog_only;
int pulse_raw_rate_in;
int pulse_raw_rate_out;
int pulse_digi_rate_in;
@ -630,6 +660,13 @@ typedef struct
dPMRVoiceFS2Frame_t dPMRVoiceFS2Frame;
//new audio filter structs
LPFilter RCFilter;
HPFilter HRCFilter;
LPFilter RCFilter8;
HPFilter HRCFilter8;
// PBFilter PBF;
char dpmr_caller_id[20];
char dpmr_target_id[20];
@ -1285,12 +1322,22 @@ void eot_cc(dsd_opts * opts, dsd_state * state); //end of TX return to CC
//Generic Tuning Functions
void return_to_cc (dsd_opts * opts, dsd_state * state);
//misc generic audio filtering for analog at 48k/1
//misc audio filtering for analog
long int raw_rms(short *samples, int len, int step);
void analog_deemph_filter(short * input, int len);
void analog_preemph_filter(short * input, int len);
void analog_dc_block_filter(short * input, int len);
void analog_clipping_filter(short * input, int len);
void init_audio_filters(dsd_state * state);
void lpf(dsd_state * state, short * input, int len);
void hpf(dsd_state * state, short * input, int len);
void lpf8(dsd_state * state, short * input, int len);
void hpf8(dsd_state * state, short * input, int len);
//from: https://github.com/NedSimao/FilteringLibrary
void LPFilter_Init(LPFilter *filter, float cutoffFreqHz, float sampleTimeS);
float LPFilter_Update(LPFilter *filter, float v_in);
void HPFilter_Init(HPFilter *filter, float cutoffFreqHz, float sampleTimeS);
float HPFilter_Update(HPFilter *filter, float v_in);
void PBFilter_Init(PBFilter *filter, float HPF_cutoffFreqHz, float LPF_cutoffFreqHz, float sampleTimeS);
float PBFilter_Update(PBFilter *filter, float v_in);
void NOTCHFilter_Init(NOTCHFilter *filter, float centerFreqHz, float notchWidthHz, float sampleTimeS);
float NOTCHFilter_Update(NOTCHFilter *filter, float vin);
//csv imports
int csvGroupImport(dsd_opts * opts, dsd_state * state);

211
src/3.c
View File

@ -321,4 +321,213 @@ uint16_t q_abs_diff(const uint16_t v1, const uint16_t v2)
//
// Wojciech Kaczmarski, SP5WWP
// M17 Project, 29 December 2023
//--------------------------------------------------------------------
//--------------------------------------------------------------------
//audio filter stuff sourced from: https://github.com/NedSimao/FilteringLibrary
//no license / information provided in source code
#define PI 3.141592653
//do we need this for test input sthit?
// #define PLOT_POINTS 100
// #define SAMPLE_TIME_S 0.001
// #define CUTOFFFREQ_HZ 50
// #define sine_freq 40
// #define sine_freq_1 500
// #define HFC 50
// #define LFC 100
void LPFilter_Init(LPFilter *filter, float cutoffFreqHz, float sampleTimeS){
float RC=0.0;
RC=1.0/(2*PI*cutoffFreqHz);
filter->coef[0]=sampleTimeS/(sampleTimeS+RC);
filter->coef[1]=RC/(sampleTimeS+RC);
filter->v_out[0]=0.0;
filter->v_out[1]=0.0;
}
float LPFilter_Update(LPFilter *filter, float v_in){
filter->v_out[1]=filter->v_out[0];
filter->v_out[0]=(filter->coef[0]*v_in) + (filter->coef[1]*filter->v_out[1]);
return (filter->v_out[0]);
}
/********************************************************************************************************
* HIGH PASS FILTER
********************************************************************************************************/
void HPFilter_Init(HPFilter *filter, float cutoffFreqHz, float sampleTimeS){
float RC=0.0;
RC=1.0/(2*PI*cutoffFreqHz);
filter->coef=RC/(sampleTimeS+RC);
filter->v_in[0]=0.0;
filter->v_in[1]=0.0;
filter->v_out[0]=0.0;
filter->v_out[1]=0.0;
}
float HPFilter_Update(HPFilter *filter, float v_in){
filter->v_in[1]=filter->v_in[0];
filter->v_in[0]=v_in;
filter->v_out[1]=filter->v_out[0];
filter->v_out[0]=filter->coef * (filter->v_in[0] - filter->v_in[1]+filter->v_out[1]);
return (filter->v_out[0]);
}
/********************************************************************************************************
* BAND PASS FILTER
********************************************************************************************************/
void PBFilter_Init(PBFilter *filter, float HPF_cutoffFreqHz, float LPF_cutoffFreqHz, float sampleTimeS){
LPFilter_Init(&filter->lpf, LPF_cutoffFreqHz, sampleTimeS);
HPFilter_Init(&filter->hpf, HPF_cutoffFreqHz, sampleTimeS);
filter->out_in=0.0;
}
float PBFilter_Update(PBFilter *filter, float v_in){
filter->out_in=HPFilter_Update(&filter->hpf, v_in);
filter->out_in=LPFilter_Update(&filter->lpf, filter->out_in);
return (filter->out_in);
}
/********************************************************************************************************
* NOTCH FILTER
********************************************************************************************************/
void NOTCHFilter_Init(NOTCHFilter *filter, float centerFreqHz, float notchWidthHz, float sampleTimeS){
//filter frequency to angular (rad/s)
float w0_rps=2.0 * PI *centerFreqHz;
float ww_rps=2.0 * PI *notchWidthHz;
//pre warp center frequency
float w0_pw_rps=(2.0/sampleTimeS) * tanf(0.5 * w0_rps * sampleTimeS);
//computing filter coefficients
filter->alpha=4.0 + w0_rps*w0_pw_rps*sampleTimeS*sampleTimeS;
filter->beta=2.0*ww_rps*sampleTimeS;
//clearing input and output buffers
for (uint8_t n=0; n<3; n++){
filter->vin[n]=0;
filter->vout[n]=0;
}
}
float NOTCHFilter_Update(NOTCHFilter *filter, float vin){
//shifting samples
filter->vin[2]=filter->vin[1];
filter->vin[1]=filter->vin[0];
filter->vout[2]=filter->vout[1];
filter->vout[1]=filter->vout[0];
filter->vin[0]=vin;
//compute new output
filter->vout[0]=(filter->alpha*filter->vin[0] + 2.0 *(filter->alpha -8.0)*filter->vin[1] + filter->alpha*filter->vin[2]
-(2.0f*(filter->alpha-8.0)*filter->vout[1]+(filter->alpha-filter->beta)*filter->vout[2]))/(filter->alpha+filter->beta);
return (filter->vout[0]);
}
void init_audio_filters (dsd_state * state)
{
//still not sure if this is even correct or not, but 48k sounds good now
LPFilter_Init(&state->RCFilter, 960, (float)1/(float)48000);
HPFilter_Init(&state->HRCFilter, 960, (float)1/(float)48000);
//still unsure if 8k sounds good or not (may be the walkie-talkie is too close to SDR, even outside)
LPFilter_Init(&state->RCFilter8, 160, (float)1/(float)8000);
HPFilter_Init(&state->HRCFilter8, 160, (float)1/(float)8000);
//PBFilter_Init(PBFilter *filter, float HPF_cutoffFreqHz, float LPF_cutoffFreqHz, float sampleTimeS);
// void NOTCHFilter_Init(NOTCHFilter *filter, float centerFreqHz, float notchWidthHz, float sampleTimeS);
}
//FUNCTIONS for handing use of above filters
//lpf
void lpf(dsd_state * state, short * input, int len)
{
int i;
for (i = 0; i < len; i++)
input[i] = LPFilter_Update(&state->RCFilter, input[i]);
}
//hpf
void hpf(dsd_state * state, short * input, int len)
{
int i;
for (i = 0; i < len; i++)
input[i] = HPFilter_Update(&state->HRCFilter, input[i]);
}
//lpf
void lpf8(dsd_state * state, short * input, int len)
{
int i;
for (i = 0; i < len; i++)
input[i] = LPFilter_Update(&state->RCFilter8, input[i]);
}
//hpf
void hpf8(dsd_state * state, short * input, int len)
{
int i;
for (i = 0; i < len; i++)
input[i] = HPFilter_Update(&state->HRCFilter8, input[i]);
}
//Generic RMS function derived from RTL_FM (RTL_SDR) RMS code (doesnt' really work correctly outside of RTL)
long int raw_rms(int16_t *samples, int len, int step) //use samplespersymbol as len
{
int i;
long int rms;
long p, t, s;
double dc, err;
p = t = 0L;
for (i=0; i<len; i+=step) {
s = (long)samples[i];
t += s;
p += s * s;
}
/* correct for dc offset in squares */
dc = (double)(t*step) / (double)len;
err = t * 2 * dc - dc * dc * len;
rms = (long int)sqrt((p-err) / len);
if (rms < 0) rms = 150;
return rms;
}

View File

@ -195,7 +195,7 @@ getFrameSync (dsd_opts * opts, dsd_state * state)
int lmin, lmax, lidx;
//assign t_max value based on decoding type expected (all non-auto decodes first)
int t_max; //maximum values allowed for t will depend on decoding type - NXDN will be 10, others will be more
int t_max = 24; //initialize as an actual value to prevent any overflow related issues
if (opts->frame_nxdn48 == 1 || opts->frame_nxdn96 == 1)
{
t_max = 10;
@ -218,6 +218,10 @@ getFrameSync (dsd_opts * opts, dsd_state * state)
{
t_max = 19; //Phase 2 S-ISCH is only 19
}
else if (opts->analog_only == 1) //shim to make sure this is set to a reasonable value
{
t_max = 24;
}
else t_max = 24; //24 for everything else
int lbuf[48], lbuf2[48]; //if we use t_max in these arrays, and t >= t_max in condition below, then it can overflow those checks in there if t exceeds t_max

View File

@ -691,6 +691,7 @@ initOpts (dsd_opts * opts)
opts->dmr_mute_encR = 1;
opts->monitor_input_audio = 0; //enable with -8
opts->analog_only = 0; //only turned on with -fA
opts->inverted_p2 = 0;
opts->p2counter = 0;
@ -1636,7 +1637,7 @@ main (int argc, char **argv)
initOpts (&opts);
initState (&state);
init_audio_filters(&state); //audio filters
InitAllFecFunction();
// CNXDNConvolution_init(); //seems to function better without initting it
@ -2103,6 +2104,7 @@ main (int argc, char **argv)
opts.dmr_mono = 0;
state.rf_mod = 0;
opts.monitor_input_audio = 1;
opts.analog_only = 1;
sprintf (opts.output_name, "Analog Monitor");
fprintf (stderr,"Only Monitoring Passive Analog Signal\n");
}

View File

@ -17,29 +17,6 @@
#include "dsd.h"
long int raw_rms(int16_t *samples, int len, int step) //use samplespersymbol as len
{
int i;
long int rms;
long p, t, s;
double dc, err;
p = t = 0L;
for (i=0; i<len; i+=step) {
s = (long)samples[i];
t += s;
p += s * s;
}
/* correct for dc offset in squares */
dc = (double)(t*step) / (double)len;
err = t * 2 * dc - dc * dc * len;
rms = (long int)sqrt((p-err) / len);
if (rms < 0) rms = 150;
return rms;
}
int
getSymbol (dsd_opts * opts, dsd_state * state, int have_sync)
{
@ -301,14 +278,9 @@ getSymbol (dsd_opts * opts, dsd_state * state, int have_sync)
sf_write_sync (opts->wav_out_raw);
}
//test running analog audio through a de-emphasis filter
analog_deemph_filter(state->analog_out, 960);
//and dc_block filter analog_dc_block_filter
analog_dc_block_filter(state->analog_out, 960);
//test running analog audio through a pre-emphasis filter
analog_preemph_filter(state->analog_out, 960);
//test running analog short sample clipping filter
analog_clipping_filter(state->analog_out, 960);
// analog audio filtering
// lpf (state, state->analog_out, 960);
hpf (state, state->analog_out, 960);
//seems to be working now, but RMS values are lower on actual analog signal than on no signal but noise
if ( (opts->rtl_rms > opts->rtl_squelch_level) && opts->monitor_input_audio == 1 && state->carrier == 0 ) //added carrier check here in lieu of disabling it above

View File

@ -66,118 +66,9 @@ void openWavOutFile48k (dsd_opts * opts, dsd_state * state)
}
}
//generic rms type function
long int gen_rms(short *samples, int len, int step)
{
int i;
long int rms;
long p, t, s;
double dc, err;
p = t = 0L;
for (i=0; i<len; i+=step) {
s = (long)samples[i];
t += s;
p += s * s;
}
/* correct for dc offset in squares */
dc = (double)(t*step) / (double)len;
err = t * 2 * dc - dc * dc * len;
rms = (long int)sqrt((p-err) / len);
//make sure it doesnt' randomly feed us a large negative value (overflow?)
if (rms < 0) rms = 150; //could also consider returning 0 and making it return the last known good value instead
return rms;
}
//super generic pre-emphasis filter with alpha value of 1 for short sample input / output
static short pstate = 1;
void analog_preemph_filter(short * input, int len)
{
int i;
for (i = 0; i < len; i++)
{
//debug
// fprintf (stderr, " I:%d;", input[i]);
input[i] = input[i] - 1 * pstate;
pstate = input[i];
//debug
// fprintf (stderr, " O:%d;;", input[i]);
//increase gain value slightly
input[i] *= 1.5;
//debug
// fprintf (stderr, " O2:%d;;;", input[i]);
}
}
//super generic short sample clipping filter, prevent short values from exceeding tolerable clip value
void analog_clipping_filter(short * input, int len)
{
int i;
for (i = 0; i < len; i++)
{
if (input[i] > 32760)
input[i] = 32760;
else if (input[i] < -32760)
input[i] = -32760;
}
}
//modified from version found in rtl_fm
static short avg = 0;
void analog_deemph_filter(short * input, int len)
{
// static short avg; // cheating...
int i, d;
// int a = (int)round(1.0/((1.0-exp(-1.0/(96000 * 75e-6))))); //48000 is rate out, but doubling it sounds 'beefier'
int a = 8; //result of calc above is 8, so save a cycle on the low end CPUs
// de-emph IIR
// avg = avg * (1 - alpha) + sample * alpha;
for (i = 0; i < len; i++)
{
d = input[i] - avg;
if (d > 0)
avg += (d + a/2) / a;
else
avg += (d - a/2) / a;
input[i] = (short)avg;
}
}
//modified from version found in rtl_fm
static int dc_avg = 0;
void analog_dc_block_filter(short * input, int len)
{
int i, avg;
int64_t sum = 0;
for (i=0; i < len; i++)
sum += input[i];
avg = sum / len;
avg = (avg + dc_avg * 9) / 10;
for (i=0; i < len; i++)
input[i] -= (short)avg;
dc_avg = avg;
}
//listening to and playing back analog audio
void edacs_analog(dsd_opts * opts, dsd_state * state, int afs, unsigned char lcn)
{
//reset static states for new channel
pstate = 1;
avg = 0;
dc_avg = 0;
int i, result;
int count = 5; //RMS has a 5 count (5 * 180ms) now before cutting off;
short analog1[960];
@ -220,7 +111,7 @@ void edacs_analog(dsd_opts * opts, dsd_state * state, int afs, unsigned char lcn
analog3[i] = sample;
}
//this rms will only work properly (for now) with squelch enabled in SDR++ or other
rms = gen_rms(analog3, 960, 1);
rms = raw_rms(analog3, 960, 1);
}
//NOTE: The core dumps observed previously were due to SDR++ Remote Server connection dropping due to Internet/Other issues
@ -274,7 +165,7 @@ void edacs_analog(dsd_opts * opts, dsd_state * state, int afs, unsigned char lcn
}
//this rms will only work properly (for now) with squelch enabled in SDR++
rms = gen_rms(analog3, 960, 1);
rms = raw_rms(analog3, 960, 1);
}
//RTL Input
@ -311,22 +202,14 @@ void edacs_analog(dsd_opts * opts, dsd_state * state, int afs, unsigned char lcn
sr += digitize (opts, state, (int)analog1[i]);
}
//test running analog audio through a de-emphasis filter
analog_deemph_filter(analog1, 960);
analog_deemph_filter(analog2, 960);
analog_deemph_filter(analog3, 960);
//and dc_block filter analog_dc_block_filter
analog_dc_block_filter(analog1, 960);
analog_dc_block_filter(analog2, 960);
analog_dc_block_filter(analog3, 960);
//test running analog audio through a pre-emphasis filter
analog_preemph_filter(analog1, 960);
analog_preemph_filter(analog2, 960);
analog_preemph_filter(analog3, 960);
//test running analog short sample clipping filter
analog_clipping_filter(analog1, 960);
analog_clipping_filter(analog2, 960);
analog_clipping_filter(analog3, 960);
// analog audio filtering
lpf (state, analog1, 960);
lpf (state, analog2, 960);
lpf (state, analog3, 960);
hpf (state, analog1, 960);
hpf (state, analog2, 960);
hpf (state, analog3, 960);
//reconfigured to use seperate audio out stream that is always 48k short
if (opts->audio_out_type == 0 && opts->slot1_on == 1)

View File

@ -1852,6 +1852,12 @@ void encodeM17STR(dsd_opts * opts, dsd_state * state)
voice2[i] *= (float)state->aout_gain/ (float)25.0f;
}
//run lpf and hpf
// lpf8 (state, voice1, 160);
// lpf8 (state, voice2, 160);
hpf8 (state, voice1, 160);
hpf8 (state, voice2, 160);
//convert out audio input into CODEC2 (3200bps) 8 byte data stream
uint8_t vc1_bytes[8]; memset (vc1_bytes, 0, sizeof(vc1_bytes));
uint8_t vc2_bytes[8]; memset (vc2_bytes, 0, sizeof(vc2_bytes));