在数字继电器中
递推DFT算法的几个实际问题
作者：Hatem A. Darwish&Magdy Fikri
翻译:06-3班测控: 张维
摘要：在本文中，作者将在数字继电器中呈现可逆性DFT算法的实际问题，重点是DFT滤波器的不稳定性导致放大错误。对估计范围和相位放大的判断也很重要。本文对新的放大错误解决方案进行了深入的测试和使用DSP为基础实现结果校正。本文将有助于继电保护工程师，特别是数字继电保护工程师。
关键词：可逆DFT，向量估计，DSP，微机继电保护
1.引言
现代数字继电器被公认为电力系统保护高效自动化的核心设备。继电器设计具有控制、测量、模拟、干预和事件记录等功能。因此，应特别注意继电器中的测量值算法，特别是在微机和DSP在强电控制的继电气设计中。
DFT基波和谐波测量广泛应用于各种参考书和工业报告中。DFT滤波器主要有两种形式：非递推算法和递推算法。与非递推算法相比，运算量的减少使递推算法得到了广泛的应用。一些技术文章强调，输入信号中的递归算法可以抑制动态干扰和系统干扰。这两种算法虽然普遍认为表现相同，但实际上只是理论上的。事实是，已经证明，DFT当过滤器用于浮点微处理时，会偏离理想的结果，实际操作中输入信号的动态范围和频域会产生一系列相位和范围放大错误。
这种现象从98年上半年开始在许多电力系统继电保护实验室使用浮点DSP出现在板研究项目中。实验结果是，即使是理想的50HZ正弦输入信号也会产生极低频率的放大误差，随着正负偏差的增加，误差也会增加。此后，虽然这个问题得到了解决，但缺乏合理的证明[3]。放大错误的原因归因于数学运算错误的积累，类似于递推滤波器中的反馈[5]-[7]。但是，递推DFT等值框图与反馈递推滤波器的内在区别削弱了该证书的有效性。本文重点关注继电安全的错误分析和补救原则。
本文已指明递推DFT几个可能的放大误差实际上解释了这些放大误差的计算和应用建议，并使用等效的Simulink模型证明了实验结果，并提出并评估了六种可选解决方案。
2.递推DFT的数学公式
为了提取输入信号的基波成分(Xj)，一个完整的递推周期可以使用以下公式[1]-[3]DFT滤波器。
其中，N是周期数，Sj和Cj是DFT正弦和余弦的表达部分。输入信号的有效值Fj和相角XX给出以下公式：
相量Fj它固定在复平面上，因此过滤器的估计值简化了电力系统相量之间的相关性。Cj和Sj相应滤波器的正弦系数代表相应的信号，因此，F’j以下公式可计算：
这里i是复数算子，只要注意Xj和Xj-N相等则F’j与F’j-1.相等。如果输入有偏差信号50Hz的偏离量(Δf)根据文献[1]相量，取样始终保持不变F’(50 Δf)与原始F关系如下：
[5]公式可证明F’(50 Δf)变化范围和Δf频域周期有限。此外，还可以观察到ΔF量和相因子ΔΦ变化。关系类型如下:
到目前为止，结果的大部分应用都是为了ΔΦ直接测量ΔF。然而，没有人注意到它产生的有限周期量ΔF，因为它相对较小。在浮点微处理器的实时操作中，这个量会被放大，以降低数字继电器的可靠性。
3.递归DFT不稳定的实验记录
递归DFT使用模拟[8][9]进行测试和使用，将测试设置参数应用于理想输入信号模式的生成DFT类似于报告[10][11]中的测试方法。包括测试设置DSP板DS1003，该板带有I/O输入输出板DS2201[12][13]。之所以选择DS1003是因为它有TMS320C40的DSP芯片，TMS320C40机器周期为6M每秒，支持浮点运算，内置40位累加器、32位寄存器和16位存储地址，其出色的运算容量可概括为[14]-[17]。一般通用I/O板DS2201有20个分布式模拟输入通道，5个12位A/D转换器，8个模拟输出端口，16个通用数字I/O口，5个PWM通道等等。A/D和D/A转换器使用10伏的参考电压，所有输入模拟信号都在DSP程序中数字化并将比例缩小到原始信号的十分之一，即峰值为10伏的输入信号将内部表示为浮点1.0。需要注意的是，浮点操作比整数操作更有可能导致舍入误差，因为整数操作的比例因数会仔细选择以避免舍入操作，而浮点操作则不同。这种现象的细节和整数与浮点数之间的关系将在4部分和5部分中解释。
(a)刚输入信号
(b)输入信号10分钟
(c)输入信号1小时
图1. Δf=0.5Hz量的放大误差实验数据
正弦信号峰值为3.3伏，频率接近49.5Hz。使用此信号后，可以使用图1中的递推DFT记录有限周期的0.5Hz估计，约0.33伏，但这个信号是0.5Hz周期将开始放大。因此，应调整预测量，相量将以一个角度旋转，产生0的比例.5Hz偏差，见图1.a和2.a。此外，误差被放大到可感知的程度，见图1.b和2.b，这些输入信号时间为10分钟，随着时间的延长有更大的放大误差。.c和2.c，这些输入时间延长到一小时。从数据2中可以注意到相位需要2s完成一次对应于0的时间.5Hz翻转频率偏差，以估计的相角变形显示相角放大误差，见图2。数据的信号输入时间为5小时，的严重变形，该数据的信号输入时间为5小时。
(a)刚输入信号
(b)输入信号10分钟
(c)输入信号1小时
图2. Δf=0.5Hz相角的放大误差实验数据
放大误差的过程非常缓慢和非线性，观察有限的可观察量需要很长时间。特别是当输入信号频率很低时，从图4和图3的对比可以看出，图4中的信号频率是图3中信号频率的四分之一，每周期只有8次，周期时间为20次s，输入信号的频率偏号的频率偏差小于0.05Hz偏离。而当0.05Hz当偏差应用于每个周期的32个信号时，误差的速度与图3捕获的速度相同，从而得出误差的速度由递推DFT数学运算的次数决定了频率的偏差f，但有限频率等于f。同样的工艺使用理想的50Hz输入信号和取样时钟也会得到相同的误差。不同时晶和制造工艺产生的数字硬件取样时钟的一般值。虽然误差很小，对不同类型的继电设施影响不大，偏差f也很小，但这种误差会产生错误的放大率。可以清楚地说，任何单一都是基于递推的DFT由于时晶偏差和频率偏差f，输入继电器可以降低可靠性，这可以通过推进DFT解释为什么仍有许多继电设计在偏离频率较大时使用最优化的非递推DFT运算。
(a)量
(b)相角
图3. Δf=0.5Hz 5小时量和相的误差
(a) Δf=0.05Hz的预估量
(b) Δf=0.05Hz的预估相角
图4. 取N=8.16小时后放大误差
和取N=32, 放大误差4小时后
4.放大误差的可能性
一个很明显的产生放大误差的原因要归功于滤波器，详见式子(5)的输入信号频率偏离。有限的周期必然导致频率偏移f，但有限周期量的放大，则根据公式(6)相反变化，因此需要解释。为了找出放大误差的原因，已经做了几次测试。测试表明，无论是从DFT递推公式、采集数据或量化数据都会产生这些误差。由于采集误差A/D产生转换、取样和条件循环；量化误差是为了配合CPU数据尾数的舍入和丢弃是由寄存器和存储器大小引起的。由于非递推，两个过程的共同作用会产生误差DFT没有放大错误，即使是相同的收集和定量设置，相应的误差也是有界变量或XXX。即使这些误差本身不会导致放大误差，它们也会推动DFT构成错误。DFT这就是这种现象的原因。其递推成分属于固有成分，不能在计算中消除舍入误差，因此该误差被积累。理论上，无论是递推还是非递推，这种误差都应在实时操作的一个周期内自行纠正。但遗憾的是，这种情况只存在于计算机没有舍入和丢弃数据尾数的情况下。当存在量化误差时，这一有限周期也被放大。在整数运算处理器中，通过设置完美的缩放系数，可以消除量化误差，因此这种现象在此之前还没有得到重视。相反，浮点动态范围内，浮点操作表明不需要缩放系数。举例来说，TMS320C最小浮点计算值为-5.8774817x10-39，最大值为3.4028235x这些小误差(小于10-39)在实时操作中会继续积累到不可预测的误差程度，因此误差是不可避免的。递推DFT在误差和相关性方面，滤波器的有效模型可以揭示更多信息。
5.递推DFT仿真模型
图5显示用途Matlab中的Simulink工具包建立的等值等值DFT模型，Matlab版本为6.5。文献(1)-(3)也用Simulink实现模型。图中的正弦信号、正弦信号1、正弦信号2是信号发生器的正弦和余弦系数。2047阶定量使用12位A/D转换器 ，延时模块Z-32代表32周期每秒输入时的整个循环时间，但Z-实时递推模块(1)和(2)代表1。乘法器1和2用于计算新信号与旧信号之间的差异和正弦和余弦系数的乘积。实虚部综合模块模拟示波器跟踪模拟，描述相量转换的量和角的变化。这里值得一提的是，无论在哪里TMS320C所有变量都是基于40里指定的数据结构[14]-[17]IEEE标准定义为只有一个浮点值。三个累加模块设置为32位可调输出，输出缩放系数为2*10-23与标准浮点数据格式的尾数相同。需要注意的是，可调浮点数不能直接送入实虚部综合模块。因此，有两个累加模块将可调浮点数转换为固定浮点数。[14]-[17]写在使用手册上CPU和DS1003的原始设置尽量减少数字的舍入运算。
图6表示递推DFT如图5所示，模型的输出信号为49.5Hz幅值为3.3V。从图6可以看出，模拟的递推DFT产生的信号与图1和图2相同，但Simulink模型产生放大误差的速度略慢于实际试验。这种差异的原因可能是两个操作系统的隐性差异，如编程语言和对A/D不同的模型量子化概括。最后，发现放大误差的产生速度取决于量化环节、可调浮点数的位数和缩放系数，减少任何环节都会加速放大误差的产生。例如，图6.c只要将缩放系数设置为2-20而不是2-23，所示的放大误差可以延迟一个小时，以保持相同的位数和量化阶。因此，图5所示的模型可以有效地研究和放大误差Simulink将模型输出设置为单个浮点值，产生32位调浮点数一样的放大误差，但前者的产生速度远低于后者。目前还只有Simulink和DS1003的数据格式与近似A/D模型的差异能解释这种现象。应注意到的是，模拟的结果也会受到Matlab版本不同和PC差异的影响，本实验所用的PC使用的是奔四1.7GHz内核。需要特别强调是，这个讨论的主题并不在Simulink模型和试验产生放大误差的速度上，但不管是因为数据结构(整数、长整数、小数、正负数)、编程语言、CPU还是A/D模型，这个误差确实被放大了。只要递推DFT公式不能稳定，放大误差将一直存在。
图5. 递推DFT的Simulink模型
(a)刚输入信号
(b)输入信号一小时
(C)输入信号七小时
图6. 用Simulink模拟的放大误差的量和相角
(a)微分模块
(c)Z平面的极点和零点
(b)积分模块
图7.递推DFT的积分模块和微分模块
6.递推DFT的不稳定性
为了探讨递推DFT公式的不稳定性，我们使用现有的图5来进行一个适合的说明。图5可以被看成像图7所示的串联微分器和串联积分器，其中积分模块有一个统一的反馈。因此，综合累加模块就不需要对数据进行截取。另一方面，如图7所示，零点和极点落在Z平面控制圆的边缘，如果这两个模块都运行良好，那么它们将不会缠上误差，DFT的输出也会比较稳定。但是，小的舍入误差使极点偏离控制圆的边缘而靠向不稳定的外侧。因此，不管什么数据格式(整数、小数、正负数等)，所有舍入误差都会因为递推DFT的不稳定性而被放大。
已经确定微分模块产生误差而积分模块使误差积累，而使用出色的浮点DSP芯片TMS320C40可以减少这一误差，因为它有较大的内存空间、寄存器空间和累加器(分别为16位、32位和40位)。但即使如此，预估的量还是因为放大误差影响而偏移，只不过使用小容量处理器会使情况更糟。除非有另一个100%Hz的信号且不存在敏感放大，否则f的有限周期被放大，正负号取决于信号是否超过标准频率50Hz。而积分器也会导致额外的fHz信号的放大，因为低频信号(f)获得较高增益而较高频信号(100 f)获得较低增益。
综上所述，可以明确有限周期的递推DFT和如图8所示带反馈的递推滤波器不同源且应用场合不同。递推DFT模块个表格与图8所示滤波器的最基本区别就是不同的反馈衰减。与递推DFT不同，待反馈的递推滤波器的极点被设计在稳定区域内，但为了配合累加环节的位数而设计了对衰减反馈的数据截取，参考文献[5]-[7]。这与最初为递推DFT放大误差设计的应用情况有很大的不同。
为了减少有限周期对放大的影响，可以采取六种解决方案：1.用固定的取样时钟；2.应用特殊滤波器；3.改变数据舍入方式；4.对整数变量的递推DFT运算进行处理；5.平行的DFT运算法则；6.最优化的DFT运算法则。
(a)Q=0.99999999，信号输入6小时后
(b)Q=0.99，刚输入信号时
(c)Q=0.99999，信号输入6小时后
图9. Q的取值对递推DFT输出的影响
7. 减小误差的推荐方法
为了减小这些误差，可以使用六种方法。每个都经过严格的审查并且对其可用性都基于对数字继电工程师的立场的观点进行了测试。每个方法都将在下面进行分别说明。
A. 应用适合的取样时钟
事实上，完美的50Hz信号用调成基频的递推DFT滤波器不会产生任何的放大误差。使用合适的取样时钟来配合频率偏移，每周期就可获得一个稳定数量的取样，这跟完美50Hz的输出信号相似。而只要适合的取样计划没有发现新的频率并进行自身的时间间隔调整，误差就只会转换的时候产生。但是，这个方法很麻烦而且并不能根除掉放大误差，因此放大误差依然以较低的速度产生。
B.应用渗漏滤波器
应用这个滤波器主要基于上面6部分的说明。因为递推DFT并不稳定，可以通过将Sj-1值与Cj-1值乘以一个比单位数小但很接近的因子Q而使积分模块的极点稍微像稳定圈内移动。所以(1)式和(2)式可以重新表示为：
修改过的递推DFT公式是图5中用Q代替单位反馈增益而得到的。从式子(7)和(8)可以发现，极点被移到稳定区域内也因此放大误差得到消除。直接将Q取值为0.99999999会对DFT的输出产生细微的影响，而放大误差依然存在，并随着时间不断积累，如图9.a。而取值0.99则会迫使输出低于3.33V，如图9.b所示，因此Q的取值是很有内涵的。经过几次尝试，发现Q取值为0.99999可以在不影响DFT输出地情况下尽量消除放大误差，如图9.c。这个结果已经经过DSP试验验证。不仅积分模块，对微分模块一起进行修改可以得到更大的改善。
尽管如此，使用渗漏滤波器还不是一个很实际的结局方案，因为它操作不方便，而且也受数限和CPU容量的限制。
C.改变数据舍入的描述
从小数格式转换到整数格式的数学操作通常有四个不同的表达方式，它们是向下舍入，更改位数，向上舍入，向最近数舍入。这些舍入表述在文献[14]中有详细说明。TMS320C40的编程语言的默认舍入方式是向下舍入，即向负无穷舍入；十进制数5.5和-5.5会分别被约成5和-6，而32位的浮点舍入运算与之相似，参见文献[14]。所以，运算过程中会产生小的误差，而递推DFT公式则对误差进行积累。不同的舍入方式也可以减小这种影响，例如更改位数。例如，用这种方法对5.5和-5.5进行舍入，就可以直接得到5和-5，也就是说，数据的细微尾数会直接被省略。事实上，当如图5的系统模型在三个累加模块使用更改位数法进行数据舍入，放大误差就会在许多情况下消失。应用向最近数舍入的方式也可以得到相似的实验结果。因此，发展研究工程师应该在递推运算涉及更改位数的舍入时保持谨慎。
但问题是，以上所说的四种舍入方法在DS1003中使用C语言和汇编语言实现的时候只能引用与浮点数向整形数的数据转换，这在处理浮点数与浮点数的数学住转换中制造了小数差别的麻烦。在浮点数运算中唯一适用的是向下舍入方法，因为它能包容放大误差。而在Simulink环境下的浮点数运算的舍入方法实际上很缺乏，现有的方法没办法与通用微控制器设备配套。幸运的是，更改位数法和向最近数舍入法可以被微控制器和DSP设备使用，虽然只是极少数。基于这种新技术的继电器可以尽可能的减少放大误差的发生。
或者使用这种整型运算递推DFT也是一个可选择的方案，因为它可以在递推DFT运算中省略掉舍入的步骤。
D.使用整型运算的递推DFT
使用这种DFT运算，它的DFT系数、输入信号样本和DFT数学运算都必须使用整型数据。设计者必须保持谨慎所有的运算结果都没有超过整型数的上限，否则会产生凑整和数据截取，所以缩放系数就必须被小心确定来避免溢出问题。这种方法可以避免放大误差和用有界波来描述取样误差。在试验中设置使用整型数的递推DFT可以产生一个稳定的量和相角，如图10所示。这个结果揭示了为什么继电器设计者都不太关注这个问题，因为到目前为止的大部分数字继电器都使用整形运算CPU，并且其缩放系数的变量都经过谨慎的推敲。
这个解决方法在面对浮点运算CPU可能没办法起到作用。因此，需要另一种平行技术来处理浮点CPU的运算。这种技术包含稳定的递推DFT运算，可以服务于浮点型和固定型CPU。
E.平行DFT技术
通过下面步骤可以将非递推和递推DFT公式平行化：
1)写出完整周期内的递推DFT公式
2)但非递推DFT的两个公式(正弦和余弦)并不在一个周期内，公式相加的周期包含两个独立的正弦和余弦波，也就是说非递推数学运算超过一个周期。
3)在周期的最后，递推DFT的正弦和余弦公式都用非递推DFT公式更新，所以，非递推DFT公式值为零。因此，这个周期内的计算的量可以基于递推的形式而在周期末尾使用非递推形式。这个最终的量的误差如图10所示。但是，数据运算的量跟非递推DFT的分布式运算相比要多。
F.优化的平行技术
鉴于以上对平行DFT技术运算量的缩减要求，(1)式和(2)式可以修改为：
这里Pj和Mj分别等于(2/N)cos(2j/N)Xj和(2/N)sin(2j/N)Xj。因此，周期信号通过相应的正弦和余弦系数被放大。这里只需要进行加法的非递推DFT的运算，因为非递推DFT的成分Pj和Mj已经有预估。优化后的平行DFT输出达到3.3V的电压和49.5Hz频率，如图10所示，图中是输入信号5小时测量。很明显，有限周期放大误差的发生几率已经降低，换做更大频率偏移的输入信号也已被证实。即使是随机轴的存储器错误优化DFT的运算法则也使测量更准确，因为误差在一个周期内已减小。
(a)信号输入10小时后的输出量
(b)信号输入10小时后的输出相角
图10. 使用方案D、E和F后的预估输出
从图10中可以看出，三个方法对放大误差有效，分别为E、D和F，开发工程师可以任意挑选一种合适的方法。在这些方法当中，本文最推荐优化的平行DFT来设计继电保护。
8. 推荐方式所需要的时间
每个方案所需要的时间和递推DFT及非递推DFT之间对比很值得一提。用DSP硬件设置分别测试各个方案，并记录下相应的执行时间。表格1显示的是当N等于32波每周期时的各个方案对应结果。可注意到，模拟通道的获得和内置功能执行外围设备都需要将近35微秒。因此，每个方案的纯工作时间都要需要减掉这个值。从表格1可以得到，递推DFT本质上比非递推DFT需要时间更少。方案E平行DFT并不比非递推DFT的内置运算耗费更多时间。最后，方案F优化平行DFT的执行时间几乎与常规递推DFT执行时间相同，但在不同的输入信号和数据格式下输出更稳定。
技术
执行时间
纯执行时间
非递推DFT
73μs
38μs
递推DFT
48μs
13μs
方案B
49μs
14μs
方案D
49μs
14μs
方案E
50μs
15μs
方案F
49μs
14μs
9. 结论
在这篇文章中分析了可能的递推DFT放大误差，并且介绍的这些放大误差的捕捉方法和发生情况。实验结果都用Simulink模型证实，并以此提出并验证了6中可供选择的减少放大误差的方案，其中推荐优化DFT。实验结果表明，随着新的应用浮点微控制器和DSP的硬件的出现，推荐的优化平行DFT应当被现代计算机继电保护所考虑。
10.参考文献
[1] A.G. Phadke, J.S. Thorp and M.G. Admiak, “A New Measurement
Technique for Tracking Voltage Phasors, Local System Frequency, and
Rate of Change of Frequency,” IEEE Trans. on Power Apparatus and
System, Vol. PAS 102, No. 5, May 1983
[2] A.G. Phadke (Chairman), W.G. H7 of the Relaying Channels, “Synchronized
Sampling and Phasor Measurements for Relaying and Control,”
IEEE Trans. on Power Delivery, Vol. 9, No. 1, Jan. 1994.
[3] H.A. Darwish, “A Close Accord on DFT Based Frequency and Phasor
Estimators for Numerical Relays,” Bulletin of Electric Research Journal
(ERJ), Shebin El-Kom, Egypt, No. 63, Oct. 2001.
[4] B. Kasztenny and EE. Rosolowski, “Two New Measuring Algorithms
for Generator and Transformer Protection,” IEEE on Power Delivery,
Vol. 13, No. 4, Oct. 1998.
[5] T. Laasko and I. Hartimo, “Noise Reduction in Recursive Digital Filters
Using High-Order Error Feedback,” IEEE Transaction on Signal
Processing, Vol. 40, No.4, May 1992.
[6] D. Williamson and S. Sridharan, “Residue Feedback in Digital Filters
Using Fractional Feedback Coefficients,” IEEE Transaction on Acoustics,
Speech, Signal Processing, Vol. ASSP-33, April 1985.
[7] T. Laakso, P. S. R. Diniz, I. Hartimo, and T. C. Macedo Jr., "Elimination
of Zero-Input and Constant-Input Limit Cycles in Single-Quantizer
Recursive Filter Structures," IEEE Transactions on Circuits and Systems,
Vol. 39, No 9, Sept. 1992.
[8] Relay Performance Testing, IEEE Power System Relaying Committee
Report, Special Publication No. 96 TP 115-0, 1996.
[9] IEEE Work group F-8 of the Relay Input Sources Subcommittee of the
IEEE Power System Relaying Committee, “Digital Simulator Performance
Requirements for Relay Testing,” IEEE Trans. on Power Delivery,
Vol. 13, No. 1, Jan. 1998.
[10] H.A. Darwish, A.I. Taalab, and T.A Kawady, "Development and Implementation
of an ANN-Based Fault Diagnosis Scheme for Generator
Winding Protection,” IEEE Trans. on Power Delivery, Vol. 16, No.2,
April 2001.
[11] H.A. Darwish, A.I. Taalab, and H. Assal, “A novel Overcurrent Relay
with Universal Characteristics,” IEEE PES Conference on Transmission
and Distribution, T&D2001, Atlanta, Oct. 28-Nov. 1, 2001.
[12] Floating-Point Processor Board DS1003, User’s Guide, dSPACE, Germany,
1996.
[13] Multi Input/OutputDS2201, User’s Guide, dSPACE, Germany, 1996.
[14] TMS320C3x/4x Floating-Point DSP Chip, User’s Guide, Texas Instruments,
USA, 1994.
[15] TMS320 Floating-Point DSP Assembly Language Tools, User’s Guide,
No. SPRU035B, Texas Instruments, USA, February 1995.
[16] TMS320 Floating-Point DSP Optimizing C-Compiler, User’s Guide,
No. SPRU034F, Texas Instruments, USA, February 1995.
[17] TMS320C4x General Purpose Applications, User’s Guide, No.
SPRU159A, Texas Instruments, USA, May 1999.
[18] Combined Phase and Earth Fault Overcurrent Relay (AU1S), The Engineering
Office for Integrated Projects (EOIP), Maadi, Cairo, Jan.
2002. Revised Jan 2005.
[19] Matlab and Simulink Software Package, Ver. 6.5, Mathworks, 2002.
展开阅读全文