本帖最后由 highgear 于 2010-8-30 08:28 编辑
(5) 问题 A 的解答
在上面的代码 CalculateHarmonic(Complex* X, int harmonic) 中可以看出dft 的各次谐波计算是独立的， 不依赖其它次谐波。而且，问题 A 不需要计算2次(100hz)，3次(150hz)等等谐波，这是 dft 的优点之一。首先，定义两个复数的结构：
typedef short int16;
typedef int int32;
typedef struct SComplex
{
int16 R;
int16 I;
} Complex;
typedef struct SComplex32
{
int32 R;
int32 I;
} Complex32;
接着， 定义两个常数以及电压电流的结构：
#define N 16 //每周期采样点数
#define LOG2_N 4 // log2(N)
struct UI {
Complex U; //电压的结果
Complex I; //电流的结果
int16 Voltage[N]; //先前的 N 个电压
int16 Current[N]; //先前的 N 个电流
Complex32 UAcc; //电压的累加器
Complex32 IAcc; //电流的累加器
int Index; //迭代索引计数器， 8-BIT MCU 可以为 char, 如果 N < 256.
Complex W[N]; //N 点的 cos, sin 系数
} ui;
初始化，cos, sin 系数数组应该事先计算好：
void UI_Init()
{
for (int i=0; i<N; i++) {
ui.W[i].R = (int16) (::cos( 2*3.1415927*i/N) * (1<<14) + 0.5); //应离线计算！！！
ui.W[i].I = (int16) (::sin(-2*3.1415927*i/N) * (1<<14) + 0.5);
ui.Voltage[i] = 0;
ui.Current[i] = 0;
}
ui.UAcc.R = 0; ui.UAcc.I = 0;
ui.IAcc.R = 0; ui.IAcc.I = 0;
ui.Index = 0;
}
下面的代码不断递推， 可以求出电压和电流的复数：
void UI_Calculate(int16 voltage, int16 current)
{
int16 d;
d = voltage - ui.Voltage[ui.Index];
ui.Voltage[ui.Index] = voltage;
ui.UAcc.R += (d * ui.W[ui.Index].R) >> (13 + LOG2_N - 16);
ui.UAcc.I += (d * ui.W[ui.Index].I) >> (13 + LOG2_N - 16);
ui.U.R = (int16) (ui.UAcc.R >> 16);
ui.U.I = (int16) (ui.UAcc.I >> 16);
d = current - ui.Current[ui.Index];
ui.Current[ui.Index] = current;
ui.IAcc.R += (d * ui.W[ui.Index].R) >> (13 + LOG2_N - 16);
ui.IAcc.I += (d * ui.W[ui.Index].I) >> (13 + LOG2_N - 16);
ui.I.R = (int16) (ui.IAcc.R >> 16);
ui.I.I = (int16) (ui.IAcc.I >> 16);
ui.Index = (ui.Index + 1) & (N-1);
}
上面的计算dft计算使用的是 16-bit, 32-bit 的定点运算，这里需要把电压和电流单位化。比如系统最大电压幅值为 Vmax = 400V，最大电流幅值 Imax = 20A, 在数字系统中统一归一化: Q15 = 2^15 = 32768. 即 Vmax，Imax在数字系统对应Q15 = 32768. 因此，演示主程序中的：
8000 ---〉8000/Q15 * Vmax = 97.66V
4000 ---〉4000/Q15 * Imax = 2.441A
至于功率，很简单, 用电压乘以电流的轭(用 j 来代替复数i, 以免混淆)：
P + jQ = U*I’ = (ui.U.R + j ui.U.I) * (ui.I.R – j ui.I.I)
P是有功功率， Q是无功功率；
功率因数为：cos(theta) = P / sqrt(P*P +Q*Q)
visual c++下的演示主程序 ：
#include "stdafx.h"
#include "Math.h"
#include "stdio.h"
#include "UI.h"
#define Magnitude(c) ((int) sqrtf(c.R*c.R + c.I*c.I))
#define PI 3.14159265f
int _tmain(int argc, _TCHAR* argv[])
{
int16 voltage, current;
Complex PQ;
UI_Init();
for (int i=0; i<1000; i++) {
voltage = (int16) (::sin(2*PI*i/N)*8000); //模拟采样电压
current = (int16) (::sin(2*PI*i/N + PI/3)* 4000); //模拟采样电流
UI_Calculate(voltage, current);
}
PQ.R = (ui.U.R * ui.I.R + ui.U.I * ui.I.I) >> 15;
PQ.I = (ui.U.I * ui.I.R - ui.U.R * ui.I.I) >> 15;
printf("Voltage: %d\n", Magnitude(ui.U));
printf("Cuurent: %d\n", Magnitude(ui.I));
printf("Power Factor: %d\n", PQ.R * 1000 / Magnitude(PQ));
::getchar();
return 0;
}
结果
Voltage: 8000
Cuurent: 3999
Power Factor: 500
|
楼主
