CN112797989B - 一种脉冲星频率参数的快速搜索方法 - Google Patents

Abstract

本发明提供一种脉冲星频率参数的快速搜索方法，具体为一种应用于星载探测器接收脉冲星信号的频率参数快速搜索方法，包括三大步骤：确定频率参数搜索的指标函数；应用动态缩减交叉熵方法进行脉冲信号频率参数大范围粗搜索，得到精确的搜索样本空间；利用小规模差分进化方法在精确样本空间中完成频率参数的快速迭代精化计算。该方案实现快速、准确的脉冲星实测信号频率参数估计，以完成实时在轨的动态信号处理过程，得到精确的脉冲信号频率及频率一阶导数的估计值。

Description

一种脉冲星频率参数的快速搜索方法

技术领域

本发明涉及航天航空技术领域，具体涉及一种脉冲星频率参数的快速搜索方法。

背景技术

X射线脉冲星导航目前已步入空间演示验证阶段，为满足脉冲星导航的实际工程应用需求，适用于在轨航天器的信号处理技术正成为研究的重点。

X射线脉冲星信号流量极弱，实际利用X射线脉冲星导航过程中，航天器需要在轨持续观测脉冲星以累积足够处理的光子，而航天器轨道运动导致的多普勒效应会使探测器接收到的脉冲星信号频率发生频移。对于动态信号处理技术，即利用预估轨道信息，其观测时间内航天器的位置信息通过轨道动力学模型外推获得，会引入模型的截断误差，且系统存在的轨道初始误差等也会被放大。若不对接收到包含轨道动态的光子TOA序列进行频率参数搜索，直接利用查询得到的频率参数信息进行轮廓折叠则会造成光子相位解算误差，导致光子发生相位偏移，后续周期的光子无法对齐，折叠轮廓发生畸变，进而影响脉冲TOA估计精度。因此，充分考虑航天器轨道动态的快速、高精度脉冲星频率参数搜索方法具有重要的实际应用价值。

脉冲星频率(周期)搜索方法可分为频域搜索方法及时域搜索方法两大类。由于光子TOA数据是非等间隔的离散序列，频域搜索方法多采用重采样或等间隔插值等近似处理方法，“频谱泄露”等问题会致使频率搜索精度较差，且经实测数据处理表明，大部分频域搜索方法的处理结果与实际情况差异明显；时域搜索方法一般具有较高精度，能够满足高精度的频率估计，但需求较大的搜索范围会影响计算效率，搜索过程指标函数的强噪声，多极值等问题也会影响频率的搜索精度。

现有技术中的搜索方法存在以下问题：(1)现有频率(周期)搜索方法大多未考虑航天器实际轨道动态或使用模拟轨道动态的仿真数据进行检验，但经实测数据处理表明，大部分频域搜索方法的处理结果与实际情况差异明显，估计精度差，计算量巨大；(2)现有频率(周期)搜索方法多仅对接收脉冲星信号的频率(周期)进行搜索，没有考虑航天器轨道动态对频率一阶导数造成的影响，导致误差吸收效果较差，处理的高精度与计算效率无法兼顾；(3)传统网格式二维搜索算法计算效率极差，且指标函数值分布杂乱，存在大量不同位置的局部取值高峰，是复杂的多极值优化问题，搜索过程中会比较容易陷入局部最优情况；(4)部分方法通过分批处理光子数据以实现降低计算量的频率(周期)搜索方法，由于每次使用的光子信息较少会使结果更容易收敛到局部最优，鲁棒性较差。

发明内容

本发明目的在于应用于星载探测器接收脉冲星信号的频率参数的快速搜索方法，实现快速、准确的脉冲星实测信号频率参数估计，以完成实时在轨的动态信号处理过程，得到精确的脉冲信号频率及频率一阶导数的估计值，具体技术方案如下：

一种脉冲星频率参数的快速搜索方法，包括以下步骤：

步骤一、确定频率参数搜索的指标函数如表达式7)：

其中：χ²为卡方检验值；k是轮廓折叠的Bin数，即封装段数；c_i是第i个封装段内的光子总数，

N_p为脉冲星周期总数，T_i为第i个封装段的中间时刻，c_j(T_i)为第j个周期折合入第i个封装段的光子数；

