Tympan / Tympan_Library

Arduino/Teensy Library for Tympan Open Source Hearing Aid
MIT License
116 stars 31 forks source link

Audio noise reduction not working #63

Open michaelliesenberg opened 2 years ago

michaelliesenberg commented 2 years ago

Hi,

i am trying to use the tympan library for my custom board. I use the TLV32ßADC6140 and everything is working fine.

I have a stereo input at 192kHz and i do an average down to 48kHz and after that i am trying to do a noice reduction but the audio output is zero. If i uncomment "#define FFT_Noise_Subtraction", everything is working fine, but using the noise subtraction algorithmus i dont see any audio signal output.

Here is my code:

#include <Arduino.h>
#include "Filter_FFT.h"
#include <stdio.h>
#include <Tympan_Library.h>
#include "AudioEffectNoiseReduction_FD_F32.h"

#define FFT_Noise_Subtraction //If i turn this off everything is working fine
#define fiter_low_high_befor_FFT
#define fiter_low_high
#define AVG
#define PHASE_GAIN_CORRECTION

elapsedMicros sincePrint;
int helptime=0;

//set the sample rate and block size
const float sample_rate_Hz = 48000.f;   //24000 or 44117 (or other frequencies in the table in AudioOutputI2S_F32)
const int audio_block_samples = 128;       //for freq domain processing choose a power of 2 (16, 32, 64, 128) but no higher than 128
const int FFT_overlap_factor = 1;         //2 is 50% overlap, 4 is 75% overlap

AudioSettings_F32 audio_settings(sample_rate_Hz, audio_block_samples);
AudioEffectNoiseReduction_FD_F32    noiseReduction_right(audio_settings);
AudioEffectNoiseReduction_FD_F32    noiseReduction_left(audio_settings);

FFT_Overlapped_F32 myFFT_Right;
FFT_Overlapped_F32 myFFT_Left;
IFFT_Overlapped_F32 myIFFT_Right;
IFFT_Overlapped_F32 myIFFT_Left;
audio_block_f32_t *myBlock_Right;
audio_block_f32_t *myBlock_Left;

audio_block_f32_t* allocate_f32_memory_help(const int num) {
  static bool firstTime=true;
  static audio_block_f32_t *data_f32;
  if (firstTime == true) {
    firstTime = false;
    data_f32 = new audio_block_f32_t[num];
  }
  return data_f32;
}

