一、全相位FFT的简介
1.全相位傅里叶变换(apFFT)
传统FFT的结果可以通过一些算法实现频谱校正,如全能重心法、比值法等等。但都是基于FFT的结果,精度有限。全相位FFT算法具有初始相位不变和有效防止频谱泄露的特性。apFFT的算法是先对数据进行全相位预处理,然后进行传统的FFT运算。
2.全相位预处理过程如下:
1,构成一个N点的汉宁窗。
2,汉宁窗对自己求卷积,得到2N-1点的卷积窗。
3,求2N-1点的卷积窗的和。
4,将卷积窗的每一项除以卷积窗的和,得到2N-1点的归一化卷积窗。
5,将数据的1:2N-1项和归一化卷积窗相乘,得到加窗的2N-1项。
6,将第1项和N+1项,第2项和N+2项 ... 第N-1项和第2N-1项相加,得到经过全相预处理的N点序列。
接下来就只需要进行FFT,就可以得到apFFT的结果。
二、CubeMX配置
1.ADC同步采样:
使用ADC主从模式采样,ADC1开启DMA传输和ADC中断(ADC2不用),DMA使用normal模式,以Word字节传输,高、低16位分别存放两个adc通道数据。再以TIM3作为触发源,进行采样。
2.其他配置:
如串口重定向、DSP库配置移植可移步四臂西瓜
三、32关键代码
1.adc启动:
__IO uint8_t over_flag=0; //采集完成标志
HAL_TIM_Base_Start(&htim3);
HAL_ADCEx_Calibration_Start(&hadc1, ADC_CALIB_OFFSET_LINEARITY, ADC_SINGLE_ENDED);
2.全相位预处理
float HanNing[Length]={0}; //汉宁窗数组
float fft_conv[2*Length-1]={0}; //FFT卷积结果
void hanning_window(float *w)
{
uint16_t n;
for (n = 0; n < Length; n++)
{
w[n] = 0.5 * (1 - cos(2 * PI * n / Length));
}
}
void convolve(float *w, float *y)
{
int n, m;
for (n = 0; n < 2 * Length - 1; n++) {
y[n] = 0;
for (m = 0; m < Length; m++) {
if (n - m >= 0 && n - m < Length) {
y[n] += w[m] * w[n - m];
}
}
}
}
// 计算序列的和
float sum(float *x, int n)
{
int i;
float s = 0;
for (i = 0; i < n; i++) {
s += x[i];
}
return s;
}
// 归一化序列
void normalize(float *x, int n)
{
float s = sum(x, n);
int i;
for (i = 0; i < n; i++) {
x[i] /= s;
}
}
hanning_window(HanNing);
convolve(HanNing,fft_conv);
normalize(fft_conv,2*Length-1);
//相位计算函数
float Phase_atan(float32_t *inputSignal,uint32_t index)
{
if(inputSignal[2*index+1] >= 0 && inputSignal[2*index] >= 0)
return 0 + atan(inputSignal[2*index+1] / inputSignal[2*index]) / PI * 180;
else if(inputSignal[2*index+1] >= 0 && inputSignal[2*index] <= 0)
return 180 + atan(inputSignal[2*index+1] / inputSignal[2*index]) / PI * 180;
else if(inputSignal[2*index+1] <= 0 && inputSignal[2*index] <= 0)
return 180 + atan(inputSignal[2*index+1] / inputSignal[2*index]) / PI * 180;
else if(inputSignal[2*index+1] <= 0 && inputSignal[2*index] >= 0)
return 360 + atan(inputSignal[2*index+1] / inputSignal[2*index]) / PI * 180;
}
void AP_FFT(void)
{
for (int i = 0; i < 2*Length-1; i++)
{
wave_2[i]=((single[i] >> 16)/65535.0f*2-1.00f)*3.57f-0.139f;
wave_1[i]=((single[i] & 0xffff)/65535.0f*2-1.00f)*3.57f-0.139f;
}
// Printf(wave_1,Length);
for(int i = 0; i < 2*Length-1; i++)
{
wave_1[i] = wave_1[i] * fft_conv[i]; //加窗
wave_2[i] = wave_2[i] * fft_conv[i]; //加窗
}
// Printf(fft_conv,Length*2-1);
for(int i = 0; i < Length-1; i++)
{
wave_1[i] = wave_1[i] + wave_1[Length+i];
wave_2[i] = wave_2[i] + wave_2[Length+i];
}
// Printf(wave_1,Length);
for (int i = 0; i < Length; i++)
{
fft_inputbuf_1[i * 2] = wave_1[i];
fft_inputbuf_1[i * 2 + 1] = 0;
fft_inputbuf_2[i * 2] = wave_2[i];
fft_inputbuf_2[i * 2 + 1] = 0;
}
arm_cfft_radix4_f32(&scfft,fft_inputbuf_1);
arm_cfft_radix4_f32(&scfft,fft_inputbuf_2);
arm_cmplx_mag_f32(fft_inputbuf_1, fft_outputbuf_1, Length);
arm_cmplx_mag_f32(fft_inputbuf_2, fft_outputbuf_2, Length);
fft_outputbuf_1[0] = 0;
fft_outputbuf_2[0] = 0;
for(int h=0;h<Length;h++)
{
fft_outputbuf_1[h]/=Length/2;
fft_outputbuf_2[h]/=Length/2;
}
// Printf(fft_outputbuf_1,Length);
arm_max_f32(fft_outputbuf_1, Length/2, &Fmax_1, &Amax_pos_1); //使用Length会让频率点在后半部分,导致频率计算错误,导致ARR->0
arm_max_f32(fft_outputbuf_2, Length/2, &Fmax_2, &Amax_pos_2); //使用Length会让频率点在后半部分,导致频率计算错误,导致ARR->0
Phase_1=Phase_atan(fft_inputbuf_1,Amax_pos_1);
Phase_1-=0.04f; //根据实际情况补偿,非必要不加
Phase_2=Phase_atan(fft_inputbuf_2,Amax_pos_2);
Phase=Phase_1-Phase_2;
// printf("Phase1=%f Phase_2=%f Phase=%f\r\n",Phase_1,Phase_2,Phase);
if(Phase>180) Phase = -360 + Phase;
else if(Phase<-180) Phase = 360 + Phase;
}
主函数:
int main(void)
{
/* USER CODE BEGIN 1 */
/* USER CODE END 1 */
/* MCU Configuration--------------------------------------------------------*/
/* Reset of all peripherals, Initializes the Flash interface and the Systick. */
HAL_Init();
/* USER CODE BEGIN Init */
/* USER CODE END Init */
/* Configure the system clock */
SystemClock_Config();
/* Configure the peripherals common clocks */
PeriphCommonClock_Config();
/* USER CODE BEGIN SysInit */
/* USER CODE END SysInit */
/* Initialize all configured peripherals */
MX_GPIO_Init();
MX_DMA_Init();
MX_TIM3_Init();
MX_USART1_UART_Init();
MX_TIM1_Init();
MX_DMA2D_Init();
MX_FMC_Init();
MX_LTDC_Init();
MX_ADC2_Init();
MX_ADC1_Init();
/* USER CODE BEGIN 2 */
RetargetInit(&huart1);
arm_cfft_radix4_init_f32(&scfft,Length,0,1);
HAL_TIM_Base_Start(&htim3); //1M频率
HAL_ADCEx_Calibration_Start(&hadc1, ADC_CALIB_OFFSET_LINEARITY, ADC_SINGLE_ENDED);
hanning_window(HanNing);
convolve(HanNing,fft_conv);
normalize(fft_conv,2*Length-1);
/* USER CODE END 2 */
/* Infinite loop */
/* USER CODE BEGIN WHILE */
while (1)
{
HAL_ADCEx_MultiModeStart_DMA(&hadc1,(uint32_t*)single,2*Length-1);
while(!over_flag);
AP_FFT();
// printf("Phase=%f\r\n",Phase);
over_flag=0;
/* USER CODE END WHILE */
/* USER CODE BEGIN 3 */
}
/* USER CODE END 3 */
}
四、与传统FFT对比
全相位处理:
传统FFT:
可以看出经全相位处理过后的相位非常稳定且精度大大提高。
五、结语
本人作为新人博主,会分享自己在备战电赛时的心得体会,希望大家多多支持,评论区留企鹅号,可分享工程。