为封装段的平均光子数

N为总光子数，

继而获得脉冲信号频率及频率一阶导数的估计值如表达式8)：

其中：

和

分别为脉冲信号频率及频率一阶导数的估计值；

步骤二、应用动态缩减交叉熵方法进行脉冲信号频率参数大范围粗搜索，得到精确的搜索样本空间，具体是：通过添加Dynamic-decreasing Operator来降低交叉熵算法(即CE算法)迭代过程的计算量，根据交叉熵算法每次迭代计算及平滑处理后得到的样本方差采用表达式11)计算样本缩减因子β，对下次迭代的样本容量按β进行缩减；

其中：γ为防止参数空间样本密度过低而设置的保险阈值；β＝min(β_l)，

和

分别为第l次和第l-1次迭代计算所得样本方差；

步骤三、利用小规模差分进化方法在精确样本空间中完成频率参数的快速迭代精化计算，输出精确的频率参数估计值。

优选的，步骤二中精确的搜索样本空间通过以下步骤获得：

步骤2.1、系统初始参数设定，具体是：定义初始均值

和方差

均值

和方差

均包含代表F0及F1的两项；精英采样率ρ，初始样本容量N′，则精英样本容量N_e表示为N′·ρ；

步骤2.2、生成样本空间，具体是：根据正态分布的概率密度

生成N′个样本F₁，...，F_N′；F_i表示频率参数F0及F1的不同组合，即生成的样本空间为二维平面范围内的N′个不同点；

步骤2.3、指标函数计算及排序，具体是：根据生成的所有样本F₁，...，F_N计算各样本的指标函数S(F_i)，并对指标函数结果进行排序；

步骤2.4、精英样本筛选，具体是：根据指标函数值的排序结果选择最优的N_e个样本作为精英样本；

步骤2.5、更新分布参数，具体是：用采集的精英样本更新分布参数，采用表达式9)计算第l次迭代的均值

和方差

步骤2.6、平滑处理，具体是：采用表达式10)对

和

进行平滑处理：

其中：α为平滑因子，其取值范围为[0,1]；

步骤2.7、样本动态缩减，具体是：根据交叉熵算法每次迭代计算及平滑处理后得到的样本方差采用表达式11)计算样本缩减因子β，对下次迭代的样本容量按β进行缩：

步骤2.8、进行判断，具体是：若满足迭代停止条件，则结束并输出

和

否则，令l＝l+1，返回步骤2.2；迭代停止条件为：

或达到规定迭代次数，ε为标准差阈值。

优选的，所述步骤三具体包括以下步骤：

步骤3.1、初始化，具体是：利用步骤二得到的

及

作为种群初始化的基点及边界范围，其中：种群基点取值即

种群边界(a，b)的取值大于

中的两项，并设定mini-batch种群规模为N_m，随机初始化样本种群，记为表达式12)：

其中：个体

为初始化的脉冲星信号频率参数组合，用于表征问题的解；i＝1，2，...，N_m；

令s为迭代次数，即后续生成的样本种群表示为

相应个体即为

步骤3.2、种群变异操作，具体是：采用表达式13)进行种群变异操作：

其中：

为变异后的种群个体；r₁，r₂，r₃∈{1，2，...，N}互不相同且与i不同；K为变异系数，取值范围为[0,2]；

步骤3.3、种群交叉操作，具体是：采用表达式14)进行种群交叉操作：

其中：

表示种群交叉后的个体；rand(0，1)为[0,1]之间的均匀分布随机数；CR为范围在[0,1]之间的交叉概率；

步骤3.4、种群自然选择及排序，具体是：利用表达式15)对试验个体即

和

的目标函数进行比较：

其中：f()表示步骤一中选定的指标函数；

步骤3.5、进行判断，具体是：若满足迭代停止条件，则结束并输出

得到精确的脉冲信号频率及频率一阶导数的估计值

否则，令s＝s+1，返回步骤3.2；迭代停止条件为

或达到规定迭代次数，ε为标准差阀值。

应用本发明的技术方案，具有以下有益效果：