Filter_FFT::Filter_FFT(int samplerate , int oversample, DATASTRUCT &datastruct): AudioStream(2,inputQueueArray)
{
  myData = &datastruct;

  dataReady=false;
  count=0;
  countAVG=0;
  SAMPLE_RATE=samplerate;
  oversampling = oversample/samplerate;
  STATE = average;

  left=NULL;
  right=NULL;

  int freq_hp = 100;//70 100
  int freq_tp = 2000;

  lowpass_average = new Dsp::FilterDesign<Dsp::Butterworth::Design::LowPass <order_lp>, 2, Dsp::DirectFormII> ;
  Dsp::Params params;
  params[0] = oversample; // sample rate
  params[1] = order_lp; // order
  params[2] = 20000; // cutoff frequency
  lowpass_average->setParams (params);

  lowpass = new Dsp::FilterDesign<Dsp::Butterworth::Design::LowPass <order_lp>, 2, Dsp::DirectFormII> ;
  params[0] = samplerate; // sample rate
  params[1] = order_lp; // order
  params[2] = freq_tp; // cutoff frequency
  lowpass->setParams (params);

  highpass = new Dsp::FilterDesign<Dsp::Butterworth::Design::HighPass <order_hp>,2, Dsp::DirectFormII>;
  params[0] = samplerate; // sample rateoutput
  params[1] = order_hp; // order
  params[2] = freq_hp;//25 // cutoff frequency
  highpass->setParams (params);

  lowpass_second= new Dsp::FilterDesign<Dsp::Butterworth::Design::LowPass <order_lp>, 2, Dsp::DirectFormII> ;
  params[0] = samplerate; // sample rate
  params[1] = order_lp; // order
  params[2] = freq_tp; // cutoff frequency
  lowpass_second->setParams (params);

  highpass_second = new Dsp::FilterDesign<Dsp::Butterworth::Design::HighPass <order_hp>,2, Dsp::DirectFormII>;
  params[0] = samplerate; // sample rateoutput
  params[1] = order_hp; // order
  params[2] = freq_hp;//25 // cutoff frequency
  highpass_second->setParams (params);

  audiodata[0] = new float[FFT_SIZE*2];
  audiodata[1] = new float[FFT_SIZE*2];
  complex_2N_buffer[0] = new float32_t[2 * FFT_SIZE];
  complex_2N_buffer[1] = new float32_t[2 * FFT_SIZE];

  AudioMemory_F32(30, audio_settings);
  int N_FFT = audio_block_samples * FFT_overlap_factor;
  noiseReduction_right.setup(audio_settings, N_FFT);
  noiseReduction_left.setup(audio_settings, N_FFT);
  myFFT_Right.setup(audio_settings, N_FFT); 
  myFFT_Left.setup(audio_settings, N_FFT); 
  myIFFT_Right.setup(audio_settings, N_FFT); 
  myIFFT_Left.setup(audio_settings, N_FFT); 

  (myFFT_Right.getFFTObject())->useHanningWindow();
  (myFFT_Left.getFFTObject())->useHanningWindow();
  (myIFFT_Right.getFFTObject())->useHanningWindow();
  (myIFFT_Left.getFFTObject())->useHanningWindow();
  //Frequency processing
  noiseReduction_right.setAttack_sec(1);//10
  noiseReduction_left.setAttack_sec(1);
  noiseReduction_right.setRelease_sec(1); //3
  noiseReduction_left.setRelease_sec(1); 
  noiseReduction_right.setMaxAttenuation_dB(0); 
  noiseReduction_left.setMaxAttenuation_dB(0); //16
  noiseReduction_right.setSNRforMaxAttenuation_dB(1.0);//1
  noiseReduction_left.setSNRforMaxAttenuation_dB(1.0);
  noiseReduction_right.setTransitionWidth_dB(6.0); 
  noiseReduction_left.setTransitionWidth_dB(6.0); 
  noiseReduction_right.setGainSmoothing_sec(0.005);
  noiseReduction_left.setGainSmoothing_sec(0.005);
  noiseReduction_right.setEnableNoiseEstimationUpdates(true); 
  noiseReduction_left.setEnableNoiseEstimationUpdates(true); 

  myBlock_Right = allocate_f32_memory_help(0);
  if (myBlock_Right != NULL) AudioStream_F32::initialize_f32_memory(myBlock_Right, 0);
    myBlock_Left = allocate_f32_memory_help(1);
  if (myBlock_Left != NULL) AudioStream_F32::initialize_f32_memory(myBlock_Left, 1);
}

void Filter_FFT::update(void)
{
  audio_block_t *left, *right;

  left = receiveWritable(0); // input 0 = left channel
  right = receiveWritable(1); // input 1 = right channel

    if(!right || !left)
        return;

  process(right, left);

}

