fft_code.zip_64点fft_EXCEL FFT_FFT excel_excel FFT_excel 计算fft

//********************************************************************** // * * // * Software License Agreement * // * * // * The software supplied herewith by Microchip Technology * // * Incorporated (the "Company") for its dsPIC controller * // * is intended and supplied to you, the Company's customer, * // * for use solely and exclusively on Microchip dsPIC * // * products. The software is owned by the Company and/or its * // * supplier, and is protected under applicable copyright laws. All * // * rights are reserved. Any use in violation of the foregoing * // * restrictions may subject the user to criminal sanctions under * // * applicable laws, as well as to civil liability for the breach of * // * the terms and conditions of this license. * // * * // * THIS SOFTWARE IS PROVIDED IN AN "AS IS" CONDITION. NO * // * WARRANTIES, WHETHER EXPRESS, IMPLIED OR STATUTORY, INCLUDING, * // * BUT NOT LIMITED TO, IMPLIED WARRANTIES OF MERCHANTABILITY AND * // * FITNESS FOR A PARTICULAR PURPOSE APPLY TO THIS SOFTWARE. THE * // * COMPANY SHALL NOT, IN ANY CIRCUMSTANCES, BE LIABLE FOR SPECIAL, * // * INCIDENTAL OR CONSEQUENTIAL DAMAGES, FOR ANY REASON WHATSOEVER. * // * * // ******************************************************************** #include <p30f3013.h> #include <dsp.h> #define PRE_PROCESSING #define BLOCK_LENGTH 64 #define LOG2_BLOCK_LENGTH 6 //2^6 = 64 //sine, cos table - pre-calculate and store in the FLASH, so we do not need to use the TwidFactorInit() to calculate the factors const fractcomplex twiddleFactors[] __attribute__ ((space(auto_psv), aligned (BLOCK_LENGTH*2)))= { 0x7FFF, 0x0000, 0x7F62, 0xF374, 0x7D8A, 0xE707, 0x7A7D, 0xDAD8, 0x7642, 0xCF04, 0x70E3, 0xC3A9, 0x6A6E, 0xB8E3, 0x62F2, 0xAECC, 0x5A82, 0xA57E, 0x5134, 0x9D0E, 0x471D, 0x9592, 0x3C57, 0x8F1D, 0x30FC, 0x89BE, 0x2528, 0x8583, 0x18F9, 0x8276, 0x0C8C, 0x809E, 0x0000, 0x8000, 0xF374, 0x809E, 0xE707, 0x8276, 0xDAD8, 0x8583, 0xCF04, 0x89BE, 0xC3A9, 0x8F1D, 0xB8E3, 0x9592, 0xAECC, 0x9D0E, 0xA57D, 0xA57D, 0x9D0E, 0xAECC, 0x9592, 0xB8E3, 0x8F1D, 0xC3A9, 0x89BE, 0xCF04, 0x8583, 0xDAD8, 0x8276, 0xE707, 0x809E, 0xF374 } ; fractcomplex FFTData[BLOCK_LENGTH] __attribute__ ((section (".ybss"), aligned (BLOCK_LENGTH*4))); fractional FFTMagnitude[BLOCK_LENGTH/2]; //to store the magnitude float THD[BLOCK_LENGTH/2]; //to store the percentage of each order harmonics(2~31 order) float THD_ALL; int data[BLOCK_LENGTH]__attribute__ ((space (xmemory)));//data to store the ADC value, here is integer, but when FFT calculation, it will be fractional unsigned char i, k; int temp_fft_data_max, temp_fft_data_min; int FFTPreProcessCnt; int main (void) { for(i=1;i<BLOCK_LENGTH/2;i++) { THD[i] = 0; } for (i=0; i<BLOCK_LENGTH; i++) data[i] = 0; //Load the ADC data to data[] array here //use WATCH window, then right click the mouse on the variable of data[] //then use the "import table" command on the pop-up menu to load the ADC data stored in the .mch file FFTPreProcessCnt = 0; #ifdef PRE_PROCESSING temp_fft_data_max = data[0]; temp_fft_data_min = data[0]; //find the maximum and minmum value in the data[] array for (k=0; k<BLOCK_LENGTH; k++) { if(data[k] > temp_fft_data_max) temp_fft_data_max = data[k]; else if(data[k] < temp_fft_data_min) temp_fft_data_min = data[k]; } //find the maximum ABSOLUTE value temp_fft_data_min = -temp_fft_data_min; if (temp_fft_data_max < temp_fft_data_min) temp_fft_data_max = temp_fft_data_min; //find the shift value if(temp_fft_data_max > 16383) //if the maximum absolute value great than 0.5 { FFTPreProcessCnt = -1; } else if(temp_fft_data_max < 8192) { if(temp_fft_data_max > 512) //if data less than 512, skip pre-processing of the data { while(temp_fft_data_max < 8192) { temp_fft_data_max = temp_fft_data_max << 1; FFTPreProcessCnt++; } } } #endif //load the ADC data into FFTData[] array for (k=0; k<BLOCK_LENGTH; k++) { if(FFTPreProcessCnt < 0) //if absolute value of the data great than 0.5, then divide by 2 FFTData[k].real = data[k] >> 1; else FFTData[k].real = data[k] << FFTPreProcessCnt; FFTData[k].imag = 0; } //FFT conversion FFTComplexIP (LOG2_BLOCK_LENGTH, &FFTData[0], (fractcomplex *) __builtin_psvoffset(&twiddleFactors[0]), (int)__builtin_psvpage(&twiddleFactors[0])); BitReverseComplex (LOG2_BLOCK_LENGTH, &FFTData[0]); //calculate the module of each order harmonic, ComputeMagnitude (&FFTData[0], &FFTMagnitude[0], BLOCK_LENGTH/2); THD_ALL = 0; for(i=2;i<BLOCK_LENGTH/2;i++) //calculate the percentage of each order harmonic { THD[i] = FFTMagnitude[i]*10000.0/FFTMagnitude[1]; //multiply with 100*100 before calculate square root THD[i] = sqrt(THD[i]); THD_ALL += THD[i]; } }