(1)本发明中的估算方法充分考虑了航天器轨道动态会对接收脉冲星信号的频率及其一阶导数产生的影响，对动态信号处理过程中产生的误差及被放大的初始误差等的吸收效果极好，在航天器存在较大初始轨道误差条件下具有良好的搜索效果，能有效提高脉冲TOA估计的精度，能吸收百米左右的速度误差，相对传统周期搜索方法，精度提升约5倍以上。

(2)本发明中的估算方法能充分利用观测时间内的接收到的全部脉冲星光子，不会出现由于光子信息较少使得结果收敛到局部最优的情况，频率参数搜索过程的鲁棒性极佳；

(3)本发明中的估算方法应用智能优化算法组合(采用dynamic-decreasing-CE算法迭代算法以及采用mini-batch-DE算法)，既扩大了初始搜索范围，同时也能完成搜索结果的快速迭代收敛。相对于传统智能优化算法、仅进行周期搜索等处理方法，其搜索效果、计算效率均有较大提升，在数据量充分的条件下，提升分别约为500％及80％左右。

除了上面所描述的目的、特征和优点之外，本发明还有其它的目的、特征和优点。下面将参照图，对本发明作进一步详细的说明。

附图说明

构成本申请的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：

图1是本发明实施例中算法组合流程图；

图2是本发明实施例中历元折叠原理示意图；

图3是本发明实施例中恢复观测轮廓与GPS模板轮廓对比图；

图4是轮廓对比相位差结果统计示意图；

图5是算法恢复轮廓平均对比相位差示意图；

图6是与一维周期搜索轮廓对比相位差结果示意图；

图7是算法搜索时间对比示意图；

图8是添加初始速度误差的组合算法结果示意图，其中：图8(a)和图8(b)分别为添加10m/s和100m/s的初始速度误差应用本发明方法所得轮廓示意图；图8(c)和图8(d)分别为添加10m/s和100m/s的初始速度误差应用传统智能优化方法所得轮廓示意图；图8(e)和图8(f)分别为添加10m/s和100m/s的初始速度误差应用传统周期搜索方法所得轮廓示意图；

图9是扩大搜索范围后算法计算时间对比示意图。

具体实施方式

以下结合附图对本发明的实施例进行详细说明。

实施例：

一种针对实测脉冲星信号频率参数的快速搜索方法，针对脉冲星信号处理过程中的频率参数搜索问题，在初始脉冲星数据库查询的频率参数不准和航天器轨道动态会对星载探测器接收到的脉冲星信号造成一定频率参数漂移的情况下，利用本实施例中所提方法(CE-DE组合算法)实现快速、准确的脉冲星实测信号频率参数估计，以完成实时在轨的动态信号处理过程，得到精确的脉冲星导航测量值，具体包括以下步骤：

第一步：确定频率参数搜索的指标函数；

第二步：应用dynamic-decreasing-CE(动态缩减交叉熵)方法进行脉冲信号频率参数大范围粗搜索，得到精确的搜索样本空间；

第三步：利用mini-batch-DE(小规模差分进化)方法在精确样本空间中完成频率参数的快速迭代精化计算。

航天器在轨接收到的脉冲星信号经过利用预估轨道进行质心修正、信号基点转换、星历表参数递推过程(脉冲星动态信号处理过程的前置部分内容参见现有技术)，得到如表达式1)所示的考虑航天器轨道动态的通用脉冲星相位传播模型：

其中：

为预估相位，

为待估相位；t表示任意光子到达时刻；φ(t)为t时刻航天器轨道动态的脉冲星相位；t₀是转换后的信号时间基点；f₀及f₁分别表示t₀时刻脉冲星的自传频率及自传频率一阶导数，可由脉冲星数据库查询获取或经过时间外推获得；m为航天器的模型线性化参数；定义

及γ_m为常矩阵，c为光速，n为脉冲星的方向矢量，n^T表示n的转秩；δr(t₀)及δv(t₀)分别表示t₀时刻航天器的位置偏差及速度偏差。此处：预估相位部分即通过航天器预估轨道信息及脉冲星可查询到数据库的信息能够计算得到的信号相位；待估相位部分则需要通过频率参数搜索方法处理得到。