void Filter_FFT::process(audio_block_t *right, audio_block_t *left)
{
#ifdef AVG

  double avg_r, avg_l;

  for(int i=0; i<AUDIO_BLOCK_SAMPLES;i=i+oversampling)
  {
    avg_r=right->data[i]/oversampling;
    avg_l=right->data[i]/oversampling;

    for(int j=1; j<oversampling;j++)
    {
      avg_r = avg_r + right->data[i+j]/oversampling;
      avg_l = avg_l + left->data[i+j]/oversampling;
    }

    array_avg[0][countAVG] = (int16_t)avg_r;
    array_avg[1][countAVG] = (int16_t)avg_l;

    countAVG++;
  }

  if(countAVG >= AUDIO_BLOCK_SAMPLES)
  {
    memcpy(&right->data[0],&array_avg[0][0],AUDIO_BLOCK_SAMPLES * sizeof(int16_t)); 
    memcpy(&left->data[1],&array_avg[1][0],AUDIO_BLOCK_SAMPLES * sizeof(int16_t));       
    countAVG=0;
  }
  else
  {
    release(right);
    release(left);
    return; 
  }

#endif

#ifdef fiter_low_high_befor_FFT

    for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
    {
      audiodata[0][i] = (float)(right->data[i]);
      audiodata[1][i] = (float)(left->data[i]);
    }

    highpass->process (AUDIO_BLOCK_SAMPLES, audiodata);
    lowpass->process (AUDIO_BLOCK_SAMPLES, audiodata);

#ifdef PHASE_GAIN_CORRECTION
    arm_scale_f32 (audiodata[0], (float)(int(myData->GAIN_MAINBOARD[0]))/1000 * (float)(int(myData->GAIN[0]))/1000, audiodata[0], AUDIO_BLOCK_SAMPLES);
    arm_scale_f32 (audiodata[1], (float)(int(myData->GAIN_MAINBOARD[1]))/1000 * (float)(int(myData->GAIN[1]))/1000, audiodata[1], AUDIO_BLOCK_SAMPLES);
    IQ_phase_correction(audiodata[0], audiodata[1], (float)(int(myData->PHASE[1]))/1000, AUDIO_BLOCK_SAMPLES);
    IQ_phase_correction(audiodata[0], audiodata[1], (float)(int(myData->PHASE[0]))/1000, AUDIO_BLOCK_SAMPLES);
#endif

    for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
    {
      right->data[i] = (int)(audiodata[0][i]);  
      left->data[i] = (int)(audiodata[1][i]);  
#ifdef FFT_Noise_Subtraction 
      myBlock_Right->data[i] = audiodata[0][i];
      myBlock_Left->data[i] = audiodata[1][i];
#endif

    }

#endif

#ifdef FFT_Noise_Subtraction 

  audio_block_f32_t *in_audio_block_right = myBlock_Right;
  audio_block_f32_t *in_audio_block_left = myBlock_Left;

  //convert to frequency domain
  myFFT_Right.execute(in_audio_block_right, complex_2N_buffer[0]);
  myFFT_Left.execute(in_audio_block_left, complex_2N_buffer[1]);

  AudioStream_F32::release(in_audio_block_right);
  AudioStream_F32::release(in_audio_block_left);
  //define some variables
  noiseReduction_right.processAudioFD(complex_2N_buffer[0], myFFT_Right.getNFFT());  //in your derived class, overwrite this!!
  noiseReduction_left.processAudioFD(complex_2N_buffer[1], myFFT_Left.getNFFT());  //in your derived class, overwrite this!!

  //rebuild the negative frequency space
  myFFT_Right.rebuildNegativeFrequencySpace(complex_2N_buffer[0]); //set the negative frequency space based on the positive
  myFFT_Left.rebuildNegativeFrequencySpace(complex_2N_buffer[1]); //set the negative frequency space based on the positive

  //call the IFFT

  audio_block_f32_t *out_audio_block_right = myIFFT_Right.execute(complex_2N_buffer[0]); //out_block is pre-allocated in here.
  audio_block_f32_t *out_audio_block_left = myIFFT_Left.execute(complex_2N_buffer[1]); //out_block is pre-allocated in here.

    for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
    {
      audiodata[0][i] = out_audio_block_right->data[i];  
      audiodata[1][i] = out_audio_block_left->data[i];        
    }
#else
  for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
  {
    audiodata[0][i] = float(right->data[i]);
    audiodata[1][i] = float(left->data[i]);
  }
#endif

#ifdef fiter_low_high

  lowpass_second->process (AUDIO_BLOCK_SAMPLES, audiodata);
  highpass_second->process (AUDIO_BLOCK_SAMPLES, audiodata);

  for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
  {        
    right->data[i] = (int)(audiodata[0][i]);  
    left->data[i] = (int)(audiodata[1][i]);  
  }
#endif

  transmit(right,0);
  release(right);

  transmit(left,1);
  release(left);

}

