Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

FFT Frequency Bins and PIC32

I am trying to derive the frequencies of an audio signal using the FFT library available for the PIC32MZ2064DAB176.

I am using the MPLAB Harmony for configuration.

For the sake of testing, two sine waves of frequencies 1002 Hz and 750 Hz are being used. This is done with the help of an online tone generator tool. I have 1002 Hz on one browser window and 750 Hz on another browser window. The output from the audio O/P jack is fed to the microcontroller ADC after a DC bias.

After performing a DC bias of 1.6 V the signal is sent to the 12-bit ADC. The maximum voltage I am expecting is a 3 V P-P, so I guess a DC bias of 1.6 V will suffice.

The signals are sampled at 48 kHz since I will need to read frequencies up to 20 kHz.

The FFT is a 1024 point FFT.

I am able to get the DC value in the 0th index of the frequency bin.

The formula being used to get the frequency value from the bin is Frequency = index * Sampling Frequency / Number of FFT Points

However, I am getting a high magnitude always in the 1st and the 2nd frequency bins for any value of the input frequency. According to my understanding, for 1002 Hz, the amplitude should be high around the 21st index of the frequency bin and for the 750 Hz signal, the amplitude should be high at around 16th index.

I am enclosing my code, the ADC Harmony configuration screenshot, the result screenshot and the signal input screenshot.

In the code, the array used for the frequency bin is "singleSidedFFT"

Any help in deriving the correct frequency value is greatly appreciated.

    /* FFT */
#define N 1024// Also change the log2N variable below!!
#define SAMPLE_FREQ 48000
#define PI 3.14

// Section: Global Data Definitions
APP_DATA appData;

/* ADC */
long count = 0;

/* FFT */
int16c  fftCoefs[N];
int16c *fftc;
int log2N = 10; 
extern const int16c twiddleFactors[];
long int freqVector[N];
int16c sampleBuffer[N]; //initialize buffer to collect samples
long int singleSidedFFT[N];


void APP_Tasks ( void )
{
    /* Check the application's current state. */
    switch ( appData.state )
    {
        /* Application's initial state. */
        case APP_STATE_INIT:
        {
            bool appInitialized = true;

            if (appInitialized)
            {
                int i;
                fftc = &fftCoefs; /* Stores the twiddle factors */

                // zero the freqVector and singleSidedFFT
                for (i=0; i<N; i++)
                {
                    freqVector = 0;
                    singleSidedFFT = 0;
                    sampleBuffer.re = 0;
                }

                // generate frequency vector this is the x-axis of your single sided fft
                for (i=0; i<N; i++)
                {
                    freqVector = i*(SAMPLE_FREQ/2)/((N/2) - 1);
                }

                /* Calculate the twiddle factors */
                DSP_TransformFFT16_setup(fftc, log2N);
                appData.state = APP_STATE_SERVICE_TASKS;

            }
            break;
        }

        case APP_STATE_SERVICE_TASKS:
        {
            /* Trigger a conversion */
            ADCCON3bits.GSWTRG = 1;

            /* Wait the conversions to complete */
            while (ADCDSTAT1bits.ARDY2 == 0);

            if (count < N)
            {
                sampleBuffer[count].re = ADCDATA2; /* fetch the result */
                sampleBuffer[count].im = 0;
                count++;
            }
            else
            {
                appData.state = APP_STATE_COMPUTE_FREQ;
                count = 0;
            }

            break;
        }

        case APP_STATE_COMPUTE_FREQ:
        {
            APP_ComputeFreq();
            appData.state = APP_STATE_SERVICE_TASKS;
            break;
        }
    }
}


void APP_ComputeFreq(void)
{
    int i;
    int16c dout[N]; //holds computed FFT 
    int16c scratch[N];

    // load complex input data into din
    DSP_TransformFFT16(dout, sampleBuffer, fftc, scratch, log2N);

    // compute single sided fft
    for(i = 0; i < N/2; i++)
    {
        singleSidedFFT = sqrt((dout.re*dout.re) + (dout.im*dout.im));
    }

    LATAbits.LATA6 = ~LATAbits.LATA6;
}

