📄 fft.c
字号:
// 32FFT.cpp : Defines the entry point for the console application.
#include "fft.h"
const int costab[32] = {16384, 16069, 15136, 13622, 11585, 9102, 6269, 3196, //定标为14位
0, -3196, -6269, -9102,-11585,-13622,-15136,-16069,
-16383,-16069,-15136,-13622,-11585, -9102, -6269, -3196,
0, 3196, 6269, 9102, 11585, 13622, 15136, 16069};
const int sintab[32] = { 0, 3196, 6269, 9102, 11585, 13622, 15136, 16069, //定标为14位
16384, 16069, 15136, 13622, 11585, 9102, 6269, 3196,
0, -3196, -6269, -9102,-11585,-13622,-15136,-16069,
-16383,-16069,-15136,-13622,-11585, -9102, -6269, -3196};
const int p[32] = { 0,16,8,24,4,20,12,28,2,18,10,26,6,22,14,30,
1,17,9,25,5,21,13,29,3,19,11,27,7,23,15,31};
long fr[32],fi[32];
void FFT();
void Calc(int *nVoltage,int *nCurrent )
{
BYTE it;
unsigned int U1,I1;
long P1,Q1;
// int P2,P3,P4,P5,P6,P7,P9,P11,P13,P15;
// int Q2,Q3,Q4,Q5,Q6,Q7,Q9,Q11,Q13,Q15;
// long Ps0,Qs0; //谐波功率
unsigned int U2,U3,U4,U5,U6,U7;
unsigned int I2,I3,I4,I5,I6,I7;
long temp;
for(it=0;it<32;it++)
{
fr[it] = nVoltage[p[it]];
fi[it] = nCurrent[p[it]];
}
FFT();
//基波电压的幅值
temp = ((fr[1] + fr[31]) >> 1) * ((fr[1] + fr[31]) >> 1) + ((fi[31] - fi[1]) >> 1) * ((fi[31] - fi[1]) >> 1);
U1=sqrt(temp);
U_base = U1;
U1R = (fr[1] + fr[31]) >> 1 ; U1I=(fi[31] - fi[1]) >> 1;
// 基波电流的幅值
temp = ((fi[1] + fi[31]) >> 1) * ((fi[1] + fi[31]) >> 1) + ((fr[1] - fr[31]) >> 1) * ((fr[1] - fr[31]) >> 1);
I1 = sqrt(temp);
I_base = I1;
I1R =(fi[1] + fi[31]) >> 1; I1I = (fr[1] - fr[31]) >> 1;
//二次谐波的电量
temp = (( fr[2] + fr[30]) >> 1) * ((fr[2] + fr[30]) >> 1) + ((fi[30] - fi[2]) >> 1) * ((fi[30] - fi[2]) >> 1);
U2 = sqrt(temp);
Un = ((fr[2] + fr[30]) >> 1) * (( fr[2] + fr[30]) >> 1) + ((fi[30] - fi[2]) >> 1) * ((fi[30] - fi[2]) >> 1);
temp = ((fi[2] + fi[30]) >> 1) * ((fi[2] + fi[30]) >> 1) + ((fr[2] - fr[30]) >> 1) * ((fr[2] - fr[30]) >> 1);
I2 = sqrt(temp);
In = ((fi[2] + fi[30]) >> 1) * ((fi[2] + fi[30]) >> 1) + ((fr[2] - fr[30]) >> 1) * ((fr[2] - fr[30]) >> 1);
//P2 = fr[2] * fi[30] + fi[2] * fr[30] ;
//Q2 = fi[30] * fi[30] + fr[30] * fr[30] - fr[2] * fr[2] - fi[2] * fi[2];
//三次谐波的电量
temp = ((fr[3] + fr[29]) >> 1) * ((fr[3] + fr[29]) >> 1) + ((fi[29] - fi[3]) >> 1) * ((fi[29] - fi[3]) >> 1);
U3 = sqrt(temp);
Un += ((fr[3] + fr[29]) >> 1) * ((fr[3] + fr[29]) >> 1) + ((fi[29] - fi[3]) >> 1) * ((fi[29] - fi[3]) >> 1);
temp = ((fi[3] + fi[29]) >> 1) * ((fi[3] + fi[29]) >> 1) + ((fr[3] - fr[29]) >> 1) * ((fr[3] - fr[29]) >> 1);
I3 = sqrt(temp);
In += ((fi[3] + fi[29]) >> 1) * ((fi[3] + fi[29]) >> 1) + ((fr[3] - fr[29]) >> 1) * ((fr[3] - fr[29]) >> 1);
//P3 = fr[3] * fi[29] + fi[3] * fr[29];
//Q3 = fi[29] * fi[29] + fr[29] * fr[29] - fr[3] * fr[3] - fi[3] * fi[3];
//四次谐波的电量
temp = ((fr[4] + fr[28]) >> 1) * ((fr[4] + fr[28]) >> 1) + ((fi[28] - fi[4]) >> 1) * ((fi[28] - fi[4]) >> 1);
U4 = sqrt(temp);
Un += ((fr[4] + fr[28]) >> 1) * ((fr[4] + fr[28]) >> 1) + ((fi[28] - fi[4]) >> 1) * ((fi[28] - fi[4]) >> 1);
temp = ((fi[4] + fi[28]) >> 1) * ((fi[4] + fi[28]) >> 1) + ((fr[4] - fr[28]) >> 1) * ((fr[4] - fr[28]) >> 1);
I4 = sqrt(temp);
In += ((fi[4] + fi[28]) >> 1) * ((fi[4] + fi[28]) >> 1) + ((fr[4] - fr[28]) >> 1) * ((fr[4] - fr[28]) >> 1);
// P4 = fr[4] * fi[28] + fi[4] * fr[28] ;
// Q4 = fi[28] * fi[28] + fr[28] * fr[28] - fr[4] * fr[4] - fi[4] * fi[4];
//五次谐波的电量
temp = ((fr[5] + fr[27]) >> 1) * ((fr[5] + fr[27]) >> 1) + ((fi[27] - fi[5]) >> 1) * ((fi[27] - fi[5]) >> 1);
U5 = sqrt(temp);
Un += ((fr[5] + fr[27]) >> 1) * ((fr[5] + fr[27]) >> 1) + ((fi[27] - fi[5]) >> 1) * ((fi[27] - fi[5]) >> 1);
temp = ((fi[5] + fi[27]) >> 1) * ((fi[5] + fi[27]) >> 1) + ((fr[5] - fr[27]) >> 1) * ((fr[5] - fr[27]) >> 1);
I5 = sqrt(temp);
In += ((fi[5] + fi[27]) >> 1) * ((fi[5] + fi[27]) >> 1) + ((fr[5] - fr[27]) >> 1) * ((fr[5] - fr[27]) >> 1);
// P5 = fr[5] * fi[27] + fi[5] * fr[27];
// Q5 = fi[27]*fi[27]+fr[27]*fr[27]-fr[5]*fr[5]-fi[5]*fi[5];
//六次谐波的电量
temp = ((fr[6] + fr[26]) >> 1) * ((fr[6] + fr[26]) >> 1) + ((fi[26] - fi[6]) >> 1) * ((fi[26] - fi[6]) >> 1);
U6 = sqrt(temp);
Un += ((fr[6] + fr[26]) >> 1) * ((fr[6] + fr[26]) >> 1) + ((fi[26] - fi[6]) >> 1) * ((fi[26] - fi[6]) >> 1);
temp = ((fi[6] + fi[26]) >> 1) * ((fi[6] + fi[26]) >> 1) + ((fr[6] - fr[26]) >> 1) * ((fr[6] - fr[26]) >> 1);
I6 = sqrt(temp);
In += ((fi[6] + fi[26]) >> 1) * ((fi[6] + fi[26]) >> 1) + ((fr[6] - fr[26]) >> 1) * ((fr[6] - fr[26]) >> 1);
// P6 = fr[6] * fi[26] + fi[6] * fr[26];
// Q6 = fi[26] * fi[26] + fr[26] * fr[26] - fr[6] * fr[6] - fi[6] * fi[6];
//七次谐波的电量
temp = ((fr[7] + fr[25]) >> 1) * ((fr[7] + fr[25]) >> 1) + ((fi[27] - fi[5]) >> 1) * ((fi[27] - fi[5]) >> 1);
U7 = sqrt(temp);
Un += ((fr[7] + fr[25]) >> 1) * ((fr[7] + fr[25]) >> 1) + ((fi[27] - fi[5]) >> 1) * ((fi[27] - fi[5]) >> 1);
temp = ((fi[7] + fi[25]) >> 1) * ((fi[7] + fi[25]) >> 1) + ((fr[7] - fr[25]) >> 1) * ((fr[7] - fr[25]) >> 1);
I7 = sqrt(temp);
In += ((fi[7] + fi[25]) >> 1) * ((fi[7] + fi[25]) >> 1) + ((fr[7] - fr[25]) >> 1) * ((fr[7] - fr[25]) >> 1);
// P7 = fr[7] * fi[25] + fi[7] * fr[25];
// Q7 = fi[25] * fi[25] + fr[25] * fr[25] - fr[7] * fr[7] - fi[7] * fi[7] ;
// P9 = fr[9] * fi[23] + fi[9] * fr[23];
// Q9 = fi[23] * fi[23] + fr[23] * fr[23] - fr[9] * fr[9] - fi[9] * fi[9] ;
// P13 = fr[13] * fi[19] + fi[13] * fr[19];
// Q13 = fi[19] * fi[19] + fr[19] * fr[19] - fr[13] * fr[13] - fi[13] * fi[13] ;
// P11 = fr[11] * fi[21] + fi[11] * fr[21];
// Q11 = fi[21] * fi[21] + fr[21] * fr[21] - fr[11] * fr[11] - fi[11] * fi[11] ;
// P13 = fr[13] * fi[19] + fi[13] * fr[19];
// Q13 = fi[19] * fi[19] + fr[19] * fr[19] - fr[13] * fr[13] - fi[13] * fi[13] ;
// P15 = fr[15] * fi[17] + fi[15] * fr[17];
// Q15 = fi[17] * fi[17] + fr[17] * fr[17] - fr[15] * fr[15] - fi[15] * fi[15] ;
P1 = fr[1] * fi[31] + fi[1] * fr[31] ;
Q1 = fi[31] * fi[31] + fr[31] * fr[31] - fr[1] * fr[1] - fi[1] * fi[1] ;
P1 +=fr[2] * fi[30] + fi[2] * fr[30] ; //二次谐波
Q1 += fi[30] * fi[30] + fr[30] * fr[30] - fr[2] * fr[2] - fi[2] * fi[2];
P1 += fr[3] * fi[29] + fi[3] * fr[29]; //三次谐波
Q1 += fi[29] * fi[29] + fr[29] * fr[29] - fr[3] * fr[3] - fi[3] * fi[3];
P1 +=fr[4] * fi[28] + fi[4] * fr[28] ;//四次谐波
Q1 += fi[28] * fi[28] + fr[28] * fr[28] - fr[4] * fr[4] - fi[4] * fi[4];
P1 += fr[5] * fi[27] + fi[5] * fr[27];//五次谐波
Q1 += fi[27]*fi[27]+fr[27]*fr[27]-fr[5]*fr[5]-fi[5]*fi[5];
P1 +=fr[6] * fi[26] + fi[6] * fr[26];//六次谐波
Q1 += fi[26] * fi[26] + fr[26] * fr[26] - fr[6] * fr[6] - fi[6] * fi[6];
P1 += fr[7] * fi[25] + fi[7] * fr[25];//七次谐波
Q1 += fi[25] * fi[25] + fr[25] * fr[25] - fr[7] * fr[7] - fi[7] * fi[7] ;
P1 += fr[9] * fi[23] + fi[9] * fr[23];//9次谐波
Q1 += fi[23] * fi[23] + fr[23] * fr[23] - fr[9] * fr[9] - fi[9] * fi[9] ;
P1 += fr[11] * fi[21] + fi[11] * fr[21];//11次谐波
Q1 += fi[21] * fi[21] + fr[21] * fr[21] - fr[11] * fr[11] - fi[11] * fi[11] ;
P1 += fr[13] * fi[19] + fi[13] * fr[19];//13次谐波
Q1 += fi[19] * fi[19] + fr[19] * fr[19] - fr[13] * fr[13] - fi[13] * fi[13] ;
//含有谐波的各个基本电度量
Un = sqrt( Un );
In = sqrt( In );
// Ps = (P1 + P2 + P3 + P4 + P5 + P6 + P7 + P9 + P11 + P13 + P15 ) >> 1;
Ps = P1>> 1;
P = P1 >>5; //(P1 + P2 + P3 + P4 + P5 + P6 + P7 + P9 + P11 + P13 + P15 ) >> 1;
// Qs = (Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q9 + Q11 + Q13 + Q15 ) >> 1;
Qs = Q1>> 1;
Q = Q1 >> 6; //(Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q9 + Q11 + Q13 + Q15 ) >> 1;
temp=(unsigned long)U1*U1 + (unsigned long)U2 * U2 +(unsigned long) U3 * U3 + (unsigned long) U4 * U4
+ (unsigned long)U5 * U5 +(unsigned long) U6 * U6 + (unsigned long)U7 * U7;
U = sqrt(temp);
temp =(unsigned long) I1 * I1 +(unsigned long) I2 * I2 +(unsigned long)I3 * I3 + (unsigned long) I4 * I4
+ (unsigned long)I5 * I5 + (unsigned long)I6 * I6 + (unsigned long)I7 * I7;
// I = sqrt(temp);
I=I1;
if(I1<5) //(I1<3)
{
I=0;
Ps=0;
Qs=0;
}
return ;
}
void FFT()
{
int it,m,j,nv,l0;
long temp_r,temp_i;
int pcos=0,s=0,d=0;
m=32; nv=1;
for (l0=5-1; l0>=0; l0--)
{
m >>= 1; nv <<= 1;
for (it=0; it<=(m-1)*nv; it=it+nv)
{
pcos = 0;
for (j=0; j<=(nv >> 1)-1; j++)
{
d = it + j + (nv >> 1);
s = it + j;
temp_r=((long)(costab[m*j]*fr[d]) - (long)(sintab[m*j]*fi[d])) >> 14; //定标为14位
temp_i=((long)(costab[m*j]*fi[d]) + (long)(sintab[m*j]*fr[d])) >> 14;
fr[d]=(fr[s]-temp_r) >> 1;
fi[d]=(fi[s]-temp_i) >> 1;
fr[s]=(fr[s]+temp_r) >> 1;
fi[s]=(fi[s]+temp_i) >> 1;
pcos += m;
}
}
}
return;
}
void CalUnballrate(int i)
{
int temp_r;
int temp_i;
unsigned long ltemp;
if(i==0)
{
temp_r = Ua_base_r - ((Ub_base_r + Uc_base_r) >> 1)
+((Uc_base_i - Ub_base_i) >> 1) + ((Uc_base_i - Ub_base_i) >> 2) + ((Uc_base_i - Ub_base_i) >> 3) - ((Uc_base_i - Ub_base_i) >> 7);
temp_i = Ua_base_i - ((Ub_base_i + Uc_base_i) >> 1)
+((Ub_base_r - Uc_base_r) >> 1) + ((Ub_base_r - Uc_base_r) >> 2) + ((Ub_base_r - Uc_base_r) >> 3) - ((Ub_base_r - Uc_base_r) >> 7);
ltemp = (unsigned long) temp_r * temp_r + (unsigned long) temp_i * temp_i;
ltemp = sqrt(ltemp);
if(ltemp <100)
{
BD_data.UnbalRateofU=0;
}
else
{
temp_r = Ua_base_r - ((Ub_base_r + Uc_base_r) >> 1)
+((Ub_base_i - Uc_base_i) >> 1) + ((Ub_base_i - Uc_base_i) >> 2) + ((Ub_base_i - Uc_base_i) >> 3) - ((Ub_base_i - Uc_base_i) >> 7);
temp_i = Ua_base_i - ((Ub_base_i + Uc_base_i) >> 1)
+((Uc_base_r - Ub_base_r )>> 1) + ((Uc_base_r - Ub_base_r) >> 2) + ((Uc_base_r - Ub_base_r) >> 3) - ((Uc_base_r - Ub_base_r) >> 7);
BD_data.UnbalRateofU = 1000 * sqrt((unsigned long ) temp_r * temp_r + (unsigned long) temp_i * temp_i)/ltemp;
if(BD_data.UnbalRateofU >1000) BD_data.UnbalRateofU =1000;
}
}
if(i==1)
{
temp_r = Ia_base_r - ((Ib_base_r + Ic_base_r) >> 1)
+((Ic_base_i - Ib_base_i) >> 1) + ((Ic_base_i - Ib_base_i) >> 2) + ((Ic_base_i - Ib_base_i) >> 3) - ((Ic_base_i - Ib_base_i) >> 7);
temp_i = Ia_base_i - ((Ib_base_i + Ic_base_i) >> 1)
+((Ib_base_r - Ic_base_r) >> 1) + ((Ib_base_r - Ic_base_r) >> 2) + ((Ib_base_r - Ic_base_r) >> 3) - ((Ib_base_r - Ic_base_r) >> 7);
ltemp = (unsigned long) temp_r * temp_r + (unsigned long) temp_i * temp_i;
ltemp = sqrt(ltemp);
if(ltemp <70)
{
BD_data.UnbalRateofI=0;
//return;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -