C中任何更快的RMS值计算?
我正在用C编写一个小型8位微控制器的软件。部分代码是读取电流互感器(ZCT)的ADC值,然后计算RMS值。 流过ZCT的电流是正弦曲线但可能会失真。 我的代码如下:
float adc_value, inst_current; float acc_load_current; // accumulator = (I1*I1 + I2*I2 + ... + In*In) double rms_current; // Calculate the real instantanous value from the ADC reading inst_current = (adc_value/1024)*2.5; // 10bit ADC, Voltage ref. 2.5V, so formula is: x=(adc/1024)*2.5V // Update the RMS value with the new instananous value: // Substract 1 sample from the accumulator (sample size is 512, so divide accumulator by 512 and substract it from the accumulator) acc_load_current -= (acc_load_current / 512); inst_current *= inst_current; // square the instantanous current acc_load_current += inst_current; // Add it to the accumulator rms_current = (acc_load_current / 512); // Get the mean square value. (sample size is 512) rms_current = sqrt(rms_current); // Get RMS value // Now the is the real RMS current
但是,它有许多浮点计算。 这给我的小型MCU增加了很大的负担。 我发现sqrt()
函数在我的编译器中不起作用。
有没有可以运行得更快的代码?
当您需要在缺少FPU的处理器上加快速度时,最好的办法是用固定点替换浮点计算。 将此与joop的建议(一个Newton-Raphson sqrt)相结合,你可以得到这样的结果:
#define INITIAL 512 /* Initial value of the filter memory. */ #define SAMPLES 512 uint16_t rms_filter(uint16_t sample) { static uint16_t rms = INITIAL; static uint32_t sum_squares = 1UL * SAMPLES * INITIAL * INITIAL; sum_squares -= sum_squares / SAMPLES; sum_squares += (uint32_t) sample * sample; if (rms == 0) rms = 1; /* do not divide by zero */ rms = (rms + sum_squares / SAMPLES / rms) / 2; return rms; }
只需通过此filter运行原始ADC样本即可。 您可以在此处添加一些位移以获得更高的分辨率,但您必须小心不要溢出变量。 我怀疑你真的需要额外的分辨率。
滤波器的输出与其输入的单位相同。 在这种情况下,它是ADC的单位:2.5 V /1024≈2.44mV。 如果您可以将此单元保留在后续计算中,则可以避免不必要的转换,从而节省周期。 如果您确实需要以伏特为单位的值(可能是I / O要求),那么您将必须转换为浮点数。 如果你想要毫伏,你可以保持整数范围:
uint16_t rms_in_mV = rms_filter(raw_sample) * 160000UL >> 16;
由于您的平方和值acc_load_current
在迭代之间变化不大,因此其平方根将几乎不变。 Newton-Raphson sqrt()函数通常仅在几次迭代中收敛。 通过每步使用一次迭代,计算被抹掉。
static double one_step_newton_raphson_sqrt(double val, double hint) { double probe; if (hint <= 0) return val /2; probe = val / hint; return (probe+hint) /2; } static double acc_load_current = 0.0; // accumulator = (I1*I1 + I2*I2 + ... + In*In) static double rms_current = 1.0; float adc_value, inst_current; double tmp_rms_current; // Calculate the real instantanous value from the ADC reading inst_current = (adc_value/1024)*2.5; // 10bit ADC, Voltage ref. 2.5V, so formula is: x=(adc/1024)*2.5V // Update the RMS value with the new instananous value: // Substract 1 sample from the accumulator (sample size is 512, so divide accumulator by 512 and substract it from the accumulator) acc_load_current -= (acc_load_current / 512); inst_current *= inst_current; // square the instantanous current acc_load_current += inst_current; // Add it to the accumulator tmp_rms_current = (acc_load_current / 512); rms_current = one_step_newton_raphson_sqrt(tmp_rms_current, rms_current); // Polish RMS value // Now the is the APPROXIMATE RMS current
笔记:
- 我将一些数据类型从
float
更改为double
(这在通用机器/桌面上是正常的)如果在微型计算机上double
非常昂贵,您可以将其更改回来。 - 我还添加了
static
,因为我不知道代码是来自函数还是来自循环。 - 我将函数设置为
static
以强制它内联。 如果编译器没有内联静态函数,则应手动内联。
希望您的项目用于测量“大”交流电压(而不是像9v电机控制那样。)如果遇到这种情况,那么您可以作弊,因为您的错误可以在可接受的范围内。
以整数表示所有数学,并使用简单的查找映射进行sqrt操作。 (您可以在启动时预先计算,如果您正在进行3阶段,那么您只需要大约600个奇数值。
这也引出了一个问题,你是否真的需要VAC RMS或其他一些功率测量? (例如,你可以用一个简单的盒子来逃避algorythm吗?)
要查找平方根,请使用微芯片中的以下应用说明来获取8位控制器
快速整数平方根
这是非常快,只能在9个循环中找到平方根
-
以2的幂为单位/乘法
可以通过仅通过位掩码操作和
+,-
更改指数来完成+,-
因此掩码/将指数提取为integer
数值然后应用偏差 。 之后add/sub
值log2(operand)
并编码回您的double
值 -
开方
应该多快和准确? 您可以在固定点上使用二进制搜索或使用sqrt(x)= pow(x,0.5)= exp2(0.5 * log2(x)) 。 再次在固定点上,它很容易实现。 您可以通过采用尾数并将其移位到您使用的值周围的某个已知指数来暂时使其成为固定点的两倍+处理偏移量或者如果您有足够的位则为
2^0
…计算
sqrt
然后存储回double
。 如果你不需要太大的精度,那么你可以保持操作数指数并直接在尾数上进行二进制搜索。
[edit1] C ++中的二进制搜索
//--------------------------------------------------------------------------- double f64_sqrt(double x) { const int h=1; // may be platform dependent MSB/LSB order const int l=0; DWORD b; // bit mask int e; // exponent union // semi result { double f; // 64bit floating point DWORD u[2]; // 2x32 bit uint } y; // fabs yf=x; yu[h]&=0x7FFFFFFF; x=yf; // set safe exponent (~ abs half) e=((yu[h]>>20)&0x07FF)-1023; e/=2; // here can use bit shift with copying MSB ... yu[h]=((e+1023)&0x07FF)<<20; // mantisa=0 yu[l] =0x00000000; yu[h]&=0xFFF00000; // correct exponent if (yf*yf>x) { e--; yu[h]=((e+1023)&0x07FF)<<20; } // binary search for (b =0x00080000;b;b>>=1) { yu[h]|=b; if (yf*yf>x) yu[h]^=b; } for (b =0x80000000;b;b>>=1) { yu[l]|=b; if (yf*yf>x) yu[l]^=b; } return yf; } //---------------------------------------------------------------------------
它返回sqrt(abs(x))
,结果与我的C ++ IDE(BDS2006 Turbo C ++)中的“math.h”实现匹配。 指数从x
值的一半开始,如果需要,对值x>1.0
进行校正。 其余的很明显
代码是用C ++编写的,但它仍然没有优化,可以肯定地改进…如果你的平台不知道DWORD
使用unsigned int
代替。 如果您的平台不支持32位整数类型,则将其块化为4 x 16位值或8 x 8位值。 如果您有64位,则使用单个值来加速该过程
不要忘记将指数也处理为11位….所以对于8位寄存器只使用2 …如果你乘以并将mantissas作为整数进行比较,则可以避免FPU操作。 乘法本身也是累积的,因此您可以使用先前的子结果。
[笔记]
对于位位置,请参阅wiki双精度浮点格式