void Filter_FFT::set_LOW_PASS_FREQUENCY(int value)
{  
  Dsp::Params params;

  params[0] = SAMPLE_RATE; // sample rate
  params[1] = order_lp; // order
  params[2] = value; // cutoff frequency
  lowpass->setParams (params);

  params[0] = SAMPLE_RATE; // sample rate
  params[1] = order_lp; // order
  params[2] = value; // cutoff frequency
  lowpass_second->setParams (params);  
}
void Filter_FFT::set_HIGH_PASS_FREQUENCY(int value)
{
  Dsp::Params params;

  params[0] = SAMPLE_RATE; // sample rateoutput
  params[1] = order_hp; // order
  params[2] = value;//25 // cutoff frequency
  highpass->setParams (params);

  params[0] = SAMPLE_RATE; // sample rateoutput
  params[1] = order_hp; // order
  params[2] = value;//25 // cutoff frequency
  highpass_second->setParams (params);  
}

void Filter_FFT::IQ_phase_correction (float32_t *I_buffer, float32_t *Q_buffer, float32_t factor, uint32_t blocksize)
{
  float32_t temp_buffer[blocksize];
  if (factor < 0.0)
  { // mix a bit of I into Q
    arm_scale_f32 (I_buffer, factor, temp_buffer, blocksize);
    arm_add_f32 (Q_buffer, temp_buffer, Q_buffer, blocksize);
  }
  else
  { // // mix a bit of Q into I
    arm_scale_f32 (Q_buffer, factor, temp_buffer, blocksize);
    arm_add_f32 (I_buffer, temp_buffer, I_buffer, blocksize);
  }
} // end IQ_phase_correction

Any idea?

chipaudette commented 2 years ago

Hi,

You have a lot happening in your code, so it's hard to know what is going wrong. Here are some thoughts:

Serial.print("MEM Cur/Pk: ");
Serial.print(AudioMemoryUsage_F32());
Serial.print("/");
Serial.print(AudioMemoryUsageMax_F32());
Serial.println();
 //convert to frequency domain
  myFFT_Right.execute(in_audio_block_right, complex_2N_buffer[0]);
  myFFT_Left.execute(in_audio_block_left, complex_2N_buffer[1]);

  AudioStream_F32::release(in_audio_block_right);
  AudioStream_F32::release(in_audio_block_left);
  //define some variables

#if 0  //ADD THIS LINE
  noiseReduction_right.processAudioFD(complex_2N_buffer[0], myFFT_Right.getNFFT());  //in your derived class, overwrite this!!
  noiseReduction_left.processAudioFD(complex_2N_buffer[1], myFFT_Left.getNFFT());  //in your derived class, overwrite this!!
#endif  //AND ADD THIS LINE

  //rebuild the negative frequency space
  myFFT_Right.rebuildNegativeFrequencySpace(complex_2N_buffer[0]); //set the negative frequency space based on the positive
  myFFT_Left.rebuildNegativeFrequencySpace(complex_2N_buffer[1]); //set the negative frequency space based on the positive

  //call the IFFT

  audio_block_f32_t *out_audio_block_right = myIFFT_Right.execute(complex_2N_buffer[0]); //out_block is pre-allocated in here.
  audio_block_f32_t *out_audio_block_left = myIFFT_Left.execute(complex_2N_buffer[1]); //out_block is pre-allocated in here.

Chip

michaelliesenberg commented 2 years ago

Hi Chip thank you for the help.

i have changed the FFT_overlap_factor to 2

const int audio_block_samples = 128
const int FFT_overlap_factor = 2 

and set

AudioMemory_F32(60, audio_settings);

Now if i uncomment

  noiseReduction_right.processAudioFD(complex_2N_buffer[0], myFFT_Right.getNFFT());  
  noiseReduction_left.processAudioFD(complex_2N_buffer[1], myFFT_Left.getNFFT());  

The audio signal is there but with artefacts and if i put the lines back, there is no audio signal coming anymore.