对于近地航天器及深空探测器，通用相位传播模型的线性化参数取值不同。近地航天器的模型线性化参数一般取值m＝2；深空航天器的模型线性化参数一般取值m＝1。本实施例以近地航天器为例，其脉冲信号相位传播模型则可表示为表达式2)：

由表达式2)可以看出，考虑轨道动态的脉冲信号相位传播模型，其频率参数f₀及f₁均会发生变化，其形式可表示为表达式3)和表达式4)：

其中：F0及F1分别为去除原始光子到达时间序列的航天器轨道效应影响后的脉冲星频率和频率一阶导数，实际处理中可使用本实施例中的频率参数搜索方法快速获得。

本实施例的搜索方法详情如下，其中涉及的算法如图1所示：

一、确定频率参数搜索的指标函数，包括以下步骤：

1、历元折合

本实施例是基于历元折合的脉冲信号频率搜索方法，在计算过程中需要恢复出每组频率参数所代表脉冲信号的经验轮廓。历元折合作为常用的脉冲轮廓恢复方法，其基本原理详见图2所示。

历元折叠过程具体如下：

假设观测周期T_obs，含有N_p个脉冲星周期，即N_p＝floor(T_obs/P)。脉冲星周期P可划分为N_b个封装段，每个封装段的长度为T_b。

具体折叠步骤为：将记录在后续周期的光子TOA(Time of Arrival)序列折合到第一个周期；计算每个封装段中的光子数c_i(i＝1，2，...，N_b)；通过归一化光子数，得到恢复轮廓

采用表达式5)表示：

其中：

c_j(T_i)为第j个周期折合入第i个封装段的光子数，c_i是第i个封装段内的光子总数，T_i为第i个封装段的中间时刻。

2、确定指标函数

在本实施例中选择卡方检验值(χ²)作为频率搜索的指标函数，其原理为：给定试验频率参数F_i，F_i表示频率参数F0及F1的不同组合，采用表达式6)并利用mod(φ′(t)，1)进行历元折叠，得到每个封装段内的光子数c_i，φ'(t)表示由试验频率参数F_i计算得到的光子TOA经验相位，由此建立检验指标函数表达式7)：

其中：k是轮廓折叠的Bin数，即封装段数；总光子数

为封装段的平均光子数

则有表达式8)，即当使指标函数值达到最大时，得到所需的结果(脉冲信号频率及频率一阶导数的估计值)：

其中：

和

分别为脉冲信号频率及频率一阶导数的估计值。

二、应用dynamic-decreasing-CE方法进行脉冲信号频率参数大范围粗搜索，具体是：

交叉熵(Cross-Entropy,CE)是主要应用于解决小概率事件估计问题的全局优化算法，其通过概率密度函数Kullback-Leibler散度自适应地产生出一系列的概率密度函数，通过寻找理论上的最佳概率密度函数，使得参数收敛到最优解的邻域内。

本实施例具体包含以下步骤：

步骤2.1、系统初始参数设定，具体是：

定义初始均值

和方差

精英采样率ρ，初始样本容量N′，则精英样本容量N_e表示为N′·ρ；由于待搜索项为频率参数F0及F1的不同组合，因此均值

和方差

均包含代表F0及F1的两项；

步骤2.2、生成样本空间，具体是：

根据正态分布的概率密度

生成N′个样本F₁，...，F_N′；F_i表示频率参数F0及F1的不同组合，即生成的样本空间为二维平面范围内的N′个不同点；

步骤2.3、指标函数计算及排序，具体是：

根据生成的所有样本F₁，...，F_N计算各样本的指标函数S(F_i)(各样本的卡方检验值)，并对指标函数结果(即卡方检验值的大小)进行排序；

步骤2.4、精英样本筛选，具体是：

根据指标函数值的排序结果选择最优的N_e个样本作为精英样本；

步骤2.5、更新分布参数，具体是：

用采集的精英样本更新分布参数，采用表达式9)计算第l次迭代的均值

和方差

步骤2.6、平滑处理，具体是：采用表达式10)对

和

进行平滑处理：

其中：α为平滑因子，其取值范围为[0,1]；