I have also tried writing a stand alone FFT function,. The result is the same. Here it is..

void APP_ComputeFreq_2(void)
{
    int16_t k, t;
    for (k = 0; k < N; k++) 
    { 
        // For each output element
        int16_t sumreal = 0;
        int16_t sumimag = 0;

        for (t = 0; t < N; t++) 
        { 
            // For each input element
            double angle = 2 * M_PI * t * k / N;
            sumreal += sampleBuffer[t].re * cos(angle) + sampleBuffer[t].im * sin(angle);
            sumimag += -sampleBuffer[t].re * sin(angle) + sampleBuffer[t].im * cos(angle);
        }
        singleSidedFFT[k] = sqrt((sumreal * sumreal) + (sumimag * sumimag));
    }
}

MPLAB Harmony ADC Config

ADC Frequency Bin

Input Signal

FFT Result 2

Thanks a lot.

like image 493
Ted Avatar asked Mar 12 '19 13:03

Ted


People also ask

What is frequency bin in FFT?

frequency bins are intervals between samples in frequency domain. For example, if your sample rate is 100 Hz and your FFT size is 100, then you have 100 points between [0 100) Hz. Therefore, you divide the entire 100 Hz range into 100 intervals, like 0-1 Hz, 1-2 Hz, and so on.

How many frequency bins are generated by the FFT?

FR = Fmax/N(Bins) With a 1024 FFT size, we divide this band into 512 bins.

How is FFT bin size calculated?

The bin width can be calculated by dividing the sample rate by the FFT length; or by dividing the bandwidth by the number of bins (which is equal to 1/2 the FFT length).

What is binning in signal processing?

Binning allows charges from adjacent pixels to be combined and this can offer benefits in faster readout speeds and improved signal to noise ratios albeit at the expense of reduced spatial resolution.


1 Answers

Ted uncovered an inconsistency between the datasheet of the microcontroller PIC32MZ Graphics (DA) Family and the version B of the specifications of 12 bit Successive Approximation Register (SAR) Analog-to-Digital Converter (ADC)

In both cases the clock source driving the sampling rate of the ADC is controlled by bits ADCSEL<1:0> of the register ADCCON3.

The datasheet, on page 452, gives the following clock sources:

11 = FRC
10 = REFCLK3
01 = System Clock (Tcy)
00 = PBCLK3

On the contrary, the specifications of the ADC, in version B, on page 14 are:

11 = System Clock (TCY)
10 = REFCLK3
01 = FRC Oscillator output
00 = Peripheral bus clock (PBCLK)

At the same point, version D of the specifications states that:

Refer to the “12-bit High-Speed Successive Approximation Register (SAR)” chapter in the specificdevice data sheet for the ADC Clock source selections.

The MPLAB Harmony ADC configurator complies with this disposition. Nevertheless, adopting the clock settings of the version B solved the sampling issue, suggesting that the family datasheet is not correct.

The sampling rate may also be affected by:

  • CONCLKDIV<5:0> ofADCCON3` : control clock divider
  • ADCDIV<6:0> of ADCxTIME : additionnal divider to define the clock source of indiviual ADC. Or alternatively ADCDIV<6:0> ofADCCON2` for shared ADC.
  • ADCxTIME<9:0> or ADCCON2<25:16>: number of clock ticks.

As the sampling was much higher than expected (625 kHz against 48kHz), the length of the frame (1024 samples = 0.0016s) was comparable to the periods of the input signal (about 1kHz). Most of the amplitude was therefore stored in the first bins of the DFT and applying a window does not solves the problem.

Once the sampling rate is corrected, the DFT features maxima corresponding to the frequencies or the input signal. These frequencies can be accurately identified by applying a window and estimating the frequency of a peak as its mean frequency wih respect to power density

like image 111
francis Avatar answered Oct 02 '22 06:10

francis