技术领域
[0001] 本
发明属于空间环境探测领域，特别涉及该领域中的一种基于
非相干散射雷达的80km-160km流星探测方法。
背景技术
[0002] 来自星际空间的复合颗粒状物质，在进入
地球大气层(一般为160km以下)时会与周边的大气分子发生剧烈摩擦，产生大量的热、光及电离现象，从而形成流星事件。流星包括前端高速运动的流星体(产生流星头回波)以及后端的长薄抛物面状的电离气体(产生流星余迹回波)。McKinley在“Meteor Science and Engineering”一书中介绍了一种利用VHF(Very High Frequency)频段流星雷达观测流星余迹的方法，其基本原理是通过接收流星余迹对入射
电磁波的菲涅尔散射回波完成
对流星的观测。这种方法的不足之处在于对测量
角十分敏感，对于单站式雷达，流星余迹与雷达波束必须近似
正交才能保证后向散射回波能够被测量到；而且由于流星雷达只能检测到初始直径小于或等于雷达
波长的流星余迹，而流星余迹的初始直径会随高度的上升而增加。有研究指出VHF流星雷达实际只能够检测出其理论探测数量3.5％的流星；另外，流星雷达极少可以观测到流星头回波，这些不足限制了对流星的观测和研究。
发明内容
[0003] 本发明所要解决的技术问题就是提供一种基于非相干散射雷达的流星探测方法。
[0004] 本发明采用如下技术方案：
[0005] 一种基于非相干散射雷达的流星探测方法，其改进之处在于，包括如下步骤：
[0006] (1)将原始基带正交两路
采样数据一分为二，一路作为流星头事件检测的输入数据，一路作为流星余迹事件检测的输入数据，两个事件检测作为独立步骤并行运行；
[0007] (2)流星事件检测：
[0008] (21)流星头事件检测：
[0009] 流星体相对测站的运动表现为距离随时间变化的二阶减速运动，流星运动公式如下：
[0010]
[0011] 式(1)中r0代表目标的起始距离，v为目标的径向速度，而a则为目标的径向
加速度，最右项为扰动项，在流星高速运动中的影响可忽略不计，流星的相干回波可以如下表示：
[0012]
[0013] 式(2)中j代表虚数，Tp为雷达脉冲重复周期，A(t)表示流星在接收机端的回波
电压值，fd为流星体的多普勒频移，若雷达波长为0.6m，其多普勒频移在数十至数百kHz之间，最右一项为流星体在不同雷达重复周期中的初始
相位；
[0014] 将式(2)表示的M个IPP的N个脉冲采样按顺序排列并补零，并假设流星体在穿过雷达脉冲时的散射功率为恒定值，那么式(2)可改写成以下的离散形式：
[0015]
[0016] 其中，ωd与ωa分别代表归一化后的多普勒相位偏移和
相位漂移，分别对应径向速度和径向减速度，P为相邻IPP之间的采样点数，对式(3)做FFT变换后可得到积累后的
功率谱，功率谱的峰值对应流星回波
能量在对应
频率下的汇聚，该过程表示如下：
[0017]
[0018] 如果式(4)中的最大值argmax(·)大于判决
门限，则认为检测出一次流星头事件，判决门限根据噪声统计设置为噪声
平均功率谱
密度 的6—10倍，其中 的计算公式为
[0019] (22)流星余迹检测：
[0020] 从后检测滤波得到的原始回波采样值计算各高度区间上的时延剖面矩阵，接着根据该矩阵计算回波
信号在各高度区域上的自相关函数，离散信号自相关函数的计算方法如下：
[0021]
[0022] 将计算得到的自相关函数进行积累，再用模糊函数进行修正，可以得到
等离子体在不同高度区域上的自相关函数，最后再对其进行快速
傅立叶变换即可得到等离子体的功率谱；时间
分辨率设为2分钟，高度分辨率为4km，从等离子体功率剖面可计算得到
电子浓度和离子成分，并且可以推算出中性大气
风速，若80km至160km区域出现电子密度高于电离层电子密度2倍的薄层，且持续一定的时间，则认为检测到流星余迹；
[0023] (3)流星参数提取：
[0024] (31)若在步骤(21)中检测到流星头事件：
[0025] 流星体的径向速度可从功率谱的最大值所对应的频率计算而来，计算方法为Vr＝-fd·λ/2，其中fd为流星体的频移；
[0026] 根据得到的多普勒频移并结合发射码型，可近似给出接收信号的时域形式：
[0027] H(t)＝Aψ(t)[cos(2fd·t)+isin(2fd·t)] (6)
[0028] 上式中ψ(t)指代发射码型，将H(t)作为匹配
滤波器对连续M个IPP的接收信号S(τ,t)进行
脉冲压缩，再将各脉压值相加完成积累，最后以速度模糊函数为分辨率，分别在fd左右各选取一定个数的频点作为式(6)新的多普勒频移值，使用更新后的
匹配滤波器依次完成脉冲压缩，比较所有脉压值的峰值，选取最大脉压峰值对应的频率值作为更新后的多普勒频移，并且从脉压时域结果中读取出流星回波信号的时延τ，进而可计算得到流星体的轨道高度，计算公式为Alt＝c·τ/2·sinθ，其中c代表光速，θ为雷达波束仰角；
[0029] 校准雷达接收天线的噪声
温度Tsys，已知系统的噪声功率为Pnoise＝κTsysB，且流星头回波的接收功率为：
[0030]
[0031] 在上式中，G为天线增益，PT为发射机
峰值功率，λ为天线波长，σm为相干散射截面积，再根据脉冲压缩后的实时
信噪比 K为脉压增益，可计算得出散射截面积σm，又已知流星头回波来自包围在流星体附近距离小于λ/4的N个电子的相干散射，用公式表示则为 其中re为电子半径，计算公式为
故在得到σm后可通过该公式计算得到流星体附近相干散射区域内的电子数；
[0032] 若流星体在连续不小于6个雷达脉冲周期内被探测到，则能够得到足够多的径向速度的观测值，使用最小二乘法对速度值在时间上进行拟合可得到速度的变化率，即流星体的减速度，根据动量方程(8)可得到流星体的BP参数，即流星体质面比，[0033]
[0034] 上式中，Γ为
空气阻力系数，在此处设为1，ρ(z)为随高度变化的空气密度，假设流星体为理想球形，则上式可改写为 其中rm为流星体的等效半径，ρm为流星体的密度，在得到BP后可进而计算得出流星体的等效尺寸以及
质量；
[0035] (32)若在步骤(22)中检测到流星余迹事件：
[0036] 探测电离层
时空间被划分成一个个宽度为距离分辨率的高度区间，通过对不同区间内回波信号计算零时延处的相关函数ACF(0)并进行解模糊运算可以求得对应高度上的平均回波功率，进而得到整个探测高度区域内的功率剖面，电子密度与回波功率的关系如下：
[0037]
[0038] 上式中C为常量，Pt为发射机功率，ΔPr(r2)为接收到距离为r处单位距离范围内的回波功率，Δr则为这个距离变化量，k为散射
波数，λD为德拜长度，计算公式为为电子—离子温度比，为等离子体参数，一般k2λD2＜＜1，假设 则式(9)可改写成：
[0039]
[0040] 式(10)的电子密度与等离子体参数无关，因此可以从回波功率值计算得出，在得到整个探测区间的电子浓度剖面后，通过与背景电离层密度进行对比，就可以直观找到流星余迹出现的时间以及轨道高度，流星事件起始和结束之间的时间差就是流星余迹在雷达波束中的整个
停留时间；
[0041] (33)若均未检测出流星头和流星余迹事件，则通过单纯的电离层回波数据更新判决门限；
[0042] (4)重复步骤(1)-(3)，对下一时刻采集的数据进行处理。
[0043] 本发明的有益效果是：
[0044] 本发明所公开基于非相干散射雷达的流星探测方法，克服了对雷达波束到达角的敏感性，拓宽了可探测流星区域的
覆盖范围，能实时检测流星事件并提取高度、径向速度、径向减速度以及等效散射截面积等流星参数。
[0045] 相比于传统的低功率VHF流星雷达，基于非相干散射雷达的流星探测方法既能够探测流星余迹，也能探测到流星头回波，扩展了对流星观测的范围。这个特性可以用于帮助人们研究流星体与电离层E层之间的耦合，以及统计每年进入到大气层中流星的质量通量(kg/year)。
[0046] 本发明所公开的方法对流星余迹的探测是基于处于热平衡态下的等离子体所激发的非相干散射，这种探测模式对雷达波束到达角不敏感，且不受流星初始直径的影响，因而能够探测高度分布更广的流星余迹，增加了对流星探测的效率。
[0047] 本发明所公开的方法可以提取多个流星参数，包括流星体的轨道高度、径向速度和径向减速度、等效散射截面积、质面比、尺寸、质量、以及流星余迹的电子浓度、漂移速度和停留时间等。另外非相干散射雷达相比普通的VHF流星雷达，具有更窄的探测波束，因而具有更高的横向分辨率，且
时间分辨率和距离分辨率均优于VHF雷达，故而求得的参数具有更好的时空
精度。
[0048] 在观测流星的同时，非相干散射雷达能够同时得到流星出现区域(电离层E层)的等离子体参数，包括电子浓度、电子温度、离子温度、等离子速度等，这有助于研究流星体与高空大气之间的耦合，对空间天气科学的研究具有重要意义。
[0049] 当流星体沿雷达波束方向运动时，本发明所公开的方法可以完整探测流星体在大气层烧蚀直至消失的整个过程，特别是测量质面比(BP)的变化，这对研究流星体如何在大气层
沉积物质具有重要意义。
附图说明
[0050] 图1是本发明
实施例1所公开方法的流程示意图；
[0051] 图2是流星回波时域采样数据的实例图；
[0052] 图3是本发明实施例1所公开方法的流星头回波检测示意图；
[0053] 图4是本发明实施例1所公开方法的流星余迹回波检测示意图。
具体实施方式
[0054] 为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图和实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
[0055] 非相干散射雷达是目前世界上用于空间环境监测和研究的最强大地基手段，具有高功率、大孔径以及低噪声等技术特点，能够直接测量电子温度、离子温度、电子密度、等离子体速度等众多电离层参数，也可用于流星体与空间碎片等空间物体的观测与研究。
[0056] 本实施例所公开的方法已在子午工程曲靖非相干散射雷达上进行了非公开应用，该雷达的主要参数为：
[0057]
[0058]
[0059] 实施例1，如图1所示，本实施例公开了一种基于非相干散射雷达的流星探测方法，包括如下步骤：
[0060] (1)将原始基带正交两路采样数据一分为二，一路作为流星头事件检测的输入数据，一路作为流星余迹事件检测的输入数据，两个事件检测作为独立步骤并行运行，流星回波时域采样数据的实例如图2所示，流星事件出现的轨道高度约为100km，信噪比约为20dB；
[0061] (2)流星事件检测：基本目的是从雷达的原始回波采样中判断流星事件的有无。一般通过在接收端检测流星头回波和/或流星余迹回波对流星进行观测。
[0062] (21)流星头事件检测：基本流程是首先将连续M个雷达脉冲周期(Inter-pulse Period–IPP，对应的时长一般十几个毫秒以内)中流星区域(80-160km)内的N个采样点(N取决于采样率)组成一个数据序列，中间按时间间隔补零。对该数据序列进行
快速傅立叶变换(Fast Fourier Transform—FFT)，得到对应的功率谱。找到功率谱的最大值，并与判决门限进行比较。若超出判决门限，则视为检测到一次流星头事件。判决门限通过恒虚警检测方法计算得到，每隔一段时间进行实时更新。在该步骤中计算功率谱的实例可参见图3，图中功率谱通过对连续雷达回波采样数据进行快速傅立叶变换得到，其中横坐标代表多普勒速度，纵坐标代表相对回波功率大小。插图中所展示的一种情况是信噪比很高时在对应频率处的流星头回波的功率谱。以下对该步骤进行具体的说明：
[0063] 流星体相对测站的运动表现为距离随时间变化的二阶减速运动，流星运动公式如下：
[0064]
[0065] 式(1)中r0代表目标的起始距离，v为目标的径向速度，而a则为目标的径向加速度，最右项为扰动项，在流星高速运动中的影响可忽略不计，流星的相干回波可以如下表示：
[0066]
[0067] 式(2)中j代表虚数，Tp为雷达脉冲重复周期，A(t)表示流星在接收机端的回波电压值，fd为流星体的多普勒频移，若雷达波长为0.6m，其多普勒频移在数十至数百kHz之间，最右一项为流星体在不同雷达重复周期中的初始相位；
[0068] 将式(2)表示的M个IPP的N个脉冲采样按顺序排列并补零，并假设流星体在穿过雷达脉冲时的散射功率为恒定值，那么式(2)可改写成以下的离散形式：
[0069]
[0070] 其中，ωd与ωa分别代表归一化后的多普勒相位偏移和相位漂移，分别对应径向速度和径向减速度，P为相邻IPP之间的采样点数，对式(3)做FFT变换后可得到积累后的功率谱，功率谱的峰值对应流星回波能量在对应频率下的汇聚，该过程表示如下：
[0071]
[0072] 如果式(4)中的最大值argmax(·)大于判决门限，则认为检测出一次流星头事件，判决门限根据噪声统计设置为噪声平均功率谱密度 的6—10倍，其中 的计算公式为
[0073] (22)流星余迹检测：当流星余迹形成并扩散数十秒后，流星余迹中的等离子体处于平衡态。此时等离子体中存在两种模式的电子
波动：分别是离子
声波和朗谬尔波。这两种模式的波具有各向同性，不受雷达波束到达角的影响。此时等离子体对入射电磁波产生非相干散射，可被雷达探测到。其基本流程是对输入数据分别进行后检测滤波，生成时延剖面矩阵，计算自相关函数，解模糊函数以及非相干积累等操作。从最终积累后的零时延处的自相关函数求得对应高度的回波功率，进而求得整个流星区域(80-160km)的功率剖面。从功率剖面并结合公式(10)可以计算得到电离层电子浓度等参数。最后将当前雷达脉冲重复周期中的电子浓度剖面与背景电离层浓度进行比较。若大于背景电离层浓度2倍或以上，则认为检测到一次流星余迹事件。背景电离层浓度从不含流星余迹的回波数据中计算得到，每隔一段时间进行实时更新。在该步骤中计算功率剖面的实例可参见图4，其中横坐标代表不同积累周期的序号，纵坐标代表轨道高度。由图可见，在第3到第10个积累周期内，在110km轨道高度上有明显高于电离层回波功率的流星余迹回波。以下对该步骤进行具体的说明：
[0074] 对探测流星余迹时用到的电离层非相干散射
信号处理方法如下所述：从后检测滤波得到的原始回波采样值计算各高度区间上的时延剖面矩阵，接着根据该矩阵计算回波信号在各高度区域上的自相关函数，离散信号自相关函数的计算方法如下：
[0075]
[0076] 将计算得到的自相关函数进行积累，再用模糊函数进行修正，消除调制方式和滤波器的影响，可以得到等离子体在不同高度区域上的自相关函数，最后再对其进行快速傅立叶变换即可得到等离子体的功率谱；时间分辨率取决于积累时长，一般设为2分钟，高度分辨率由模糊函数决定，一般为4km左右，从等离子体功率剖面可计算得到电子浓度和离子成分等物理参数，并且可以推算出中性大气风速，若80km至160km区域出现电子密度高于电离层电子密度2倍的薄层，且持续一定的时间，则认为检测到流星余迹；
[0077] (3)流星参数提取：若在步骤(2)中检测出流星(流星头/流星余迹)事件，则执行步骤(3)，从回波中提取出流星的相关参数：从流星头回波的数据中能够直接测量得到流星体出现的轨道高度以及径向速度。对径向速度随时间的变化值进行线性拟合可以反演得到流星体的径向减速度。根据流星体进入大气层时的动量方程可以计算得到流星体的弹道参数(ballistic parameter-BP)，即流星体的质量与表面积比的比率(质面比)。若假设流星体形状为球形且来源为
小行星，则可以进一步计算得到流星体的等效半径以及质量。
[0078] (31)若在步骤(21)中检测到流星头事件(流星体经过雷达波束)，则可从功率谱中得到流星体的多普勒频移。通过公式(6)以该频移一定范围内的频点值作为变量生成匹配
滤波器组，分别与回波数据进行脉冲压缩。以最大响应值对应的频率更新多普勒频移值，可从其计算得到流星体的径向速度。另外还可通过时延计算出流星体的轨道高度。根据回波信噪比可计算得到产生相干散射的电子数。如果在连续6个或以上的回波数据中计算得到流星体的径向速度，可利用最小二乘拟合计算得到流星体的减速度。通过动量方程(8)可以计算得到流星体的质面比，进而可求得等效直径和质量，最后保存上述结果。以下对该步骤进行具体的说明：
[0079] 流星体的径向速度可从步骤(2)得到的功率谱的最大值所对应的频率计算而来，计算方法为Vr＝-fd·λ/2，其中fd为流星体的频移；
[0080] 根据得到的多普勒频移并结合发射码型，可近似给出接收信号的时域形式：
[0081] H(t)＝Aψ(t)[cos(2fd·t)+isin(2fd·t)] (6)
[0082] 上式中ψ(t)指代发射码型，将H(t)作为匹配滤波器对连续M个IPP的接收信号S(τ,t)进行脉冲压缩，再将各脉压值相加完成积累，最后以速度模糊函数为分辨率，分别在fd左右各选取一定个数的频点作为式(6)新的多普勒频移值，使用更新后的匹配滤波器依次完成脉冲压缩，比较所有脉压值的峰值，选取最大脉压峰值对应的频率值作为更新后的多普勒频移，并且从脉压时域结果中读取出流星回波信号的时延τ，进而可计算得到流星体的轨道高度，计算公式为Alt＝c·τ/2·sinθ，其中c代表光速，θ为雷达波束仰角；
[0083] 对流星头回波等效散射截面积的测量需要校准雷达接收天线的噪声温度Tsys，已知系统的噪声功率为Pnoise＝κTsysB，且流星头回波的接收功率为：
[0084]
[0085] 在上式中，G为天线增益，PT为发射机峰值功率，λ为天线波长，σm为相干散射截面积，再根据脉冲压缩后的实时信噪比 K为脉压增益，可计算得出散射截面积σm，又已知流星头回波来自包围在流星体附近距离小于λ/4的N个电子的相干散射，用公式表示则为 其中re为电子半径，计算公式为
故在得到σm后可通过该公式计算得到流星体附近相干散射区域内的电子数；
[0086] 若流星体在连续不小于6个雷达脉冲周期内被探测到，则能够得到足够多的径向速度的观测值，使用最小二乘法对速度值在时间上进行拟合可得到速度的变化率，即流星体的减速度，在得到流星体的减速度后，根据动量方程(8)可得到流星体的BP参数，即流星体质面比，
[0087]
[0088] 上式中，Γ为空气阻力系数，在此处设为1，ρ(z)为随高度变化的空气密度，假设流星体为理想球形，则上式可改写为 其中rm为流星体的等效半径，ρm为流星体的密度，根据之前的研究结果，非相干散射雷达测量到的流星主要来自小行星体，其密度为1g/cm3。因此在得到BP后可进而计算得出流星体的等效尺寸以及质量；
[0089] (32)若在步骤(22)中检测到流星余迹事件，可以从电子密度剖面中得到流星余迹的轨道高度以及出现和消失的时间。从流星余迹轨道高度随时间的变化可以计算出流星余迹相对雷达波束的漂移速度，最后保存上述结果。以下对该步骤进行具体的说明：
[0090] 从流星余迹的回波功率能够直接测量得到流星余迹的电子浓度，轨道高度以及停留时间。通过非线性反演还能推算得到流星余迹的漂移速度(大气风速)。下面介绍电子密度的测量方法：探测电离层时空间被划分成一个个宽度为距离分辨率的高度区间，通过对不同区间内回波信号计算零时延处的相关函数ACF(0)并进行解模糊运算可以求得对应高度上的平均回波功率，进而得到整个探测高度区域内的功率剖面，电子密度与回波功率的关系如下：
[0091]
[0092] 上式中C为常量，Pt为发射机功率，ΔPr(r2)为接收到距离为r处单位距离范围内的回波功率，Δr则为这个距离变化量，k为散射波数，λD为德拜长度，计算公式为为电子—离子温度比，为等离子体参数，一般k2λD2＜＜1，假设 则式(9)可改写成：
[0093]
[0094] 式(10)的电子密度与等离子体参数无关，因此可以从回波功率值计算得出，在得到整个探测区间的电子浓度剖面后，通过与背景电离层密度进行对比，就可以直观找到流星余迹出现的时间以及轨道高度，流星事件起始和结束之间的时间差就是流星余迹在雷达波束中的整个停留时间；
[0095] (33)若均未检测出流星头和流星余迹事件，则通过单纯的电离层回波数据更新判决门限；
[0096] (4)重复步骤(1)-(3)，对下一时刻采集的数据进行处理。