步骤2.7、样本动态缩减，具体是：根据交叉熵算法(即CE算法)每次迭代计算及平滑处理后得到的样本方差采用表达式11)计算样本缩减因子β，对下次迭代的样本容量按β进行缩：

其中：γ为防止参数空间样本密度过低而设置的保险阈值，取值范围为[0,1]；β＝min(β_l)，

和

分别为第l次和第l-1次迭代计算所得样本方差；

步骤2.8、进行判断，具体是：若满足迭代停止条件，则结束并输出

和

否则，令l＝l+1，返回步骤2.2；迭代停止条件为：

或达到规定迭代次数，ε为标准差阈值；由于交叉熵算法过程计算量较大，并结合精化搜索空间要求，本案例设定为3次CE迭代。

三、应用mini-batch-DE方法完成精确样本空间中的频率参数快速搜索，具体是：

差分进化(Differential Evolution，DE)是一种基于群体进化的算法，具有记忆个体最优解和种群内信息共享的特点，即通过种群内个体间的合作与竞争来实现对优化问题的求解，其本质是一种基于实数编码的具有保优思想的贪婪遗传算法。基本思想是：对当前种群进行变异和交叉操作，产生另一个新种群；然后利用基于贪婪思想的选择操作对这两个种群进行一对一的选择，从而产生最终的新一代种群。

本实施例具体包括以下步骤：

步骤3.1、初始化，具体是：

利用步骤二得到的

及

作为种群初始化的基点及边界范围，其中：种群基点取值即

种群边界(a，b)的取值略大于

中的两项，并设定mini-batch种群规模为N_m，随机初始化样本种群，记为表达式12)：

其中：个体

为初始化的脉冲星信号频率参数组合，用于表征问题的解；i＝1，2，...，N_m；令s为迭代次数，即后续生成的样本种群表示为

相应个体即为

步骤3.2、种群变异操作，通过差分操作实现个体变异，具体是：

采用表达式13)进行种群变异操作：

其中：

为变异后的种群个体；r₁，r₂，r₃∈{1，2，...，N}互不相同且与i不同；K为变异系数，取值范围为[0,2]；

步骤3.3、种群交叉操作，目的为随机选择个体，具体是：采用表达式14)进行种群交叉操作：

其中：

表示种群交叉后的个体；rand(0，1)为[0,1]之间的均匀分布随机数；CR为范围在[0,1]之间的交叉概率；

步骤3.4、种群自然选择及排序，具体是：利用表达式15)对试验个体即

和

的目标函数进行比较：

其中：f()表示步骤一中选定的指标函数；

步骤3.5、进行判断，具体是：若满足迭代停止条件，则结束并输出

否则，令s＝s+1，返回步骤3.2；迭代停止条件为

或达到规定迭代次数，ε为标准差阀值，本案例中ε＝[10^-8，10^-12]；经多次测试，DE迭代结果在20次以内能够达到稳定状态，即当前停止条件也可设定为s＝20。

即为所求的精确的脉冲信号频率及频率一阶导数的估计值。

本实施例采用2017年8月31日至9月3日时间范围内慧眼HXMT卫星的Crab脉冲星实际观测数据，选取信噪比较高的HE高能探测器探测数据进行测试，详情如下：

预估轨道根数如表1所示：

表1预估轨道根数

根数	数值
		轨道历元(UTC)	2017/8/31 10:00:00
半长轴a	6949847.97m
		偏心率e	5.59e-4
轨道倾角i	42.944711°
		升交点赤经Ω	234.100076°
近地点幅角ω	71.209481°
		平近点角M	350.218648°

初始频率参数是Jodrell Bank发布的历元为2017年8月27日的Crab频率参数，如表2所示。

表2初始频率参数

根数	数值
		星历参考历元[MJD]	57992.165480818206561
F0[Hz]	29.639022542325657056
		F1[Hz s-1]	-3.6867173882290316336e-10
F2[Hz s-2]	8.8497789817619221739e-21

利用GPS轨道信息进行光子数据质心修正，然后直接利用初始计时参数递推到首光子到达时间处的递推计时参数进行历元折合，将得到的GPS模板轮廓视作模板进行对比。

应用本发明方案恢复的观测轮廓与GPS模板轮廓的对比结果如图3所示，折叠Bin数为1000。