Here is the Init in constructor:

  int N_FFT = audio_block_samples * FFT_overlap_factor;
  noiseReduction_right.setup(audio_settings, N_FFT);
  noiseReduction_left.setup(audio_settings, N_FFT);
  myFFT_Right.setup(audio_settings, N_FFT); 
  myFFT_Left.setup(audio_settings, N_FFT); 
  myIFFT_Right.setup(audio_settings, N_FFT); 
  myIFFT_Left.setup(audio_settings, N_FFT); 

  (myFFT_Right.getFFTObject())->useHanningWindow();
  (myFFT_Left.getFFTObject())->useHanningWindow();
  (myIFFT_Right.getFFTObject())->useHanningWindow();
  (myIFFT_Left.getFFTObject())->useHanningWindow();
  //Frequency processing
  noiseReduction_right.setAttack_sec(1);//10
  noiseReduction_left.setAttack_sec(1);
  noiseReduction_right.setRelease_sec(1); //3
  noiseReduction_left.setRelease_sec(1); 
  noiseReduction_right.setMaxAttenuation_dB(0); 
  noiseReduction_left.setMaxAttenuation_dB(0); //16
  noiseReduction_right.setSNRforMaxAttenuation_dB(1.0);//1
  noiseReduction_left.setSNRforMaxAttenuation_dB(1.0);
  noiseReduction_right.setTransitionWidth_dB(6.0); 
  noiseReduction_left.setTransitionWidth_dB(6.0); 
  noiseReduction_right.setGainSmoothing_sec(0.005);
  noiseReduction_left.setGainSmoothing_sec(0.005);
  noiseReduction_right.setEnableNoiseEstimationUpdates(true); 
  noiseReduction_left.setEnableNoiseEstimationUpdates(true); 

  myBlock_Right = allocate_f32_memory_help(0);
  if (myBlock_Right != NULL) AudioStream_F32::initialize_f32_memory(myBlock_Right, 0);
    myBlock_Left = allocate_f32_memory_help(1);
  if (myBlock_Left != NULL) AudioStream_F32::initialize_f32_memory(myBlock_Left, 1);

And here is the code inside update for noise removal:


 audio_block_f32_t *in_audio_block_right = myBlock_Right;
  audio_block_f32_t *in_audio_block_left = myBlock_Left;

  //convert to frequency domain
  myFFT_Right.execute(in_audio_block_right, complex_2N_buffer[0]);
  myFFT_Left.execute(in_audio_block_left, complex_2N_buffer[1]);

  AudioStream_F32::release(in_audio_block_right);
  AudioStream_F32::release(in_audio_block_left);

  //define some variables
  noiseReduction_right.processAudioFD(complex_2N_buffer[0], myFFT_Right.getNFFT());  //in your derived class, overwrite this!!
  noiseReduction_left.processAudioFD(complex_2N_buffer[1], myFFT_Left.getNFFT());  //in your derived class, overwrite this!!

  //rebuild the negative frequency space
  myFFT_Right.rebuildNegativeFrequencySpace(complex_2N_buffer[0]); //set the negative frequency space based on the positive
  myFFT_Left.rebuildNegativeFrequencySpace(complex_2N_buffer[1]); //set the negative frequency space based on the positive

  //call the IFFT

  audio_block_f32_t *out_audio_block_right = myIFFT_Right.execute(complex_2N_buffer[0]); //out_block is pre-allocated in here.
  audio_block_f32_t *out_audio_block_left = myIFFT_Left.execute(complex_2N_buffer[1]); //out_block is pre-allocated in here.

    for (int i=0; i<AUDIO_BLOCK_SAMPLES; i++) 
    {
      audiodata[0][i] = out_audio_block_right->data[i];  
      audiodata[1][i] = out_audio_block_left->data[i];   
    }

Do you think i have to change the noise settings ? Does it record the noise at beginning or it subtracts the noise with the parameters on the struct? And is it normal that once i do FFT and back iFFT i get some artefacts?