Hot!Problem with FFT power specteum

Author
akira_ss
New Member
  • Total Posts : 2
  • Reward points : 0
  • Joined: 2018/08/07 00:06:53
  • Location: 0
  • Status: offline
2018/08/07 17:32:34 (permalink)
0

Problem with FFT power specteum

Hi, sorry for my bad English.
I can get FFT spectrum in fftbuffer, but can not get power spectrum generated by the function SquareMagnitudeCplx().
Is there something wrong with my code, or program parameter wrong? 
i use dsPIC33FJ64GP802.
Someone please tell me.
 
#include <xc.h>
#include <dsp.h>
#include <timer.h>
#include <p33FJ64GP802.h>
#include <p33Fxxxx.h>
#define LOG2N 7
#define FFT_SIZE 128 //=2^7
#define NUM_TAPS 138
unsigned int Index;
unsigned int flag;
fractional SigIn[FFT_SIZE];
fractional FilterOut[FFT_SIZE];
fractcomplex fftBuffer[ FFT_SIZE ] __attribute__ ((space(ymemory), aligned (FFT_SIZE * 2 * 2))); // Place fftBuffer in Y-memory aligned to the size the array in bytes
fractional window[FFT_SIZE]; // Hamming window
fractional powerspec[FFT_SIZE]; // Power Spectrum
fractcomplex twiddleFactors[ FFT_SIZE/2 ]__attribute__((space(xmemory)));
fractional _YBSS(128) delay[FFT_SIZE];
const fractional lpfCoeffs[NUM_TAPS] __attribute__ ((space(auto_psv), aligned (FFT_SIZE * 2 * 2))) ={
1424,-2012,-1898,131,4372,8848,10772,8848,4372,131,-1898,-2012,1424
}; //filter coefficient
FIRStruct LPF; //LPF fc1=1000 fc2=2000
void initADC(void)
{
 AD1PCFGL = 0x03C0; //AN0-3 analog
 AD1CON1bits.SSRC = 2; //Internal counter ends sampling and starts conversion (auto-convert)
 AD1CON1bits.AD12B = 1; //12-bit, 1-channel ADC operation
 AD1CON1bits.ADON = 1; //ADC module is operating
    AD1CON1bits.ASAM = 1;
 AD1CON2 = 0x00; //ADCON2 no scan
 AD1CON3bits.SAMC = 0; //Autosample time
 AD1CON3 = 0b0000000000001001; // Tad={Tcy(ADCS+1)}/2>334ns, Then ADCS>18.7, Tad=10*Tcy
    AD1PCFGL=0xFFFF;
    //AD1PCFGH=0xFFFF;
    AD1CSSL = 0b0000000000000000;
    AD1PCFGLbits.PCFG1 = 0; //port pin in analogmode, port read input disabled, ADC samples pin voltage
    IFS0bits.AD1IF = 0; // Clear the A/D interrupt flag bit
    IEC0bits.AD1IE = 0; // Do Not Enable A/D interrupt
}

//return ADC result 32uS@Tcy=36.85MHz
unsigned int readADC()
{
 AD1CHS0 = 0x0001; //input channel in this AN1 is input
 AD1CON1bits.SAMP = 1; // ADC Sample Enable bit
 while(!AD1CON1bits.DONE); //wait until done ADC
 AD1CON1bits.DONE=0; //status = 0
 return ADC1BUF0; //result
}

// Timer subroutine
void __attribute__((interrupt, auto_psv)) _T3Interrupt(void)
{
    IFS0bits.T3IF = 0;
 while(!IFS0bits.AD1IF);
    SigIn[Index++] = readADC();
    //FIR (FFT_SIZE, &FilterOut[Index++], &SigIn[Index++], &HPF);
    if(Index == FFT_SIZE)
    {
        T3CONbits.TON = 0;
        flag = 1;
    }
    IFS0bits.AD1IF = 0;
    IFS0bits.T3IF = 0;
}


