一种从电磁频谱中提取频率源数量及基频值的方法与流程
本发明涉及一种从电磁干扰发射频谱中提取谐波分量及设备中频率源数量及基频值的方法，属于电磁兼容技术领域。
背景技术：
随着电子技术的进步和电子信息产业的蓬勃发展，电气电子设备数量呈爆炸式增长，这些设备在工作时会产生一些有用或无用的电磁发射，这些发射可能导致其它电子设备性能降级、传输错误、控制失灵，甚至严重威胁人员的健康与安全，危害极大，因此要对设备的电磁干扰进行控制。
对于电磁兼容工程师而言，要对设备的电磁干扰进行控制，需要首先查找电磁干扰源，然后提取其电磁发射特征进行分析，之后再采取相应的电磁兼容性措施。
目前排查干扰源的方法是对设备采取电磁兼容性测试，根据相应标准(如国军标GJB-151A《军用设备和分系统电磁发射和敏感度要求》)中规定的发射极限值(极限值是指根据相应标准规定，被测设备的测试值不应超过的最大值。如GJB-151A中第5.3.15.2极限规定：“电场辐射发射不应超过图RE102-1～RE102-3所示的值。”)，判断设备是否符合电磁兼容性要求。当发现设备的电磁发射超过限值时，电磁兼容工程师一般只能使用个人经验法或近场探头检测法来判断产生干扰的原因。经验法对工程师的专业技能要求较高，且存在很大的不确定性。近场探头检测法的成本很高，耗费时间长，且通常需要对设备的状态进行变更(如打开机箱或断开某些电路板)才能进行，效果并不理想。
在大量的电磁干扰测试中发现，被测设备的电磁发射频谱虽然复杂，但相当一部分的案例中，具有稳定频率的窄带干扰信号及其由于非线性产生的谐波在被测频谱中占据主导地位。这些发射频谱中的超标谐波虽然占据的频谱很宽，但是由稳定的频率源(如晶振或稳定频率的开关电路)产生的，所以称为谐波分量。如果可以根据这些超标的数据，经过分析确定被测设备中频率源的数量及其基频值，便可以确定该类谐波分量的来源，进而可以分析其发射特征，进而采取针对性较强的措施，具有重要的意义。
目前尚未有相关文献报导。
技术实现要素：
本发明技术解决问题：克服现有技术的不足，针对电磁干扰发射频谱中的主要贡献者——谐波分量，提供一种从电磁频谱中提取频率源数量及基频值的方法，使用自相关-方差分析综合法对谐波分量进行分析，提取被测设备中基频数量及基频值的方法。
本发明技术解决方案：一种从电磁频谱中提取频率源数量及基频值的方法，其步骤如下：
首先采集被测设备的发射频谱数据，并使用多分辨率分析法得到其中与频率源紧密相关的谐波分量；然后计算谐波分量的自相关序列，从自相关序列中提取出备选基频值，同时将谐波分量中该备选基频值对应的谐波置0；然后使用方差分析法检验该备选基频值的有效性：如有效，则将结果记录为有效基频值，否则将其删除；然后重新计算谐波分量的自相关序列，重复上述步骤直到确定不再包含新的频率源为止；最后剔除有效基频值中重复的结果后计数并输出，若有效基频值与某个元器件的工作频率相等，则可以找到频谱数据中谐波产生的根源，可以指导被测设备电磁兼容性设计和整改。
本发明使用自相关-方差分析综合法对谐波分量进行分析，提取被测设备中基频数量及基频值的过程包括下列步骤：
步骤一：读取并显示测数据，选择待分析的频段
使用测试频率序列AFtest＝{f0(1),f0(2),...,f0(n0)}为横坐标，测试幅度值序列AVtest＝{vt0(1),vt0(2),…,vt0(n0)}为纵坐标，做出被测设备的幅度频谱图。同时，根据被测设备的特征，调用该被测设备适用的极限值。
之后，根据测试结果与标准极限值的关系，选择关心的频段(如超标频段)进行分析。选择的频段频率序列为Ftest＝{f(1),f(2),...,f(n)}，幅度值序列为Vtest＝{vt(1),vt(2),...,vt(n)}。
步骤二：剔除宽带分量，提取测试数据中谐波分量Vharm
观察绝大多数电磁干扰测试数据，发现数据基本由两部分组成，其一即为所谓的谐波分量，是指测试数据中具有稳定频率的窄带干扰信号及其由于非线性产生的部分。其二，电磁干扰频谱整体具有一定的趋势，这种趋势常常是由于干扰设备中存在宽带的器件引起的，故称为“宽带分量”。
本发明使用基于小波变换的多分辨率分析法将数据中与频率源无关的宽带分量剔除，得到频谱数据中与频率源紧密相关的谐波分量。
步骤201：获取频谱数据的宽带分量
对测试结果Vtest进行离散小波变换，得到测试结果Vtest在相应分辨率下的近似分量和细节分量。对第5层的近似分量进行重建，得到测试结果的宽带分量：Vtrend＝{vtr(1),vtr(2),...,vtr(n)}。
步骤202：获取频谱数据的谐波分量
将测试结果序列的宽带分量Vtrend从原测试结果序列Vtest中剔除，得到频谱数据中与频率源紧密相关的谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}，即进行如下处理：
vh(i)＝vt(i)-vtr(i) (1)
Vharm＝{vh(1),vh(2),...,vh(n)}为一个序列，Vtrend＝{vtr(1),vtr(2),...,vtr(n)}和Vtest＝{vt(1),vt(2),...,vt(n)}同样为序列，这里是针对Vtest＝{vt(1),vt(2),...,vt(n)}序列中的每一个元素和Vtrend＝{vtr(1),vtr(2),...,vtr(n)}序列中的每个元素进行相减处理。
步骤三：确定频谱数据中的所有频率源的基频值
为了存放频率源的基频值，定义数组fBase(n_Base)用于存放获取的频率源的基频值，其初始值为0。而n_Base为整型变量，其初始值为1。
本发明使用一种结合了自相关变换和方差分析检验的方法，确定对频谱数据中的频率源的数量及所有频率源的基频值。
步骤301：将谐波分量Vharm进行自相关变换
对于谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}，使用公式(2)求其自相关函数为：
τ为自相关函数的自变量。
步骤302：使用自相关函数提取备选基频值
求自相关函数R(τ)除R(0)之外的最大值点的位置，为了便于理解，假设该最大值出现在τ＝P处。实际测试中，两个测试频点之间的频率间隔为Δf＝f(2)-f(1)，则可以使用式(3)确定其最可能存在的备选基波值为：
fB＝Δf·P (3)
公式(3)中，fB表示频率源的备选基频值，Δf为待处理频段中两个测试频点间的频率间隔，P为自相关函数中除R(0)之外的最大值点的位置。
步骤303：对备选基频值的谐波进行置零处理
对谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}的不同值，按照式(4)进行如下处理：
步骤304：使用方差分析法对备选基频值进行有效性检验
将R(τ)按P的间隔进行分组，得到m组不同的数据(矩阵的每一列视作一组)，如式(5)所示：
公式(5)中，xji表示矩阵X中第j行第i列的数据，且xjs＝R((j-1)·P+i)，m是R(τ)中的周期数量，满足m×P≤n。
使用公式(6)对矩阵X的每行进行归一化处理，令
yji＝xji/xj1,j＝1,2,...,m；i＝1,2,...,P (6)
yji表示矩阵Y中第j行第i列的数据，矩阵Y称为归一化矩阵。
根据式(7)、(8)、(9)分别求得Y的组内均值组内方差组间方差
使用式(10)构造统计量F：
则统计量F服从自由度为(P-1,n-P)的F-分布。选定显著性水平α(如α＝0.05)，查表或使用计算机计算出F分布表中Fα(P-1,n-P)的值，记为Fα。
若F＞Fα，则说明假设在置信度为1-α下成立，即该备选基频值有效，并继续进行步骤305；否则，说明自相关函数R(τ)中不存在准周期P，该备选基频值为“伪基频”，将其剔除，并继续进行步骤四。
步骤305：将有效基频值写入基频数组中
将步骤304得到的有效基频值fB存入基频数组fBase(n_Base)中，之后变量n_Base自动加1，然后重复进行步骤301、步骤302、步骤303和步骤304，直到步骤304检验发现自相关函数R(τ)中不再存在准周期，进行步骤四。
步骤四：输出频率源数量、基频值并与器件频率匹配确定干扰源
步骤401：对提取结果进行后处理，并输出结果
使用式(11)数组fBase(n_Base)中的任意两个数据进行检验，若式(11)的条件满足，则说明两个数据为重复的，将两个数据合并为两个数据的平均值。
完成上述步骤后，对有效基频值数组进行计数，将数组fBase(n_Base)中数据的数量作为频率源数量输出，将fBase(n_Base)中的每一个数据作为频率源的基频值输出。
步骤402：使用基频值与器件频率匹配，确定干扰源
将频率源的基频值依次与设备中电子元器件的工作频率进行比对：若从测试数据中提取得到的频率源的基频值与设备中已知的频率源频率相吻合，则可以找到频谱数据中谐波产生的根源，并确定不同频率源对整个发射频谱曲线的贡献大小，指导操作人员依次采取干扰抑制措施，直至满足相应标准极限值要求。
整套方法的实施流程框图见附图4所示。
本发明与现有技术相比的优点在于：
(1)本发明根据实际电路产生干扰的机理不同，使用基于小波变换的多分辨率分析的方法，提取得到与频率源紧密相关的谐波分量，将所有的谐波分量的噪声基底统一到了0附近，可以提高针对性、增强数据特征，提高了频率源提取的精度和速度；
(2)本发明综合使用自相关函数法与方差分析的方法，一方面克服了自相关函数法不能对数据中是否存在周期性进行检验的劣势，另一方面也克服了传统方差分析法计算量过大的劣势，可以准确高效地提取得到设备中可能的频率源个数及频率源基频值；
(3)本发明易于编程实现，可以编写为专用软件，读入测试频谱数据后，直接输出频率源个数及频率源基频值，使得排查干扰源的整个过程更加方便快捷。
(4)本发明得到的提取结果可以对应到测试频谱数据中，物理意义明确，可扩展性强，对电磁发射频谱的特性认知具有重要作用。
附图说明
图1为典型电磁干扰测试设备配置图；
图2为本发明的分析步骤流程框图；
图3为电磁干扰测试得到的幅度频谱图、相应标准极限值对比图；
图4为对测试结果进行小波变换的过程图；
图5为测试结果及其宽带分量的对比图；
图6为对谐波分量进行自相关变换得到的自相关函数R(τ)的曲线图。
具体实施方式
下面结合附图和实例对本发明做进一步的详细说明。
典型的试验配置如附图1所示，电磁干扰测试一般可以分为“传导发射测试”及“辐射发射测试”两类，其中传导发射测试使用电流探头、线性阻抗网络及测量接收机作为主要仪器，对电源导线或天线端子上可能的传导发射信号进行测量。而辐射发射测试是指使用不同类型的天线及测量接收机作为主要仪器，对设备可能辐射到空间中的辐射发射信号进行测量。上述测试过程还需要数据记录装置(计算机)，测量接收机分别与接收天线(或线性阻抗网络)和计算机采用有线方式连接。
在本发明中，测试接收机接收到的信息可以通过有线的方式传输到数据记录计算机。计算机的最低配置为CPU 2GHz，内存2GB，硬盘180GB；操作系统为windows 2000/2003/XP。在本发明中使用的计算机中还安装有电磁兼容性检测软件EMC2003。将计算机与测试接收机相连，启动软件后，计算机按照标准规定的频段和方法，控制测量接收机在一定的频率范围内扫描，得到干扰设备产生辐射的频率(单位Hz)，以及每个频率上对应的幅度测量值(根据标准规定，不同类型测试的单位可能为dBμA、dBμV或dBμV/m)。
测量接收机接收到的信息存入计算机的数据库并可导出为纯文本文档。纯文本文档的内容分为两列：左侧列为扫描的不同频点值序列AFtest，单位为Hz；右侧为对应频点的测量幅度值序列AVtest，单位则可能为dBμA、dBμV或dBμV/m中的一种；两列内容之间每行均使用英文半角的逗号(，)相隔。
本发明采用改进的自相关函数法提取电磁频谱中基频数量及频率的方法包括下列步骤，具体的步骤流程框图见图2所示：
步骤一：读取并显示测试数据，选择待分析的频段
步骤101：读取测试数据文件，做出测试数据的幅度频谱图
测试数据文件中的数据可以读取为两个序列：左侧一列的数据可以读取为AFtest＝{f0(1),f0(2),...,f0(i),...,f0(n0)}，右侧一列的数据可以读取为AVtest＝{vt0(1),vt0(2),...,vt0(i),...,vt0(n0)}，两个序列分别包含n0个数据，n0表示测试点的数量。
所述AFtest＝{f0(1),f0(2),...,f0(i),...,f0(n0)}表示测试频率点依次组成的序列，f0(1)表示测试中的第1个频率点的值、f0(2)表示测试中的第2个频率点的值、f0(i)表示测试中的第i个频率点的值、f0(n0)表示测试中的第n0个频率点的值，i表示序列中任意一个点的标识号，n0为测试点的数量。AVtest＝{vt0(1),vt0(2),...,vt0(i),...,vt0(n0)}表示测试幅度值依次组成的序列，与AFtest的序列值一一对应，即vt0(1)表示在频点f0(1)上测得的幅度值、vt0(2)表示在频点f0(2)上测得的幅度值、vt0(i)表示在频点f0(i)上测得的幅度值，vt0(n0)表示在频点f0(n0)上测得的幅度值，i表示测试中任意一个幅度值点的标识号。
AFtest＝{f0(1),f0(2),...,f0(n)}，AFtest为横坐标，AVtest＝{vt0(1),vt0(2),...,vt0(n)}，AVtest为纵坐标，做出被测设备的幅度频谱图。同时，根据被测设备的特征，调用该被测设备适用的相应极限值，如图3所示。
在附图3中，横坐标为频率值(使用对数坐标)，纵坐标为幅度值，实线代表干扰设备测试得到的频谱幅度曲线，虚线代表根据该设备特征确定使用的GJB-151A中规定的CE102的极限值。实线超过虚线的部分即为频谱超标范围。
步骤102：选择待分析的频段
根据测试结果与标准极限值的关系，为降低计算量，选择关心的频段(如超标频段)进行分析。选择的频段频率序列为Ftest＝{f(1),f(2),...,f(i),...,f(n)}，幅度值序列为Vtest＝{vt(1),vt(2),...,vt(i),...,vt(n)}。
所述Ftest＝{f(1),f(2),...,f(i),...,f(n)}为待分析频段中的测试频率点依次组成的序列，f(1)表示频段中的第1个频率点的值、f(2)表示频段中的第2个频率点的值、f(i)表示频段中的第i个频率点的值、f(n)表示频段中的第n个频率点的值，i表示序列中任意一个点的标识号。Vtest＝{vt(1),vt(2),...,vt(i),...,vt(n)}表示对应待分析频段的幅度值依次组成的序列，与Ftest的序列值一一对应，即vt(1)表示在频点f(1)上测得的幅度值、vt(2)表示在频点f(2)上测得的幅度值、vt(i)表示在频点f(i)上测得的幅度值、vt(n)表示在频点f(n)上测得的幅度值，n为待分析频段中频率(幅度)点数量。
步骤二：剔除宽带分量，提取测试数据中谐波分量Vharm
本发明使用基于小波变换的多分辨率分析法将数据中与频率源无关的宽带分量剔除，剩余的即为频谱数据中与频率源紧密相关的谐波分量。
宽带分量是指电磁干扰频谱整体具有的趋势，这种趋势常常是由于干扰设备中存在宽带的器件造成的，故称为宽带分量。
小波分析是一种包含尺度伸缩和时间平移的双参数的函数分析方法。多分辨率分析是使用小波分析把信号在一系列不同层次的空间上进行分解的方法。这种信号分解的能力能将各种交织在一起的不同频率组成的混合信号分解成不相同频率的子信号，进而可以只分析关心的子信号。
步骤201：获取测试数据的宽带分量
图4给出了对测试结果Vtest进行离散小波变换的方法(该方法参考文献Mallat S.G.A theory for multiresolution signal decomposition:the wavelet representation[J].Pattern Analysis and Machine Intelligence，IEEE Transactions on，1989，11(7):674-693.)，即将测试结果序列Vtest，使用经过特定的数字低通和高通滤波器后分解为第一层的近似分量和细节分量，再对第一层的近似分量使用上述方法继续进行分解，得到第二层的近似分量和细节分量，依次类推可以得到测试结果序列Vtest在不同分辨率下的近似分量和细节分量。使用最后分解得到的近似分量进行小波逆变换重建，可以得到原序列中频率较低的部分，该部分为测试结果的宽带分量Vtrend＝{vtr(1),vtr(2),...,vtr(i),...,vtr(n)}。图5给出了测试数据的宽带分量，其中较细的虚线为原幅度频谱图，较粗的实线为其宽带分量。
所述的宽带分量Vtrend＝{vtr(1),vtr(2),...,vtr(i),...,vtr(n)}中vtr(1)表示宽带分量在待分析频段中第1个频率点的值、vtr(2)表示宽带分量在待分析频段中第2个频率点的值、vtr(i)表示宽带分量在待分析频段中第i个频率点的值、vtr(n)表示宽带分量在待分析频段中第n个频率点的值。
步骤202：获取频谱数据的谐波分量
将测试结果序列的宽带分量Vtrend从原测试结果序列Vtest中剔除，可以得到频谱数据中与频率源紧密相关的谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}，即：
vh(i)＝vt(i)-vtr(i)，i＝1,2,3,...,n (12)
所述的公式(12)中，vh(i)表示谐波分量在待分析频段中第i个频率点的值，vt(i)表示测试数据在待分析频段中在第i个频率点的值，vtr(i)表示宽带分量在待分析频段中在第i个频率点的值，i为0～n的正整数。
步骤三：确定频谱数据中的频率源的数量及所有频率源的基频值
为了存放频率源的基频值，定义浮点型数组fBase(n_Base)，其初始值为0。而n_Base为整型变量，其初始值为1。
本发明使用一种结合了自相关变换和方差分析检验的方法，确定对频谱数据中的频率源的数量及所有频率源的基频值。
步骤301：将谐波分量Vharm进行自相关变换
对于谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}，使用公式(13)求其自相关函数为：
公式(13)中，vh(i)表示谐波分量在待分析频段中第i个频率点的值，i表示待分析频段中任意一个幅度点的标识号；表示对vh(i)和vh(i+τ)从i＝1到i＝n依次求和，求和的结果即为自相关函数R(τ)；τ为自相关函数的自变量。
图6为对谐波分量进行自相关变换得到的自相关函数R(τ)的曲线图。该图的横坐标为R(τ)的自变量τ，纵坐标为R(τ)的值。
步骤302：使用自相关函数提取优势基频值
求自相关函数R(τ)中除R(0)之外的最大值点的位置，为了便于理解，假设该最大值出现在τ＝P处。实际测试中，两个测试频点之间存在频率间隔为Δf＝f(2)-f(1)，则可以使用式(14)确定最有可能的备选基波值为：
fB＝Δf·P (14)
所述公式(14)中，fB表示频率源的备选基频值，Δf为待处理频段中两个测试频点间的频率间隔，P为自相关函数中除R(0)之外的最大值点的位置。
步骤303：对备选基频值的谐波进行置零处理
对谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}的不同值，按照式(15)进行如下处理：
所述公式(15)中，等号左侧的vh(i)表示经过处理后的谐波分量在待分析频段中第i个频率点的值，等号右侧表示需要按照i的取值不同来确定vh(i)的取值：当0.99mP＜i＜1.01mP时，vh(i)＝0；当i为其它值时，vh(i)的取值不变。另外，m＝1,2,3,...,[n/P]表示m的取值范围为大于0且不大于n/P的正整数。
步骤304：使用方差分析法对备选基频值进行有效性检验
将R(τ)按P的间隔进行分组，得到m组不同的数据(矩阵的每一列视作一组)，如式(16)所示：
公式(16)中，xji表示矩阵X中第j行第i列的数据，且xjs＝R((j-1)·P+i)，m是R(τ)中的周期数量，满足m×P≤n。
使用公式(17)对矩阵X的每行进行归一化处理，令
yji＝xji/xj1,j＝1,2,...,m；i＝1,2,...,P (17)
yji表示矩阵Y中第j行第i列的数据，矩阵Y称为归一化矩阵。
根据式(18)、(19)、(20)分别求得Y的组内均值组内方差组间方差
式(18)中，表示归一化矩阵Y中第i列的均值，yji表示归一化矩阵Y中第j行第i列的数据，表示对归一化矩阵Y的第i列的数据yji从j＝1到j＝m依次求和，j表示归一化矩阵Y中任意一行的标识号，i表示归一化矩阵Y中任意一列的标识号，m是R(τ)中的周期数量，亦是归一化矩阵Y的行数。
式(19)中，yji表示归一化矩阵Y中第j行第i列的数据，表示归一化矩阵Y中第i列的均值(由式(18)求得)，表示先对从j＝1到j＝m依次求和，再对结果从i＝1到i＝P依次求和，得到的结果称为组内方差j表示归一化矩阵Y中任意一行的标识号，i表示归一化矩阵Y中任意一列的标识号，m是R(τ)中的周期数量，亦是归一化矩阵Y的行数，P表示R(τ)中除R(0)之外的最大值点的位置，亦是归一化矩阵Y的列数。
式(20)中，表示归一化矩阵Y中第i列的均值(由式(18)求得)，表示归一化矩阵Y中所有数据的平均值(由求得)。表示将乘以m后，再对结果从i＝1到i＝P依次求和，得到的结果称为组间方差而表示对从i＝1到i＝P依次求和后再除以P可以得到i表示归一化矩阵Y中任意一列的标识号，m是R(τ)中的周期数量，亦是归一化矩阵Y的行数，P表示R(τ)中除R(0)之外的最大值点的位置，亦是归一化矩阵Y的列数。
使用式(21)构造统计量F：
统计量F服从自由度为(P-1,n-P)的F-分布。选定显著性水平α(如α＝0.05)，查表或使用计算机计算出F分布表中Fα(P-1,n-P)的值，记为Fα。
若F＞Fα，则说明假设在置信度为1-α下成立，即该备选基频值有效，记作有效基频值，继续进行步骤305；否则，说明自相关函数R(τ)中不存在准周期P，该备选基频值为“伪基频”，将其剔除，并继续进行步骤四。
本步骤中，F-分布是统计学中常见的一种样本分布形式。显著性水平α是指使用该分析检验结果是否有效时，可能犯错误的概率为α。而1-α为置信度，表示是指对待定参数估计的统计可信程度。
步骤305：将有效基频值写入基频数组中
将步骤302得到的有效基频值fB存入基频数组fBase(n_Base)中，之后变量n_Base自动加1，然后重复进行步骤301、步骤302、步骤303和步骤304，直到步骤304检验发现自相关函数R(τ)中不再存在准周期，进行步骤四。
步骤四：输出频率源数量、基频值并与器件频率匹配确定干扰源
步骤401：对提取结果进行后处理，并输出结果
使用式(22)对基频数组fBase(n_Base)中的任意两个数据进行检验，若式(22)的条件满足，则说明两个数据为重复的，将两个数据合并为两个数据的平均值，以此反复，直到基频数组fBase(n_Base)中的任何两个数据不再满足式(22)为止。
所述式(22)中，fBase(i)和fBase(j)分别表示基频数组fBase(n_Base)的第i个和第j个数据，且i≠j。
完成上述步骤后，对有效基频值数组进行计数，将数组fBase(n_Base)中数据数量作为频率源数量输出，将fBase(n_Base)中的每一个数据作为频率源的基频值输出。
步骤402：使用基频值与器件频率匹配，确定干扰源
将频率源的基频值依次与设备中电子元器件的工作频率进行比对：若从测试数据中提取得到的频率源的基频值与设备中已知的频率源频率相吻合，则可以找到频谱数据中谐波产生的根源，并确定不同频率源对整个发射频谱曲线的贡献大小，指导操作人员依次采取干扰抑制措施，直至满足相应标准极限值要求。
下面通过一个具体的实施例对本发明提供的使用自相关-方差分析综合法提取电磁频谱中频率源数量及基频值的方法做出进一步阐述。步骤如下：
对被测设备进行GJB-151A/GJB-152A规定的RE102试验，并将试验结果输出为测试数据文件。下面需要按照附图2所示的流程图，对频率源的数量及基频值进行求解说明：
步骤一：将电磁发射测试中得到的某被测设备的测试数据文件中两列数据分别读入MatlabR2010a，通过图形绘制命令将其显示出来，并按照该被测设备的特点，选择GJB-151A中5.3.15.2节RE102限值的基准曲线。显示的结果如图3所示。根据测试结果与极限值的关系，从图3中可以看到30MHz～200MHz全频段曲线超标严重，则选择待分析频段为30MHz～200MHz。
步骤二：使用图4所表示的流程，使用Matlab R2010a编程实现对待分析频段的测试结果Vtest进行小波变换，得到各层的近似分量和细节分量。对第5层的近似分量进行小波逆变换重建，得到图5中显示的宽带分量Vtrend的提取结果，图5中较细的虚线为原幅度频谱图，较粗的实线为其宽带分量。将测试结果序列的宽带分量Vtrend从原测试结果序列Vtest中剔除，得到频谱数据中与频率源紧密相关的谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}。
步骤三：使用公式将谐波分量Vharm＝{vh(1),vh(2),...,vh(n)}进行自相关变换，计算得到其自相关函数R(τ)，得到的结果如图6所示，图6是谐波分量的自相关函数R(τ)的曲线，横轴是自变量τ，纵轴是R(τ)的值。
然后求自相关函数R(τ)中除R(0)之外的最大值点的位置，得到其最大值的位置为τ＝P＝55。
由于两个测试频点之间存在频率间隔为Δf＝f(2)-f(1)＝0.2MHz，则最有可能的备选基波值为：fB＝11MHz。将谐波分量Vharm中fB的所有倍频置0。
将R(τ)按P的间隔进行分组，得到m组不同的数据，记作矩阵X，矩阵X的部分取值如式(23)所示：
对矩阵X的每行进行归一化处理，令yji＝xji/xj1,j＝1,2,...,m；i＝1,2,...,P，得到归一化矩阵Y。分别求Y的组内均值组内方差组间方差
构造统计量得到F＝4.6642，选定显著性水平α＝0.05，使用计算机计算出F分布表中Fα(P-1,n-P)的值为1.4540，由于F＞Fα，说明结果有效，将有效基频值fB存入基频数组fBase(n_Base)中，之后变量n_Base自动加1。
由于结果有效，需要重新开始步骤3，第2次求得备选基频值fB为13MHz，统计量F＝3.5506，Fα＝1.4117，由于F＞Fα，说明结果有效，则重新开始步骤3。第3次求得备选基频值fB为17MHz，统计量F＝3.6127，Fα＝1.3528，由于F＞Fα，说明结果有效，则重新开始步骤3。第4次求得备选基频值fB为2MHz，统计量F＝2.1225，Fα＝2.3928，由于F＜Fα，说明备选基频值无效，则开始步骤四。
步骤四：检验有效基频值数组中无重复数据，则对有效基频值数组进行计数，说明被测设备中有3个频率源，频率源的基频值分别为：11MHz，13MHz和17MHz。将分析得到的基频值与元器件依次与设备中电子元器件的工作频率进行比对，由于频率源的基频值分别为11MHz，13MHz和17MHz。通过对比电路，发现电路中存在11MHz，13MHz和17MHz的晶振，则将电磁干扰源定位在上述三个设备上，可以直接对三个设备采取针对性措施。
提供以上实施例仅仅是为了描述本发明的目的，而并非要限制本发明的范围。本发明的范围由所附权利要求限定。不脱离本发明的精神和原理而做出的各种等同替换和修改，均应涵盖在本发明的范围之内。
- 还没有人留言评论。精彩留言会获得点赞！