地震是最具有破坏力的自然灾害之一，20世纪以来，高达180万人死于地震, 经济损失高达数千亿美元^[1]。地震发震前总伴随有孕震过程，如何利用更加灵敏的仪器记录整个地震前兆的完整信息对于发展地震学科、研究地震深层次本质因素具有重要意义。目前常用的地震测量仪器如伸缩仪、应变仪以及加速度计等都是基于陶瓷式、电磁式以及摆式等结构，这些方案在实际应用中存在着种种问题：如陶瓷式地震计通常不具有非常高的灵敏度；电磁式地震计容易受到强电场等外界环境干扰，且经常出现很大的零漂现象，其抗雷击，抗工频串扰能力也略显不足；机械摆式的地震计动态范围小，工作带宽低，对于大震级的地震信号(180 dB)以及高频(1 kHz以上)剧烈震动信号很难响应^[2]。

光纤地震计是地震监测手段的一个新的发展方向，将光纤作为地震计的核心元件，其原理是利用光纤中传输的光的参数变化进行地震测量，先天具备灵敏度高、抗电磁干扰能力强等优势。另外光纤地震计耐高温、耐腐蚀、造价低廉，易大规模长距离铺设，甚至可以埋入温度较高的地层深处，既可以更加快速、准确地获得第一手地震信息，又同时可以对地层应力、磁场、地震波等物理量进行长期监测。

光纤地震计是利用干涉仪输出光信号参数的变化计算震动、位移、加速度等物理量。根据测量参数选取的不同分为幅度型、频率型以及相位型，其中相位型干涉仪具有最高的灵敏度与稳定性^[3-6]。常用的相位计算方法为相位生成载波(PGC)算法，该算法于1982年被美国海军实验室提出，因其具有实现简单、抗干扰能力强等特点，所以被广泛工程化应用^[7]。但是该算法属于开环检测方法，在布设环境参数变化很大、器件不稳定的情况下，容易引入较大的相位噪声与谐波失真，此时需要人为外加校准信号或者定期重置系统静态工作点^[8]，很大地影响了地震计的测量实时性。

文中提出了一种基于COSTAS环的光纤地震计解调算法，实际上是对传统PGC算法增加一个闭环^[9-11]相位反馈功能，实现对地震信号的实时测量。该算法依托FPGA与DSP数字信号处理平台，通过对光纤干涉信号进行正交锁相，消除检测系统因环境变化以及自身状态漂移引入的相位干扰，提高测量精度以及信噪比。

1 光纤地震计解调原理

光纤地震计中的干涉仪基本结构如图 1所示。

激光器将光注入到干涉仪中，系统采集光电探测器输出的干涉信号并计算干涉信号中的相位变化。

1.1 相位生成载波算法

PGC算法是相位型光纤干涉仪常用的解调方法。如图 1所示，系统产生调相波信号cosω₀t加载到压电陶瓷环上，此时干涉仪中产生固定相位变化，在光电探测器两端进行差分接收，得到干涉信号为

式中：P为接收到信号的光强值；A为干涉信号直流分量大小；B为干涉信号交流幅度；C为干涉仪中相位调制深度；cosω₀t为调相波；φ_s(t)为所求信号的相位变化，PGC算法即利用干涉信号中的多种频率分量求解相位变化量φ_s(t)。将式(1)中的干涉信号进行贝塞尔函数展开，可以得到干涉信号频谱成分：

式中：J_k(C)为贝塞尔函数系数；k为贝塞尔函数阶数。由式(2)可知，干涉信号中包含不同频次分量，此时需要利用COSTAS环结构对所求信号的不同频次分量进行提取。

1.2 COSTAS环原理

对干涉信号进行解调时，需要将被测信号的不同频率分量提取出来，其方法是在解调过程中提供一个与调相波信号同频同相的本地载波信号，获取这个本地载波信号的方法即载波同步过程。本系统采用COSTAS环进行载波同步，完成对不同频次分量的信号提取，其基本功能框图如图 2所示。

如图 2所示COSTAS环的主要功能是提取输入信号P₁中的频率为ω₀的信号分量成分，同时消除该频段的相移变化量Δφ。根据式(2)所示贝塞尔展开，取输入信号P的一阶分量P₁为

式中：cos (ω₀t+Δφ)为同步载波；Δφ为系统未知相移；J₁(C) sinφ_s(t)为所求信号。压控振荡器产生正交的2路信号即I、Q信号(正弦、余弦)，分别与输入信号P₁做相乘运算，得到混频后的信号为

混频所得信号分别经过低通滤波器后得到

