一种基于互相关算法的涡流阻抗求解方法及装置与流程

本发明涉及涡流阻抗技术领域，更具体的说是涉及一种基于互相关算法的涡流阻抗求解方法及装置。

背景技术：

当具有导电性能的工件置于涡流探头磁场时，便会在涡流线圈中产生特定的阻抗值，由相同材料制成的所有试样在探头中便会产生相同的阻抗，据文献记载，其阻抗变化曲线是被测物体性质和检测仪器的仪器特性函数。被测物体性质包括：1)电导率；2)被检物体尺寸；3)磁导率；4)工件缺陷等。仪器特性包括：1)检测线圈中交流场的频率；2)检测线圈的大小和形状；3)检测线圈与被检物体间的距离(提离)等。硬度是基于试样的初始幅值磁导率来检测的，在保证测试线圈的仪器特性、电导率及尺寸标样相同、检测线圈与被检工件的位置固定等条件下，标定出合格曲轴的阻抗值，将待测试样的阻抗值同已知试样的阻抗值进行比较，形成一个公差带，用公差带即可判断待测曲轴硬度区间，公差带范围决定硬度分选的分辨率。这样，基于涡流技术的硬度检测实际上变成一种在试样阻抗值区间进行比较的分选方法，分选的关键是研究出一种可靠性好、准确度高、实时性强的阻抗求解技术。

有学者将机器学习应用到涡流检测特征提取，如葛亮利用支持向量机(SVM)技术实现了窄搭接焊缝涡流信号的特征提取与缺陷识别，缺陷判别准确率达到80％以上，但该算法需要巨大的GPU开销，实时性和成本不能得到较好的保证，对于定性分析来说，准确率仍有欠缺；Sasi B研究了基于线性代数法的双频涡流检测技术，但该方法运算量较大，处理周期较长、实时性较差；范孟豹基于相敏检波器(PSD)及低通滤波器(LPF)技术，从阻抗二维信息中分离出幅度和相位，给出了基于分立元器件的涡流检测信号X-R正交分解电路原理图，奠定了阻抗信号定量分析的基础，但该电路仅针对单频激励的特征信息提取，对需要多个频率才能反映试件特性的情况有较大局限性；高军哲提出基于锁相放大器的阻抗信息求解方法，采用模拟乘法器和模拟低通滤波器，设计了锁相放大电路，受限于实现方案中采用8阶巴特沃斯滤波器作为低通滤波器，其-3dB截止频率仅为30Hz，输出幅值约有33ms周期的抖动，在测试结果的重复性上需要提升。

因此，如何提供一种准确度高、实时性强的涡流阻抗求解方法及装置是本领域技术人员亟需解决的问题。

技术实现要素：

有鉴于此，本发明提供了一种基于互相关算法的涡流阻抗求解方法及装置，在阐述互相关算法原理的基础上，基于Xilinx FPGA硬件平台，给出了互相关算法的核心组件的数字化设计，因算法仅需单片FPGA及少量接口器件即可实现，利于检测设备的小型化和便携化；经互相关算法得到的结果与LCR测试结果比对，误差保持在±0.5％内,算法准确性高；互相关算法约2.4s完成8组频率响应下的单只工件分选，算法实时性强，适合在线检测；同相及正交支路分量输出趋近直流的程度是互相关算法精度的体现，同相支路的直流输出仍有低频抖动，后续工作的开展从成形滤波后引入滑动平均滤波算法入手，让更多周期的分量数据参与均值运算，则直流分量计算准确度预计会进一步提高。

为了实现上述目的，本发明提供如下技术方案：

一种基于互相关算法的涡流阻抗求解方法，具体步骤包括如下：

步骤一：得到涡流检测线圈的阻抗信号，将所述阻抗信号与高斯噪声信号叠加；

步骤二：引进与阻抗信号同频的同相参考信号，同时引进与所述同相参考信号相差90°的同频信号正交参考信号；

步骤三：计算叠加后的信号分别与同相参考信号、正交参考信号的乘积；

