基于STM32 DSP庫函數(shù)的電力諧波分析,輸入時域信號采樣值,進行Q31 FFT計算,在有頻譜泄露、柵欄效應(yīng)的前提下計算基波和諧波的頻率、幅度、相位。函數(shù)void spectrum_with_cfft(q31_t * x_X, const int NPT,const float df,float Magnitude_harmonic[5],float Angle_harmonic[5],float f_harmonic[5])實現(xiàn)頻譜分析的功能,附帶測試信號的生成代碼。
- //基于cfft的頻譜分析,RAM空間復(fù)雜度8*NPT
- //頻譜分析方法文獻:基于加漢寧窗的FFT高精度諧波檢測改進算法
- //reference: STM32的DSP庫FFT函數(shù) https://wenku.baidu.com/view/9cbc1a94eff9aef8951e061b.html
- //輸入?yún)?shù):
- //x_X: of size 2*NPT,
- // 作為輸入的ADC采樣數(shù)據(jù)時,in the format of {adc0,0,adc1,0,adc2,0,........}
- // 計算的頻域數(shù)據(jù)也存儲在x_X中, 格式是{X[0]的實部,X[0]的實部,X[1]的實部,X[1]的虛部,X[2]的實部,X[2]的虛部,......}
- //df: 頻率分辨率,等于采樣頻率除以采樣點數(shù)
- //Magnitude_harmonic: 計算的基波及諧波的電壓有效值
- //Angle_harmonic: 計算的基波及諧波的相位角
- //f_harmonic: 計算的基波及諧波的頻率
- void spectrum_with_cfft(q31_t * x_X, const int NPT,const float df,float Magnitude_harmonic[5],float Angle_harmonic[5],float f_harmonic[5])
- {
- int i;//時域信號的下標(biāo)
- int k;//頻域信號的下標(biāo)
-
- //譜線鄰域分析用的變量
- int k0;//諧波譜線的理論下標(biāo)
- const int harf_neighborhood=8;//正負15Hz鄰域, 2/df
- float mag_neighborhood[harf_neighborhood*2+1];
- float ang_neighborhood[harf_neighborhood*2+1];
- int idx_harmonic;//諧波的下標(biāo)
- float max_mag;
- int ka,km;
- float alpha_m;
- float dkm;
-
- //加直流偏置
- for(i=0;i<NPT;i++)
- {
- //這里因為單片機的ADC只能測正的電壓 所以需要前級加直流偏執(zhí)
- //加入直流偏執(zhí)后,需要在軟件上減去2048即一半,達到負半周期測量的目的(需要根據(jù)具體情況來進行配置)
- //x[2*i] = ((signed short)(adc_buf[i]-2048)) << 16;
- //左移位數(shù)需要小于32-12,以避免定點FFT計算時溢出。
- x_X[2*i] = ((signed short)(x_X[2*i]-2048)) << 8;
- x_X[2*i+1]=0;
- }
-
-
- //數(shù)據(jù)加窗,漢寧窗
- for(i=0;i<NPT;i++)
- x_X[2*i] = x_X[2*i]* 0.5*(1-cos(2*PI*i/(NPT+1)));
-
-
- arm_cfft_q31(&arm_cfft_sR_q31_len2048,x_X,0,1);
-
- //計算的頻率、幅度、相位
- for (idx_harmonic=0;idx_harmonic<5;idx_harmonic++)
- {
- float magnitude_XH5_km,angle_XH5_km;
- if (idx_harmonic==0)
- {
- //基波
- k0=(idx_harmonic+1)*50/df;//諧波的序號的理論值
- //計算諧波鄰域內(nèi)的幅度譜、相位譜
- for (k=k0-harf_neighborhood;k<=k0+harf_neighborhood;k++)
- {
- float Xreal_k=x_X[2*k]/60.-(x_X[2*(k+1)]+x_X[2*(k-1)])/90.+(x_X[2*(k+2)]+x_X[2*(k-2)])/360.;
- float Ximag_k=x_X[2*k+1]/60.-(x_X[2*(k+1)+1]+x_X[2*(k-1)+1])/90.+(x_X[2*(k+2)+1]+x_X[2*(k-2)+1])/360.;
- Xreal_k/=256;//與(adc_value-2048)) << 8 對應(yīng)
- Ximag_k/=256;//與(adc_value-2048)) << 8 對應(yīng)
- mag_neighborhood[k-k0+harf_neighborhood]=sqrt(Xreal_k*Xreal_k+Ximag_k*Ximag_k);//先不乘以2。實際上單邊譜cn是雙邊譜Fn的2倍
- ang_neighborhood[k-k0+harf_neighborhood]=atan2(Ximag_k,Xreal_k)/PI*180;
- }
-
- //計算中心頻率、幅度、相位
- //找到幅度極大譜線ka和臨建次大譜線kb,取km=min(ka,kb)
- max_mag=0;
- for (k=k0-harf_neighborhood;k<=k0+harf_neighborhood;k++)
- {
- if (mag_neighborhood[k-k0+harf_neighborhood]>max_mag)
- {
- ka=k;
- max_mag=mag_neighborhood[k-k0+harf_neighborhood];
- }
- }
- if (mag_neighborhood[ka-1-k0+harf_neighborhood]<mag_neighborhood[ka+1-k0+harf_neighborhood])
- {
- km=ka;
- }
- else
- {
- km=ka-1;
- }
- alpha_m=mag_neighborhood[km-k0+harf_neighborhood]/mag_neighborhood[km+1-k0+harf_neighborhood];
- dkm=(4-3*alpha_m)/(1+alpha_m);
- magnitude_XH5_km=mag_neighborhood[km-k0+harf_neighborhood];
- angle_XH5_km=ang_neighborhood[km-k0+harf_neighborhood];
- }
- else
- {
- float Xreal_km,Ximag_km;
- float temp=f_harmonic[0]/df*(idx_harmonic+1);
- km=floor(temp);
- dkm=temp-km;
- Xreal_km=x_X[2*km]/60.-(x_X[2*(km+1)]+x_X[2*(km-1)])/90.+(x_X[2*(km+2)]+x_X[2*(km-2)])/360.;
- Ximag_km=x_X[2*km+1]/60.-(x_X[2*(km+1)+1]+x_X[2*(km-1)+1])/90.+(x_X[2*(km+2)+1]+x_X[2*(km-2)+1])/360.;
- Xreal_km/=256;//與(adc_value-2048)) << 8 對應(yīng)
- Ximag_km/=256;//與(adc_value-2048)) << 8 對應(yīng)
- magnitude_XH5_km=sqrt(Xreal_km*Xreal_km+Ximag_km*Ximag_km);//先不乘以2。實際上單邊譜cn是雙邊譜Fn的2倍
- angle_XH5_km=atan2(Ximag_km,Xreal_km)/PI*180;
- }
- f_harmonic[idx_harmonic]=(km+dkm)*df;//中心頻率
- //校正諧波參數(shù)
復(fù)制代碼
【必讀】版權(quán)免責(zé)聲明
1、本主題所有言論和內(nèi)容純屬會員個人意見,與本論壇立場無關(guān)。2、本站對所發(fā)內(nèi)容真實性、客觀性、可用性不做任何保證也不負任何責(zé)任,網(wǎng)友之間僅出于學(xué)習(xí)目的進行交流。3、對提供的數(shù)字內(nèi)容不擁有任何權(quán)利,其版權(quán)歸原著者擁有。請勿將該數(shù)字內(nèi)容進行商業(yè)交易、轉(zhuǎn)載等行為,該內(nèi)容只為學(xué)習(xí)所提供,使用后發(fā)生的一切問題與本站無關(guān)。 4、本網(wǎng)站不保證本站提供的下載資源的準(zhǔn)確性、安全性和完整性;同時本網(wǎng)站也不承擔(dān)用戶因使用這些下載資源對自己和他人造成任何形式的損失或傷害。 5、本網(wǎng)站所有軟件和資料均為網(wǎng)友推薦收集整理而來,僅供學(xué)習(xí)用途使用,請務(wù)必下載后兩小時內(nèi)刪除,禁止商用。6、如有侵犯你版權(quán)的,請及時聯(lián)系我們(電子郵箱1370723259@qq.com)指出,本站將立即改正。
|
|