光纤声传感反正切解调的采样率与信号幅度关系*
2023-11-23安秉文吴先梅张金英
安秉文 吴先梅 张金英
(1 中国科学院声学研究所 声场声信息国家重点实验室 北京 100190)
(2 中国科学院大学 北京 100049)
(3 北京理工大学光电学院 北京 100081)
0 引言
声波在医疗成像、无损检测及地球物理勘探等领域中有着重要应用[1]。在声传感领域中,光纤光栅因其尺寸小、灵敏度高、抗电磁干扰、传输损耗小及易复用等优点[2-4],得到了广泛应用[5-7]。在光纤光栅声传感系统的几种不同解调方法中,基于3× 3 耦合器的干涉解调方法[8]从实现难易度、成本及解调精度上都表现优异,其常结合微分-交叉相乘(Differential cross-multiplying,DCM)算法[9-10]和反正切解调算法[11-13]进行解调。
当前声传感领域关于反正切解调算法的研究多在于改进算法[14-18]来改善3×3 耦合器的对称性和光电探测器的一致性问题,而关于幅度的研究较少。陈家熠等[19]研究了反正切解调所需的幅度下限,表明声信号引起的相位变化幅值需大于π。在实际应用反正切算法解调光纤光栅传感的声波信号时,若声场幅度较大,为得到不明显失真的解调结果,需要设置远超奈奎斯特采样定理对声波采样要求的采样率,且采样率随声场幅度增加也要相应增加,这表明声场幅度与采样率之间存在着某种约束关系。虽然也存在使用双波长激光源合成波长从而降低采样率的解决方法[20],但大幅提升了装置数量及解调复杂度,在实际中难以应用。张楠[21]研究了一种反正切解调算法在外差系统下的幅度动态范围,但由于外差法带来了额外的调制频率,不适用于医学超声等高频情况,且缺乏对不同相位展开方式的分析。
本文首先介绍了光纤光栅反正切解调原理及相位展开的数学模型,接着用正弦信号和高斯脉冲信号揭示了3 种不同相位展开方式的反正切算法所需采样率与信号幅度的线性关系及内在原因,并通过计算和实验进行了分析和验证。
1 光纤光栅反正切解调数学模型
以基于3×3 耦合器的光纤光栅声传感系统为例,介绍反正切算法的数学模型。
图1 为光纤光栅声传感系统结构图,光纤光栅传感器在声信号作用下中心波长发生偏移,借助基于3×3耦合器的马赫-曾德干涉仪(Mach–Zehnder interferometer,MZI)转化为相位变化,进一步经过三路光电探测器后被采集,输出的三路信号为
图1 光纤光栅声传感系统Fig.1 Fiber Bragg grating acoustical sensing system
式(1)中,Vi(i=1,2,3)为输出的第i路信号,B、C分别为三路输出电信号的直流分量和交流信号的幅度;s0为初始相位;s=Δλ为包含声传感信息的相位,其中λ为光栅的中心波长,d为马赫-曾德干涉仪的臂长差,Δλ为声信号引起的波长偏移,与声信号呈线性关系[22]:
其中,P11、P12为光纤的光弹系数,n为有效折射率,E为杨氏模量,ν为泊松比,P为声压;βi=2π(i-1)/3为3×3耦合器带来的相位差。
由于声信号幅度P、光栅中心波长偏移Δλ和相位s之间是线性关系,s的变化代表了声信号的变化,因此解调出s即可得到光纤光栅传感的声信号。
对式(1)进行三角函数和差化积运算,即可得到反正切解调出的信号上标表示是计算结果:
实际中采集到的信号是光电探测器的三路输出信号,为探究解调结果与原信号的差别,假定s已知。先给定信号s作为输入信号,再通过公式(1)计算出三路输出信号Vi,最后使用公式(3) 计算直接进行反正切解调的结果ˆs。图2以正弦信号s=20 sin(2πt)为例展示了上述过程,并以此说明在反正切算法中相位展开的必要性。
图2 相位展开前正弦信号的反正切解调Fig.2 Arctangent demodulation before phase unwrapping
图2(c)~(e)展示的三路输出信号与图2(a)的输入信号之间除波形的周期性、对称性外较难直观找到联系,这也是使用反正切或DCM 算法进一步解调的原因。图2(b)为直接使用公式(2)的反正切解调结果,其幅度被限定在反正切函数的值域(-π/2,π/2),因此出现相位缠绕现象,输入s=α+nπ (n为整数)时会解调出ˆs=α,α ∈(-π/2,π/2)。需进行相位展开还原真实的输入信号。
下面介绍3 种典型的相位展开方式。方式1为反正切算法最初被提出时使用的相位展开方式[11],即根据反正切值是否位于[-π/4,π/4],及前一反正切值的正负判断是否需在反正切值上加减π/2。方式2[12]不再根据反正切值本身是否属于不同区间,而是根据前后反正切值的差值是否属于[-π/2,π/2]而加减π。方式3[13]则是借助复平面的概念拓展至4 个象限,再根据两次返回值的差值绝对值是否属于[-π,π]而加减2π。本文将以上3 种相位解调方式称为3种反正切解调算法。
在信号幅度较小的情况下,这3 种解调算法均能够得到正确的解调结果。但是当信号幅度较大时,3种算法都可能出现解调结果失真。随着采样率增加,有的解调算法能够得到正确的解调结果,而有的算法仍然需要继续增加采样率才能得到正确的解调结果,这表明信号幅度和采样率之间的关系与相位展开方式有关,因此下文将进行不同相位展开方式下的采样率与信号幅度关系分析。
2 采样率与信号幅度关系分析
在声传感信号中,窄带的正弦信号和宽带的脉冲信号具有一定代表性,本节将对这两种信号的反正切解调算法进行分析,对比3 种不同相位展开方式下的解调结果,并与DCM 算法[10]进行比较,最后用脉冲信号实验进行验证。
2.1 正弦信号
给定正弦信号s=Asin(ωt),其中A为输入信号的幅度,ω=2πf为输入信号的角频率,f为信号频率。
由于公式(1)中的直流分量B可通过减去均值消除,交流信号的幅度C可进行幅度归一化,因此都取1,不再单独研究二者的影响;引入采样率fs与信号频率f的比值Kj=fs/f表示一周期内采样的点数,并记为相对采样率,下标j=1,2,3 表示第j种反正切解调方式。在接下来的分析中不失一般性地取f=1。输入信号仍用图2(a)所示的信号。
为误差判据,用于分析解调结果的失真度。计算了在采样率刚刚满足误差判据时3 种反正切解调方式及DCM 算法的解调结果与输入信号间的误差,如图3 所示,其中图3(a)、图3(c)、图3(e)和图3(g)依次为方式1~方式3 和DCM 算法解调出的信号,图3(b)、图3(d)、图3(f)和图3(h)依次为解调结果与原信号的误差|s-ˆs|。图3 中3 种反正切解调方式的相对采样率为K1=147、K2=80、K3=40,DCM 算法使用的相对采样率KDCM=160 与K1接近。计算时发现3 种反正切解调方式的误差均在3.55×10-15,而即使KDCM=2000时,其最大误差为0.003>0.01%×A=0.002。这是因为DCM 算法的求导、积分等步骤建立在微元假设上,因此误差较大,而反正切算法则不是。在计算时间上3种反正切解调方式均在10-4s量级,DCM算法为0.11 s,表明DCM在实时性上的不足。
图3 正弦信号的解调结果及误差Fig.3 Demodulation results and error of sine signal
为进一步探究采样率与信号幅度之间的关系,计算了信号幅度A ∈[1,1000]时3 种反正切解调方式在满足判据时所需的最小相对采样率,DCM 算法由于无法满足误差判据,不进行计算。得到相对采样率与输入信号幅度的关系如图4(a)所示,此时的最大误差如图4(b)所示,图中点线、划线和实线分别对应方式1、方式2和方式3。
图4 正弦信号不同反正切解调算法所需相对采样率及误差Fig.4 Relative sample rates and error of sine signal under different arctangent algorithms
由图4(a)看出3种方式的相对采样率与输入信号幅度存在着明显的线性关系,线性拟合后得到方式1~方式3 所需的相对采样率分别约为K1=8A、K2=4A、K3=2A。图4(c)表明在满足误差判据时,3 种方式解调出的信号与原信号之间的误差远小于误差判据,这体现了反正切解调算法具有较高的精度。
采样率与信号幅度相关的原因在于相位展开方式的判据是固定的数值,而输入信号是变化的。具体来说,3 种相位展开方式是根据相邻两次返回值的差值或是所属区间的关系来移动坐标基准、展开相位的,因此输入信号斜率最大处前后的导数(方式1、方式2 和方式3)或函数值(方式1)决定了相对采样率的要求。
以方式1 为例说明。正弦函数的导数是余弦函数,其在t=0 时导数最大,此时st=0=0。下一时刻:
当K1=8A时,若A较大,则:
正好处于方式1 的[-π/4,π/4]边缘,对应着相对采样率需要为输入信号幅度的8 倍。方式2 和方式3的判据分别为±π/2和±π,因此所需的相对采样率分别为信号幅度的4倍和2倍。
需要注意的是,如图4(a)局部放大后的图4(b),在输入信号幅度较小时,应设置一个采样率的下限以满足奈奎斯特定理。奈奎斯特定理要求采样率至少为信号频率的2 倍,即KN=2,实际中要采用更大的采样率,因此这里设置了相对采样率的下限KL=24。也出于这一考虑,在进行单个信号的分析时选择的信号幅度为A=20,此时3 种反正切方式需要的最低采样率高于设置的采样率下限,即K3=2A >KL。
图5 计算了图2(a)作为输入信号时两种不同采样率下各种算法的解调结果。图5(a)、图5(c)、图5(e)和图5(g)为K=2A时的解调结果,图5(b)、图5(d)、图5(f)和图5(h)为采样率增加一倍后即K=4A时的解调结果。从图5 可见,图5(a)~(c)的解调结果严重失真,波形各异,比较图5(b)与图5(a)可知,采样率更高的解调结果并未更接近输入信号,这表明存在一个采样率的阈值,使得低于该值时解调结果不可信。但比较图5(g)、图5(h)和图3(g),可发现DCM 算法的解调结果则随着相对采样率的增加逐渐接近输入信号,这是因为其越来越满足微元假设,DCM 在采样率较低时的解调结果可信度更好。由图5 中两种采样率的解调结果可知,3 种反正切解调方式中,方式1 对采样率的要求最高,方式2次之,方式3最低。
图5 不同采样率时的正弦信号解调结果Fig.5 Demodulation results of sine signal under different sample rates
2.2 高斯脉冲信号
为模拟在超声检测中常用的脉冲信号,进行了高斯脉冲信号的计算分析,其中心频率与实验用换能器中心频率一致,为f=2.5 MHz。为与前文统一引入相对采样率Kj=fs/f,其中f为高斯脉冲的中心频率。计算分析了不同带宽下的高斯脉冲信号解调结果,均与正弦信号基本一致,即相对采样率和输入信号幅度之间存在明显线性关系,且3 种反正切算法在计算时间和解调精度上都较DCM 算法表现更佳,这里出于篇幅不再赘述。
与正弦窄带信号不同的是,带宽对解调结果有一定影响。图6 展示了不同带宽下的输入高斯脉冲信号波形,图6(a)~(e)分别表示带宽为中心频率20%、40%、60%、80%和100%的情况。
图6 不同带宽的高斯脉冲信号Fig.6 Gaussian pulse signal of different bandwidths
图7 给出了不同带宽高斯脉冲信号所需的相对采样率与信号幅度的倍数与带宽之间的关系,图7 中菱形实线、圆圈划线和三角点线分别表示反正切解调方式1、方式2 和方式3,标注了每种方式在不同带宽下的相对采样率与幅度之间的拟合系数。
图7 高斯脉冲信号不同带宽下的相对采样率与幅度关系Fig.7 Relationship between relative sample rate and signal amplitude of Gaussian pulse signal under different bandwidths
由图7 中曲线变化趋势可以看到,随着带宽的增加,3 种解调方式所需的相对采样率均降低,且在高斯脉冲信号的带宽较窄时所需的相对采样率与正弦信号要求一致。
计算了不同带宽下高斯脉冲信号的导数,将导数最大值s′(t0)和函数值s(t0) 绘制在图8,t0为对应的时刻,时间步长为Δt=8.33×10-9s。三角划线和左侧纵坐标轴用来描述导数最大值,圆圈实线和右侧纵坐标轴描述最大导数对应时刻的函数值。
图8 高斯脉冲信号不同带宽下的导数最大值及其函数值Fig.8 Maximum derivatives of Gaussian pulse signal and their function values under different bandwidths
由图8 可解释高斯脉冲的相对采样率要求。不同带宽的导数最大值s′(t0)均在约3× 108,s(t0+Δt)-s(t0)=s′(t0)Δt ≈2.4,无论函数值s(t0)为多少,前后两次反正切值不同时属于[-π/4,π/4],方式1 需对当前及之后数据进行修正,加减。而=6与2π相近,约是3 种方式的判据π/4、π/2 和π 的8 倍、4 倍和2倍。此外,图8 中导数最大值随带宽的增加而减小的规律解释了所需相对采样率随带宽增加而减小的现象。
3 实验验证
使用脉冲信号发生器和2.5 MHz的压电换能器在水中发射声波,光纤光栅在距离声源60 mm 处接收声波。分别用625 MS/s,1250 MS/S和2500 MS/s的采样率采样,则3 种采样率将依次仅满足方式3、方式2~方式3,及满足方式1~方式2和方式3 的采样率要求。以625 MS/S的采样结果为例,将三路采集到的输出信号展示在图9。
图9 实验采集的三路输出信号Fig.9 Three outputs in experiment
下面依次分析不同采样率时不同解调方式下的解调结果。图10 展示了采样率为625 MS/s时3 种反正切解调方式及DCM 算法的解调结果,图10(a)~(d)分别对应反正切解调方式1、方式2、方式3 和DCM 算法,图中标注了解调结果的峰峰值Avv(峰峰值是信号幅度的2倍)。
图10 625 MS/s 采样率时的解调结果Fig.10 Demodulation results of four methods under 625 MS/s
图10(a)~(b)均出现了信号前后直流分量解调结果不一致的“相位跳变”现象,且与图10(c)的幅度相比可以看出信号未解调完全。DCM 算法虽未出现这一现象,但其幅度低于3种反正切解调方式。
从计算时间上看,3 种反正切解调方式的计算时间均在数十毫秒,而DCM算法则超过20 min,无法满足高频大数据量输入信号的实时解调,这与前文仿真的结论一致。因此在后文1250 MS/s 和2500 MS/s的解调中不再与DCM算法比较。
图11 计算了采样率1250 MS/s 和2500 MS/s时3种反正切解调方式的解调结果。
图11 1250 MS/s 和2500 MS/s 采样率时的解调结果Fig.11 Demodulation results of three arctangent algorithms under 1250 MS/s and 2500 MS/s
图11(a)出现了“相位跳变”现象,图11(b)~(f)解调结果的峰峰值均为146.2,可以相信信号的相位幅度约为73,而相对采样率与幅度的比值为≈3.4,位于2 倍和4 倍之间,故3 次实验按采样率由小到大依次满足了方式3、方式2~方式3和方式1 方式3 的采样率要求。方式3 对采样率的要求最低,因此在实际中应用较多。不同采样率下解调结果的结论与前文计算分析一致,也验证了前文计算的正确性。
4 结论
本文基于光纤光栅3×3 干涉解调的反正切算法及相位展开的数学模型,数值计算分析了3 种反正切算法下采样率与信号幅度之间的关系,并与DCM算法进行了比较,在实验中进行验证。得到结论如下:
(1) 揭示了采样率和信号幅度相关的原因在于相位展开的判据,且频率固定时3 种反正切解调方法均与信号幅度之间呈线性。以是否跨越[-π/4,π/4]及函数值符号为判据的反正切解调方式1 所需的采样率为fs=8Af以相邻两次返回函数值差值是否属于[-π/2,π/2]和[-π,π]为判据的方式2 和方式3 所需的采样率分别为fs=4Af和fs=2Af。采样率与信号幅度之间线性函数的系数取决于信号斜率最大处的导数(方式1~方式3)及函数值(仅方式1)。
(2) 对于高斯脉冲信号,反正切解调的采样率要求随着带宽增加而略有降低,这是因为导数的最大值随带宽增加而下降。从倍数关系的数值来看,实验时基本上仍需向上取整来按照结论(1)的采样率要求保证信号不失真。
(3) 与3 种反正切解调方式相比,DCM 算法无论从实时性(采样率要求、计算时间)还是计算精度的角度都表现都更逊色一些,但在采样率不足以满足反正切解调时不失为一个好的选择。当用硬件实现积分、微分等步骤时或许可以实现实时解调。
(4) 从仿真结果及实验结果来看,借助了复平面的反正切解调方式3 对采样率的要求最低,在应用前可预先根据先验知识或水听器测试估算传感器处声场的大致幅度,再根据所用光纤光栅传感器的灵敏度综合计算出所需采样率。
值得注意的是,借助复平面相位展开的方式3需要输入实数,否则可能错误判断象限,造成解调结果的偏差甚至错误。当出于干涉仪相位的不稳定、输出信号幅度差异等因素使用椭圆拟合时,由于拟合后的三路信号为复数信号,此时应选择其他方式解调,如在实验时采样率加倍或后处理时插值后再使用方式2解调。