步骤四：进行滤波通过三角函数运算可得到线圈阻抗的幅值和相位值。

优选的，在上述的一种基于互相关算法的涡流阻抗求解方法中，所述步骤一中涡流检测线圈的阻抗信号为：

式中Vsig为阻抗信号的幅值，为阻抗信号的频率，φ0为线圈阻抗信号的相位；阻抗信号s(t)与高斯噪声信号n(t)叠加后的表达式x(t)为：

优选的，在上述的一种基于互相关算法的涡流阻抗求解方法中，所述步骤二中，同相参考信号r1(t)与阻抗信号同频，φr为同相参考信号的相位值，r1(t)为：

正交参考信号r2(t)是与r1(t)相差90°的同频信号，φr为正交参考信号的相位值，则r2(t)为：

优选的，在上述的一种基于互相关算法的涡流阻抗求解方法中，所述步骤三中，叠加后的信号x(t)与同相参考信号r1(t)的乘积记为Vpsd1(t)：

叠加后的信号x(t)与正交参考信号r2(t)的乘积记为Vpsd2(t)：

优选的，在上述的一种基于互相关算法的涡流阻抗求解方法中，所述步骤四中，将步骤三中的乘积通过截止频率仅为0.5Hz的低通滤波后有：

进一步，将参考信号的幅度值Vr设置成单位电压，系数1/2为系统增益，令相角差φsig＝φ0-φr，通过三角函数运算可得到线圈阻抗的幅值Vsig和相位值φsig，则有：

一种基于互相关算法的涡流阻抗求解装置，包括：DDS核、电桥差动回路、功率放大器、乘法器、第一滤波器、第二滤波器和CORDIC电路；其中，DDS核提供余弦激励信号，同时提供余弦信号作为同向参考信号，提供正弦信号作为正交参考信号；激励信号经D/A输出后与参考线圈、检测线圈组成电桥差动回路；所述电桥差动电路输出线圈阻抗信号；阻抗信号经放大后进入A/D转化成数字量，进行互相关算法处理；互相关运算方向分为同相支路和正交支路，两支路在FPGA内同步进行，阻抗信号与同相参考信号和正交参考信号经乘法器，同步通过第一滤波器及第二滤波器低通滤波处理，后得到阻抗同向和正交幅值，最后经CORDIC电路运算得到阻抗幅值Amp1和相角值Angle1。其中，注：菱形图左下代表检测线圈，右下代表参考线圈；左上代表与线圈检测线圈匹配的电阻，右上代表与参考线圈匹配的电阻。

优选的，在上述的一种基于互相关算法的涡流阻抗求解装置中，所述第一滤波器和所述第二滤波器采用CIC+HB+FIR结构进行低通滤波器的设计，先对阻抗倍频信号进行梳状滤波和半带滤波处理，降低倍频信号的数据速率，再用较低阶数的有限长单位冲击响应滤波器成形处理，滤除经梳状滤波和半带滤波后阻带内的频谱混叠，满足最后一级滤波的要求，进而准确地获取阻抗幅度信息。

优选的，在上述的一种基于互相关算法的涡流阻抗求解装置中，经CIC滤波后的同相支路信号进入同相支路第1级HB滤波，依次经过同相支路第2级、第3级、第4级HB滤波、第5级FIR滤波输出；经CIC滤波后的正交支路信号进入同相支路第1级HB滤波，依次经过正交支路第2级、第3级、第4级HB滤波、第5级FIR滤波输出。

优选的，在上述的一种基于互相关算法的涡流阻抗求解装置中，所述CORDIC电路采用Xilinx Artix-7内置的CORDIC IP核完成阻抗幅度值Vsig和相位值φsig，IP核求解模式选用定向模式。

经由上述的技术方案可知，与现有技术相比，本发明公开提供了一种基于互相关算法的涡流阻抗求解方法及装置，在阐述互相关算法原理的基础上，基于Xilinx FPGA硬件平台，给出了互相关算法的核心组件的数字化设计，因算法仅需单片FPGA及少量接口器件即可实现，利于检测设备的小型化和便携化；经互相关算法得到的结果与LCR测试结果比对，误差保持在±0.5％内,算法准确性高；互相关算法约2.4s完成8组频率响应下的单只工件分选，算法实时性强，适合在线检测。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据提供的附图获得其他的附图。

