matlab如何仿真递推型dft算法,数字化继电器中递推DFT算法的几点实际问题.doc

数字化继电器中

递推DFT算法的几点实际问题

作者：Hatem A. Darwish&Magdy Fikri

翻译：测控06-3班 张维

摘要：在本文中作者将会呈现数字化继电中可逆DFT算法的现实问题，并着重于DFT滤波器的不稳定性导致的放大错误。在估计幅度和相位放大的判断方面也很重要。本文针对新的放大错误的解决方案进行深入测试，，利用DSP为基础实现结果校正。本文将对继电保护工程师尤其是关注与数字化继电保护工程师有帮助。

关键词：可逆DFT，向量估计，DSP，微机继电保护

1.引言

现代数字化继电器是目前被公认的用于电力系统保护高效化和自动化的核心装置。继电器被设计有包括控制、测量、模拟、干预和事件记录的功能。因此，继电器中的测量数值算法应当尤其被关注，特别是在以微机和DSP实现的强电控制的继电气设计中。

DFT滤波器在各种参考书和工业报告中都被广泛用于基波和谐波测量。DFT滤波器主要有两种形式：非递推算法和递推算法。与非递推算法相比，运算量的减少使递推算法得到普遍的使用。有部分技术文章已经强调了一下在输入信号中递归算法对于动态干扰和系统干扰有抑制作用。虽然普遍认为两种算法都有相同的表现，但其实只是理论上的。事实是，已经证明，DFT滤波器被用于浮点微处理时会偏离理想结果，实际操作中输入信号的动态范围和频域会产生一系列相位和幅度的放大错误。

98年上半年开始这种现象在很多的电力系统继电保护实验室运用浮点DSP板研究项目中出现。实验结果是，即使是理想的的50HZ正弦输入信号也会产生一个极低频率的放大误差，并随着频率的正负偏离增大，误差也随之增加。此后这个问题虽然被解决，但却缺乏一个合理的证明[3]。在[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’(50)关系为：

[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分钟，随着时间的延长将得到更大的放大误差，参见图1.c和2.c，这些输入时间延长到一小时。从数据2中可注意到相位需要2s时间去完成一次相应于0.5Hz频率偏离的翻转，而相角的放大误差将以估计的相角变形来显示，见图2。而从图3中可以观察到比较严重的相位和量的变形，这个数据的信号输入时间为5小时。

(a)刚输入信号

(b)输入信号10分钟

(c)输入信号1小时

图2. Δf=0.5Hz相角的放大误差实验数据

放大误差的产生过程非常缓慢且非线性，需要花费很长的时间才能观察到有限的可观察的量，尤其当输入信号频率很低时，从图4和图3的对比中可以看出，图4中信号频率为图3信号频率的四分之一，仅为8次每周期，周期时间同为20s，输入信号比图3中信号频率偏离小，为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的递推构成才是这种现象的原因。它的递推成分属于固有成分，不能在计算中消除舍入误差，因此这个误差就被累积。在理论上来说，不管递推还是非递推，这种误差都应当在实时操作的一个周期内被自我修正。但很可惜，这种情况只有在计算机不存在舍入和丢弃数据尾数时才存在。而当存在量化误差时，这种有限周期也被放大。通过设置完美的缩放系数，量化误差在整数运算处理器中可以被消除，因此在此之前这个现象还没有得到关注。相反的，在大的动态范围内，浮点运算说明缩放系数并不被需要。举例来说，TMS320C40的浮点计算最小值为-5.8774817x10-39，最大值为3.4028235x1038，这些小误差(小于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代表的是实时递推模块(1)和(2)。乘法器1和2被用作以计算新信号和旧信号之差与正弦、余弦系数的乘积。而实虚部综合模块则是用描述相量转化的量和角的变化来模拟示波器追踪仿真。在这里很值得一提的是不管在TMS320C40里指定什么样的数据结构[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.

展开阅读全文