式(3)、(4)均为锁相输出结果，但是由于地震计中的电路系统的相频特性不同，会带有误差项Δφ，将二者相乘后得到

此时，调整环路滤波器相位参数Δφ使输出信号P₆=0，即令sin2Δφ=0；因为在Δφ很小的情况下存在sin2Δφ≈Δφ，所以Δφ=0，则本地载波信号与调相波信号同频同相，取P₅信号为COSTAS环输出结果，可以得到

此时完成了对干涉信号P的一阶分量提取。根据算法需要，设计不同参数的COSTAS环分别完成对不同频率分量的信号提取。

2 光纤地震计系统设计

光纤地震计由硬件测量系统与相位解调算法组成，硬件部分为光纤光路与数字信号处理系统，相位解调算法依托FPGA与DSP数字信号处理系统，在该处理系统中完成对干涉相位信号的计算。

2.1 光纤地震计硬件结构

光纤地震计的硬件结构如图 3所示。

光纤地震计硬件主要包括光源、干涉仪、探测器、数字信号处理系统以及计算机5个部分。

系统工作过程中，光源将光注入至干涉仪中，数字信号处理系统产生的调相波信号cosω₀t，通过干涉仪中的压电陶瓷环加载到光纤上，使整个干涉仪有固定频率的相位变化。

此时在干涉仪输出端得到式(1)所示干涉信号。探测器实质为数字信号处理系统的模拟前端，负责将干涉的光信号转变为电信号并做适当的信号调理，包括光电转换的带宽限制、增益倍数的适量调节、转换后高频滤波等过程。

探测器输出的电信号由16位模数转换器采集输入数字信号处理系统中。该系统包括FPGA处理前端与DSP相位求解部分，同时嵌入COSTAS环与PGC相位解调算法，通过这2种算法将采集的干涉信号进行相位求解，最终结果输入至计算机中。

2.2 光纤地震计解调算法实现

经过模数转换后的干涉信号在FPGA与DSP数字信号处理系统中完成相位计算，其实现过程如图 4所示。

模数转换信号如式(1)所示，经过贝赛尔展开如式(2)所示，所得干涉信号包含k阶分量成分。通过2个COSTAS环分别对频率为cos ω₀t以及cos 2ω₀t的2种分量成分提取，得到信号为

式中P_I，P_Q为COSTAS环载波同步之后的正弦分量、余弦分量。将式(5)、(6)相除，得到被测信号相位的正切值：

此时使用反正切算法，即通过查表的方法即可求得干涉仪相位变化φ_s(t)。

3 算法仿真

根据式(3)、(4)所示，在不引入COSTAS环进行闭环锁相反馈时，系统解调准确度与相位误差Δφ有关。在不引入COSTAS环时，被测信号φ_s(t)的正弦分量、余弦分量式(5)、(6)则修正为

此时求得被测信号值为：

式中M为常数因子，设解调误差因子为α，则有

实际输出误差为

对该误差关系在零点处做泰勒展开，得到相位误差为

此时，对应解调结果中的谐波失真值THD为正确结果φ_s(t)与输出额外高阶分量的功率比值，可表示为

仿真参数设置为采样率f_s=2 MHz, 载波频率ω₀=20 kHz，取相位误差斜率较大处时，COSTAS环反馈进行校正仿真结果。

如图 5所示，通过COSTAS环对该相位进行补偿，此时对应解调结果的频域如图 6所示。未使用COSTAS环反馈前，解调结果谐波失真THD≈20.08%，经过COSTAS环反馈调节后，对应信号幅度明显提高，解调结果的THD≈0.52%。

4 光纤地震计实验

光纤地震计布设至黑龙江省哈尔滨市二龙山国家地震台(45.744 169 N、127.420 135 E)进行地震监测与数据记录，实际布设场景如图 7所示。

光纤地震计在地震台处进行全天候监测，地震台位于山体基岩相连处，可以良好地传导地脉动信号。

如图 8所示为光纤地震计记录马鲁古海附近地震，坐标为2.19 N、126.68 E，震级为6.3里氏。光纤地震计可以清晰地看到地震前兆信号P波信号波形，在1 min后有强烈地震波S波信号，X、Y、Z轴均能捕捉到高频峰值。

5 结束语

文中提出一种基于COSTAS环的光纤地震计解调算法，在传统PGC算法基础上增加相位反馈环路，用于校正电路以及压电陶瓷环产生的非线性相移。该方法有效抑制因相位变化产生的谐波失真，提高解调精度。该算法内嵌于FPGA与DSP数字信号处理系统，应用于光纤地震计中，实现对干涉信号相位计算，同时完成地震波前兆监测、地震波波形记录功能。