图1附图为本发明的结构示意图；

图2附图为本发明实施例DDS参考信号与激励源状态机跳转图；

图3附图为本发明实施例参考信号和预置激励源实测曲线；

图4附图为本发明实施例同相支路(a)和正交支路(b)乘法器输出信号实测曲线；

图5附图为本发明实施例CIC滤波器幅频响应图；

图6附图为本发明实施例同相支路(a)和正交支路(b)CIC滤波器输出信号曲线；

图7附图为本发明实施例HB滤波器响应曲线；

图8附图为本发明实施例FIR滤波器响应曲线；

图9附图为本发明实施例联合滤波频域(a)和时域(b)波形仿真曲线；

图10附图为本发明实施例联合滤波实测输出曲线；

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明实施例公开了一种基于互相关算法的涡流阻抗求解方法及装置，在阐述互相关算法原理的基础上，基于Xilinx FPGA硬件平台，给出了互相关算法的核心组件的数字化设计，因算法仅需单片FPGA及少量接口器件即可实现，利于检测设备的小型化和便携化。

一种基于互相关算法的涡流阻抗求解装置包括：

算法平台：Xilinx Artix-7FPGA开发板、AD9226模/数接口板、AD9767数/模接口板；辅助硬件：功率放大电路、增益调节电路、含线圈的电桥电路；测试仪器：泰克示波器2024C、是德示波器3024A、安捷伦函数发生器33500B、LCR测试仪TH2816B。为便于对算法组件进行FPGA设计，并对比不同频率下线圈阻抗的变化情况。

DDS参考信号与线圈激励源：

激励源采用Xilinx Artix-7DDS(直接数字频率合成器)IP核实现，通过控制DDS相位角输出27Hz，73Hz，194Hz，520Hz，1.4KHz，3.7KHz，10KHz，27KHz等8组预置正弦信号，相位角的变化由状态机控制，图2为状态机跳转图。以27Hz频率信号的产生为例，阐述状态机的跳转过程。rdy为状态机准备就绪信号，高电平有效，状态机s1～s8状态分别对应8个预置频率信号，系统时钟f＝50MHz，相位角范围Φ∈[-π，π]，相位角采用32位数据量化，含3位量化整数位及29位量化小数位，量化后相位角步进量N为累加次数，当NR-π＜π时，相位角以-π为起始点，以R为值进行累加，累加过程就是25Hz信号在[-π，π]内一个周期输出波形的过程；当NR-π≥π时，状态机从状态s1跳转至状态s2，对67Hz频率下的R进行累加，准备67Hz频率波形的输出。同理，其他状态的跳转条件与s1跳转至s2状态相似，根据状态机的循环跳转，可产生8组扫频频率信号。

状态机控制的8组预置频率信号经双通道D/A输出，用泰克示波器测得结果如图3所示，1通道表示正弦信号，2通道表示余弦信号，且相互正交。1通道信号仅作为检测系统的正交参考信号，2通道信号作为检测系统的同向参考信号的同时，还作为线圈的激励源直接输出到检测线圈。经测试，参考信号和预置激励源的频率误差满足设计指标的0.01％。

乘法器：

乘法器是互相关算法的核心，互相关算法需两路乘法器，两乘法器的被乘数分别为同相参考信号及正交参考信号，乘数端为检测线圈的阻抗信号。线圈阻抗信号经A/D采样后，经乘法器后的两支路乘法器的结果如图4，因参考信号和阻抗信号同频，所以经乘法器后相当于把阻抗信号倍频，仅幅度减小了1/2倍，乘积后同相和正交支路仍保持严格的正交关系，但正交支路幅值位于零点以上。

低通滤波器：