选择观测时长在1800s以上的数据段对本发明方案及应用其它传统智能优化方法进行测试，增大相位分辨率，将轮廓折叠Bin数扩大到10000，分别对本发明方案和CE、DE、PSO(粒子群优化算法)方法性能进行对比分析，具体测试结果统计如图4所示。

平均统计结果如表3所示，应用交叉熵算法的动态信号处理精度最高，本发明方案综合了CE及DE算法，精度次之，与交叉熵算法的轮廓对比相位误差平均约为0.0003，且随着数据量的增加呈逼近趋势。

表3相位差统计结果

算法	平均相位差
		CE	0.00028
DE	0.00057
		PSO	0.00059
本发明方案	0.00047

图5给出了四种算法的对比相位差精度随观测数据长度变化的对比结果，可以看出，四种算法的相位差结果精度均会随着观测时间的增加而提高；本发明方案及应用CE的方法在数据量较少情况下具有较好的鲁棒性，对比相位差结果稳定，仍然能够得到高精度的脉冲TOA估计结果。

如图6所示，利用观测时间内100组不同的数据进行测试得到的与传统周期搜索方法处理的相位差结果对比。仅进行周期搜索恢复的轮廓与GPS对比轮廓的相位差多数会大于2e-3，且其回复轮廓的显著性会降低约1％左右，图中传统周期搜索方法恢复轮廓的相位差统计均值约为2.3e-3，本发明方案轮廓恢复的平均相位差约为4.7e-4，相对的精度提升约为500％。粗略计算，动态信号处理过程中传统周期搜索方法与本发明方案进行处理相比在导航方面会多造成10km左右的额外位置误差。

在计算时间方面，本发明方案中的交叉熵算法部分利用Dynamic-decreasingOperator能有效降低迭代过程的计算量。四种算法计算时间随数据长度的变化对比如图7所示。本发明方案在观测时间较短的情况下其计算时间不具优势，但随观测时间的增加变化较小，速度优势逐渐显现，在处理时间超过1000s的观测数据时计算速度优势明显。本发明方案在搜索精度上略逊于交叉熵算法，但计算速度具有极大优势。如图7所示，在数据长度约1800s时，本发明方案计算耗时约为140s左右，相比于交叉熵算法700s左右的耗时，计算效率提升约500％；相比于DE和PSO算法250s左右的耗时，计算效率提升约80％。

本实施例还在较为极端条件下进行算法测试，如下：

1、测试条件：通常情况下，由于航天器星载系统计算能力限制，为了在满足X射线脉冲星导航的实时性要求前提下得到足够精确的导航测量值解算，脉冲星周期搜索的范围一般较窄，这就造成了在动态信号处理过程中，若航天器预估轨道误差较大，探测到光子序列的实际频率会超出设定的搜索范围，搜索得到结果就是设置的边界，即频率搜索过程出现的截断误差，因此，会造成相当程度的相位偏差，无法准确估计脉冲TOA。

在动态信号处理过程中给预估轨道添加10m/s及100m/s的初始误差，分别对扩大搜索范围后的本发明方案，应用传统智能优化算法的二维频率参数搜索方法及传统周期搜索方法的轮廓恢复效果进行分析。

2、测试结果：

结果如图8所示，针对可能出现的航天器预估轨道初始误差较大，导航过程误差发散等情况，应用本发明方案可以有效吸收脉冲星频率参数发生的漂移量。

其它传统周期搜索方法及应用传统智能优化算法的二维搜索算法在扩大相应范围后计算代价将急剧增加，具体统计结果如图9所示。其它方法在能达到本发明方案的近似脉冲TOA估计精度的情况下，计算代价急剧增大。综上所述，本发明方案在精度及计算效率均有较为明显的优势。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (3)

1.一种脉冲星频率参数的快速搜索方法，其特征在于，包括以下步骤：

步骤一、确定频率参数搜索的指标函数如表达式7)：

其中：χ²为卡方检验值；k是轮廓折叠的Bin数，即封装段数；c_i是第i个封装段内的光子总数，

N_p为脉冲星周期总数，T_i为第i个封装段的中间时刻，c_j(T_i)为第j个周期折合入第i个封装段的光子数；