int main(void) {
    int i, j, k;
    
    /// Initialize Clock
 OSCCON = 0x00300; // PRIPLL??
 CLKDIV = 0x0100; // prescaler set1:2
 PLLFBD = 0x0026; // 40times?8MHz÷2×40÷2 = 80MHz
    RCONbits.SWDTEN = 0; //watchdog timer off
    
    initADC();
    OpenTimer3(T3_ON & T3_GATE_OFF & T3_PS_1_1 &T3_SOURCE_INT, 1600-1);
    ConfigIntTimer3(T3_INT_PRIOR_5 & T3_INT_ON); /// T3setting:T3= 10kHz = 16000kHz / 1600
     
    flag = 0;
    Index = 0;
    
    TwidFactorInit(LOG2N, twiddleFactors, 0);
    HammingInit (FFT_SIZE, window);
    FIRStructInit( &LPF, NUM_TAPS, (fractional *)lpfCoeffs, PSVPAG, delay );
    FIRDelayInit (&LPF);
    
    DISICNT = 0x0000; //enable interruput
    while(1)
    {
        //while(!flag);
        flag = 0;
        fftBuffer[0].real = 0;
        fftBuffer[1].real = 0;
        for(j = 0; j < FFT_SIZE; j++)
        FIR (1, &FilterOut[j], &SigIn[j], &LPF);
        VectorWindow(FFT_SIZE, FilterOut, FilterOut, window);
        for(i = 2; i < FFT_SIZE; i++) {
            fftBuffer[i].real = FilterOut[i];
            fftBuffer[i].imag = 0;
        }
        
        FFTComplexIP(LOG2N, fftBuffer, twiddleFactors, COEFFS_IN_DATA);
        BitReverseComplex(LOG2N, fftBuffer );
        SquareMagnitudeCplx( FFT_SIZE , fftBuffer, powerspec);

 
 
 
#1

3 Replies Related Threads

    DarioG
    Allmächtig.
    • Total Posts : 54081
    • Reward points : 0
    • Joined: 2006/02/25 08:58:22
    • Location: Oesterreich
    • Status: offline
    Re: Problem with FFT power specteum 2018/08/08 01:36:36 (permalink)
    0
    I use this code (excerpt)

    // Sampling Control
    #define Fosc        FCY
    #define Fcy            (Fosc/2)
    #define Fs           44100
    #define SAMPPRD    (Fcy/Fs)-1
    //#define NUMSAMP     512        //http://www.microchip.com/forums/fb.ashx?m=699633
    // bandwidth=170Hz => 256 samples per FFT
    #define LOG2_BLOCK_LENGTH     8    // = Number of "Butterfly" Stages in FFT processing
    #define FFT_BLOCK_LENGTH    (1<<LOG2_BLOCK_LENGTH)     // = Number of frequency points in the FFT
    #define SAMPLING_RATE        Fs    // = Rate at which input signal was sampled
                                        // SAMPLING_RATE is used to calculate the frequency
                                        // of the largest element in the FFT output vector


     
    #ifdef USA_DMA   // visto che non c'è il DMA, la sistemazione dei dati ADC la faccio in IRQ!
          p_real=&opBuffWindow1[0].real;
          // sta roba funziona così: le Window sono grosse il doppio, metà viene riempita consecutivamente
          // da IRQ o DMA con l'ADC; poi questo blocco prende ogni valore e lo sdoppia in real:complex (usando così l'intero buffer))
          p_cmpx=&opBuffWindow1[0];
          for(i=0; i<FFT_BLOCK_LENGTH; i++) {     // The FFT function requires input data
                            // to be in the fractional fixed-point range [-0.5, +0.5]
            (*p_real) >>= 1 ;        // So, we shift all data samples by 1 bit to the right.
            p_real++;            // Should you desire to optimize this process, perform
            }                    // data scaling when first obtaining the time samples
                            // Or within the BitReverseComplex function source code

          p_real=&opBuffWindow1[(FFT_BLOCK_LENGTH/2)-1].real;    // Set up pointers to convert real array
          p_cmpx=&opBuffWindow1[FFT_BLOCK_LENGTH-1]; // to a complex array. The input array initially has all
                  // the real input samples followed by a series of zeros

          for(i=FFT_BLOCK_LENGTH; i>0; i--) {   // Convert the Real input sample array
                            // to a Complex input sample array  
            (*p_cmpx).real = (*p_real--);    // We will simpy zero out the imaginary  
            (*p_cmpx--).imag = 0x0000;    // part of each data sample
            }
    #endif
    //        m_Led1Bit = 1;        // CHECK Timing!        
        
    //http://www.microchip.com/forums/m397450.aspx
    //http://www.microchip.com/stellent/idcplg?IdcService=SS_GET_PAGE&nodeId=1406&dDocName=en023598
                FFTComplexIP(LOG2_BLOCK_LENGTH, opBuffWindow1, (fractcomplex *)__builtin_psvoffset(&twiddleFactors[0]), (int)__builtin_psvpage(&twiddleFactors[0]));
            
          BitReverseComplex(LOG2_BLOCK_LENGTH, opBuffWindow1);
          SquareMagnitudeCplx(FFT_BLOCK_LENGTH, opBuffWindow1, &opBuffWindow1[0].real);
          VectorCopy(96, lcdData, &opBuffWindow1[0].real); // ogni Bin vale 44100/FFT_BLOCK_LENGTH ossia 170Hz circa
          // full range il valore è circa 2100... 1KHz
    //      m_Led1Bit = 0;        // CHECK Timing!    500uS @50MHz 8.11.17    

        // Compute the frequency (in Hz) of the largest spectral component
          VectorMax(FFT_BLOCK_LENGTH/2, &opBuffWindow1[0].real, &peakFrequencyBin);
          peakFrequency = peakFrequencyBin*(SAMPLING_RATE/FFT_BLOCK_LENGTH);

     
    Note the "adjust" part, that I use with DMA - in case.
     
    And of course you should have the constant array filled with proper values.

    GENOVA :D :D ! GODO
    #2
    akira_ss
    New Member
    • Total Posts : 2
    • Reward points : 0
    • Joined: 2018/08/07 00:06:53
    • Location: 0
    • Status: offline
    Re: Problem with FFT power specteum 2018/08/10 00:16:29 (permalink)
    3 (1)
    Thank you for the reply, dario.
     
    I changed my code, but still can't get power spectrum.
    Power spectrum's value return all 0. The fft buffer's value is okay and I can get fft specrtrum.
     
    I also tried FFT sample code by microchip. In the sample code, I can get the power spectrum.
    So, my ADC code or interrupt code maybe do something bad. But I can't find out where is wrong.
    #3
    DarioG
    Allmächtig.
    • Total Posts : 54081
    • Reward points : 0
    • Joined: 2006/02/25 08:58:22
    • Location: Oesterreich
    • Status: offline
    Re: Problem with FFT power specteum 2018/08/10 02:37:44 (permalink)
    0
    Ah ok, so FFT analysis is ok but Powers not... strange...

    GENOVA :D :D ! GODO
    #4
    Jump to:
    © 2018 APG vNext Commercial Version 4.5