低通滤波器设计是互相关算法的最关键部分，其作用是将倍频后的阻抗信号整形为定幅的直流输出。为保证检测系统的信噪比，要求低通滤波器的通频带很窄，滤波器的Q(品质因素)值很高，直接设计较为困难。信噪比改善SINR表达式为(11)，B表示激励源信号的噪声带宽，BN表示检测系统的等效噪声带宽。

激励源经乘法器后最大频率fmax＝50KHz，若噪声带宽B＝500Hz，要求达到SINR＝100的情况下，需令BN＝0.05Hz，要求滤波器Q值高达10⁶；进一步地，令A/D采样频率为62.5KHz、通带截止频率fp＝0.5Hz、阻带截止频率fc＝2Hz、通带波纹dev＝0.01dB、阻带衰减为80dB时，利用Matlab工具设计出来的滤波器高达209101阶，这对于任何工程实践说，几乎是不可能完成的任务。基于以上原因，采用CIC+HB+FIR结构进行低通滤波器的设计，先对阻抗倍频信号进行CIC(梳状滤波)和HB(半带滤波)处理，降低倍频信号的数据速率，再用较低阶数的FIR(有限长单位冲击响应滤波器)成形处理，滤除经CIC和HB后阻带内的频谱混叠，满足最后一级滤波的要求，进而准确地获取阻抗幅度信息。

CIC滤波器：

CIC滤波器是一种具有线性相位的特殊FIR滤波器，基于零点相消，已经被证明是在高速抽取中非常有效的单元，由于CIC滤波器只有加法器、寄存器和积分器，抽取过程没有乘法运算只有加法运算，便于节省运算量，实现形式简单，缺点是阻带衰减不明显，需要采用多级级联的方式来满足阻带衰减。其冲击响应h(n)为:

其中，M为滤波器长度，按照降采样后频段内不产生频率混叠、占用硬件资源少且兼顾运算速度的原则，采用5级级联结构，128倍抽取，经CIC抽取后，阻抗倍频信号的数据率从62.5KHz降低到488Hz，设计后的CIC滤波器在Artix-7的CIC滤波器IP核中的频率响应如图5，同相和正交支路的设计参数相同。

阻抗倍频信号经CIC滤波后的结果如图6，经CIC滤波后，同相和正交支路对5组高频阻抗倍频信号幅度衰减十分明显，小于488Hz的3个低频阻抗倍频信号几乎没有衰减，8组预置频率锯齿波化，数据速率得到明显降低。按照设计目标得到的实测结果和FPGA中的实际幅频响应曲线相吻合，各级滤波器组件设计指标实现良好。

半带滤波器及成形滤波器：

将半带滤波器和成形滤波器看成联合滤波系统进行整体设计，其设计过程为：联合滤波系统的初始数据率f＝488Hz，通带上限边缘频率fp＝0.4Hz，阻带下限边缘频率fc＝2Hz，系统采用5级抽取，即抽取因D＝F0/Fk＝488/15.25＝2⁵，抽取后数据率降为15.25Hz，前4级为HB半带滤波器，最后1级为FIR成形滤波器，前4级半带滤波器的通阻带容限devi＝1/4000，FIR成形滤波器通带阻带容限dev＝0.001。联合滤波系统部分Matlab代码如下：

fpi＝zeros(1，D-1)；

for i＝1:D-1

fpi(i)＝fp*(2^i)/F0；

end

b1＝firhalfband('minorder'，fpi(1)，devi)；L1＝length(b1)

b2＝firhalfband('minorder'，fpi(2)，devi)；L2＝length(b2)

b3＝firhalfband('minorder'，fpi(3)，devi)；L3＝length(b3)

b4＝firhalfband('minorder'，fpi(4)，devi)；L4＝length(b4)

fc＝[fp fs]；mag＝[1 0]；devk＝[devi，dev]；

[n，wn，beta，ftype]＝kaiserord(fc，mag，devk，Fout*2)；

fpm＝[0fc(1)*2/Fout fc(2)*2/Fout 1]；

magpm＝[1 1 0 0]；

b5＝firpm(n，fpm，magpm)；L5＝length(b5)