为封装段的平均光子数

N为总光子数，

继而获得脉冲信号频率及频率一阶导数的估计值如表达式8)：

其中：

和

分别为脉冲信号频率及频率一阶导数的估计值；

步骤二、应用动态缩减交叉熵方法进行脉冲信号频率参数大范围粗搜索，得到精确的搜索样本空间，具体是：通过添加Dynamic-decreasing Operator来降低交叉熵算法迭代过程的计算量，根据交叉熵算法每次迭代计算及平滑处理后得到的样本方差采用表达式11)计算样本缩减因子β，对下次迭代的样本容量按β进行缩减；

其中：γ为防止参数空间样本密度过低而设置的保险阈值；β＝min(β_l)，

和

分别为第l次和第l-1次迭代计算所得样本方差；

步骤三、利用小规模差分进化方法在精确样本空间中完成频率参数的快速迭代精化计算，输出精确的脉冲信号频率及频率一阶导数的估计值。

2.根据权利要求1所述的脉冲星频率参数的快速搜索方法，其特征在于，步骤二中精确的搜索样本空间通过以下步骤获得：

步骤2.1、系统初始参数设定，具体是：定义初始均值

和方差

均值

和方差

均包含代表F0及F1的两项；精英采样率ρ，初始样本容量N′，则精英样本容量N_e表示为N′·ρ；

步骤2.2、生成样本空间，具体是：根据正态分布的概率密度

生成N′个样本F₁，...，F_N′；F_i表示频率参数F0及F1的不同组合，即生成的样本空间为二维平面范围内的N′个不同点；

步骤2.3、指标函数计算及排序，具体是：根据生成的所有样本F₁，...，F_N计算各样本的指标函数S(F_i)，并对指标函数结果进行排序；

步骤2.4、精英样本筛选，具体是：根据指标函数值的排序结果选择最优的N_e个样本作为精英样本；

步骤2.5、更新分布参数，具体是：用采集的精英样本更新分布参数，采用表达式9)计算第l次迭代的均值

和方差

步骤2.6、平滑处理，具体是：采用表达式10)对

和

进行平滑处理：

其中：α为平滑因子，其取值范围为[0,1]；

步骤2.7、样本动态缩减，具体是：根据交叉熵算法每次迭代计算及平滑处理后得到的样本方差采用表达式11)计算样本缩减因子β，对下次迭代的样本容量按β进行缩：

步骤2.8、进行判断，具体是：若满足迭代停止条件，则结束并输出

和

否则，令l＝l+1，返回步骤2.2；迭代停止条件为：

或达到规定迭代次数，ε为标准差阈值。

3.根据权利要求2所述的脉冲星频率参数的快速搜索方法，其特征在于，所述步骤三具体包括以下步骤：

步骤3.1、初始化，具体是：利用步骤二得到的

及

作为种群初始化的基点及边界范围，其中：种群基点取值即

种群边界(a，b)的取值大于

中的两项，并设定mini-batch种群规模为N_m，随机初始化样本种群，记为表达式12)：

其中：个体

为初始化的脉冲星信号频率参数组合，用于表征问题的解；i＝1，2，...，N_m；

令s为迭代次数，即后续生成的样本种群表示为

相应个体即为

步骤3.2、种群变异操作，具体是：采用表达式13)进行种群变异操作：

其中：

为变异后的种群个体；r₁，r₂，r₃∈{1，2，...，N}互不相同且与i不同；K为变异系数，取值范围为[0,2]；

步骤3.3、种群交叉操作，具体是：采用表达式14)进行种群交叉操作：

其中：

表示种群交叉后的个体；rand(0，1)为[0,1]之间的均匀分布随机数；CR为范围在[0,1]之间的交叉概率；

步骤3.4、种群自然选择及排序，具体是：利用表达式15)对试验个体即

和

的目标函数进行比较：

其中：f()表示步骤一中选定的指标函数；

步骤3.5、进行判断，具体是：若满足迭代停止条件，则结束并输出

得到精确的脉冲信号频率及频率一阶导数的估计值

否则，令s＝s+1，返回步骤3.2；迭代停止条件为

或达到规定迭代次数，ε为标准差阈值。

Images