b1、b2、b3、b4分别对应前4级HB滤波器的系数，系数相同，其系数在Xilinx Artix-7滤波器IP核中幅频响应如图7；b5代表第5级FIR滤波器系数，其幅频响应如图8。红色曲线代表malab设计的理想滤波器曲线，蓝色曲线代表经过量化后的实际响应曲线。

为考察联合滤波系统整体性能，对其进行仿真。仿真工具：matlab。仿真数据：0.4Hz、25Hz、50Hz、67Hz的单频正弦信号，高斯白噪声序列。仿真步骤：1)单频正弦信号叠加高斯白噪声序列；2)叠加信号依次经半带滤波器和成形滤波器；3)用matlab画出滤波前和滤波后的频域波形和时域波形。仿真结果：如图9，原始混合信号只剩下0.4Hz的单频正弦信号，其他频率及高斯噪声序列被有效滤除。

经CIC滤波后的同相支路信号进入同相支路第1级HB滤波，依次经过同相支路第2级、第3级、第4级HB滤波、第5级FIR滤波输出；分别为10.a,10.c,10.e,10.g,10.i；经CIC滤波后的正交支路信号进入同相支路第1级HB滤波，依次经过正交支路第2级、第3级、第4级HB滤波、第5级FIR滤波输出，分别为10.b,10.d,10.f,10.h,10.j；从图10.i、10.j知知：同相支路信号经最后一级FIR成形滤波后，幅度无限趋近于0；正交支路信号经最后一级FIR成形滤波后，输出信号趋于平滑光整，幅度为原始信号的1/4。

CORDIC电路：

坐标旋转数字计算机CORDIC是解决计算机数学计算的一种高速算法，

运算代替乘法运算，使得矢量的旋转和定向计算不需要反三角函数、乘法、开平方等运算，适合于FPGA器件。同相支路第5级FIR滤波器输出作为式(9)、(10)的Vpsd11，正交支路第5级FIR滤波器输出作为式(9)、(10)的Vpsd21，采用Xilinx Artix-7内置的CORDIC IP核完成阻抗幅度值Vsig和相位值φsig，IP核求解模式选用定向模式。但由于最后一级成形滤波器数据更新率仅为15.25Hz，获取一次数据需约60ms的刷新时间，而进行一次幅度和相位角的CORDIC运算需要40个周期，运用脉冲沿检测及流水线技术，将数据域从15.25Hz时钟切换到610Hz时钟，仅需1个15.25Hz的时钟周期开销。

实时性分析：阻抗信号经乘法器运算流水线延迟为1个时钟周期，经CIC滤波处理的延迟为22个时钟周期，经HB滤波处理的延迟为10个时钟周期，经FIR滤波处理的延迟为41个时钟周期，其乘法器、CIC滤波、4级HB滤波、FIR滤波的工作频率分别为62.5KHz、62.5KHz、488Hz、244Hz、122Hz、61Hz、30.5Hz，CORDIC运算产生1个时钟周期开销，求和所有组件用时时间，经计算，互相关检测算法用时约2.4s。

算法结果与LCR测试对比：

当参考线圈内放置工件时，DDS产生的8组预置频率信号的每组持续时间由DDS平均设置成1s，记录基于互相关算法下的检测线圈阻抗幅值和相位值，其阻抗幅值测试对比结果如表1，表1中b1,c1,d1,e1为发明实施例幅度测试结果，表1中b2,c2,d2,e2为采用LCR仪幅度测试结果。

参考线圈放置标定气门而检测线圈空载时，输出是一个稳恒直流，可作为算法基准值；检测线圈内放置硬度不同的工件时，输出是正比于线圈阻抗

信号幅值的直流，二者的差值正好反映硬度区间的变化值，这充分证明了互相关算法用于硬度分选的工程适用性。硬度场分布相同的气门因线圈激励频率不同而引起的幅值变化差异较大，互相关算法非常容易能识别出来；硬度场分布差异微小的工件所引起的阻抗变化为1mVrms或10mVrms量级时，互相关算法仍能优异地区分其差异性。

本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。