CN103412298A - 一种自动获取船舶螺旋桨变速旋转时间区间方法 - Google Patents

Abstract

本发明涉及一种自动获取船舶螺旋桨变速旋转时间区间方法，属于水声信号处理技术领域。本发明基于船舶水下辐射的螺旋桨空化噪声受螺旋桨转速调制的特点，及螺旋桨匀速旋转和变速旋转时调制包络信号的特点与区别，运用图像处理的方法对经恒虚警处理的调制包络谱瀑布图上代表其叶频分量的谱线的识别与判断，将螺旋桨匀速旋转的时间段与变速旋转的时间段予以区分，从而获得螺旋桨变速旋转的起止时刻。本发明解决了无人值守的被动水声探测系统无法自动获取螺旋桨加速或减速旋转的时间段，进而无法获取螺旋桨完全处于变速旋转过程中辐射噪声数据的问题，未增加额外的硬件设备，从而为被动水声探测系统进一步正确提取螺旋桨旋转速度变化特征提供了保障。

Description

一种自动获取船舶螺旋桨变速旋转时间区间方法

技术领域

本发明涉及水声信号处理技术领域，具体涉及一种自动获取船舶螺旋桨变速旋转时间区间方法。

背景技术

船舶在进出港口时，通常会经历航速变化等复杂的运动过程。加速离港或减速进港是它们的典型运动状态。螺旋桨旋转速度变化的快慢，与船舶的类别有着一定的联系，若能够获取螺旋桨处于加速或减速旋转过程中所辐射的空化噪声数据，则可以通过基于分数阶傅里叶变换等的信号处理的方法从中提取螺旋桨旋转速度的变化特征，从而可用于进一步提升港口水域被动水声探测系统区分不同类型船舶的能力。

然而，对于无人值守的被动水声探测系统而言，采用现有的水声信号处理流程，不能够自动识别螺旋桨处于变速旋转的过程。螺旋桨完全处于加速或减速旋转状态下的船舶辐射噪声数据。由于螺旋桨空化噪声幅度在匀速旋转状态下受周期信号调制，而在加速或减速旋转状态下受chirp周期信号调制，因此在提取螺旋桨转速特征时，应该将螺旋桨变速旋转时接收的信号与匀速旋转时接收的信号分开处理，并且应采用不同的特征提取方法。若不进行区分而直接运用基于分数阶傅里叶变换等的方法处理包络解调后输出的信号，由于无法在整个观测时间内实现调频率的匹配，因此得到的特征提取结果会产生非常大的估计误差，不能获取准确螺旋桨旋转速度变化特征值。因此，如何从观测时间内接收到的船舶辐射噪声信号中，运用恰当的信号处理流程，正确识别出螺旋桨变速旋转的起止时间，是正确提取螺旋桨旋转速度变化特征的先决条件。

发明内容

(一)要解决的技术问题

本发明要解决的技术问题是：如何针对无人值守的被动水声探测系统无法自动获取船舶螺旋桨完全处于加速或减速旋转过程时的船舶水下辐射噪声信号问题，提出一种自动获取船舶螺旋桨变速旋转时间区间的方法。

(二)技术方案

为了解决上述技术问题，本发明提供了一种自动获取船舶螺旋桨变速旋转时间区间方法，包括以下步骤：

S1、对采样后的观测信号r(n)依次进行带通滤波、绝对值包络检波，并消除直流分量，得到包络信号，其中，所述绝对值包络检波的过程为：

$E_0\left(l\right)=\frac{1}{M}\underset{j=1}{\overset{M}{\mathrm{\∑}}}\left|r(\mathrm{Mn}+j)\right|---\left(1\right)$

其中M是绝对值检波的平均点数，满足：

M≤int(F_s/f_h)                                      (2)

式(2)中F_S为采样频率，f_h是带通滤波器的最高频率，int表示取整函数，在绝对值检波后，采样率降为F_S/M；

消除直流分量的过程为：

$E\left(l\right)=E_0\left(l\right)-\frac{1}{L}\underset{l=1}{\overset{L}{\mathrm{\∑}}}{E}_0\left(l\right)---\left(3\right)$

式中L是E₀(l)的点数；

S2、以宽度为T_W的滑动窗对E(l)进行截取，然后依次进行2^K点快速傅里叶变换和频域恒虚警率处理，K为正整数，经频域恒虚警率处理后超过阈值的点保留，低于阈值的点置零，再对输出结果进行灰度转换，形成一行调制包络谱分析结果，频域恒虚警率处理的阈值Th由公式(4)确定：

Th＝m+dσ                                     (4)

式(4)中d是由恒虚警概率决定的检测因子，m为无船舶目标时经步骤S1处理后频谱幅度值的均值，σ为无船舶目标时经步骤S1处理后频谱幅度值的方差；

S3、以步进时间t_D将滑动窗进行滑动，重复步骤S2，获得多行调制包络谱分析结果，并形成瀑布图的形式，从而获得包络谱瀑布图，记为G(x，y)，其中x＝X₁，X₁+1，…，X_H和y＝1，2…Y分别代表时间和频率轴上的变量：

$\left\{\begin{array}{c}{X}_1=\mathrm{int}(2^{K-1}M{F}_l/F_s)\\ X_H=\mathrm{int}(2^{K-1}M{F}_h/F_s)\\ H=X_H-X_1+1\\ Y=\mathrm{int}((T-T_W)/t_D)\end{array}---\left(5\right)\right.$

式(5)中(F_l，F_h)为步骤S1的处理过程所分析的频率范围，T为观测信号r(n)的总观测时间；

S4、初始化标志位：Y1＝0，Y2＝Y，W1＝0，W2＝Y；

S5、对G(x，y)进行0角度的Radon变换，并将Radon变换的横坐标s转换成频率号x，转换公式为：

$x=\mathrm{int}\left(\frac{2^{K-1}M(F_h+F_l)}{F_S}\right)+s---\left(6\right)$

从而得到随频率变化的0角度Radon变换输出R(x)；

S6、对R(x)的峰值进行分析：

若R(x)＝0，变速时间区间识别过程结束；

当R(x)≠0时，确定R(x)最大值的位置，即峰值位置的频率号x_M，并计算该峰值位置两侧相临的连续超过该峰值一半的频率点数N_P，以获得谱峰宽度：

$W_R=\frac{N_P{F}_S}{2^{K}M}---\left(7\right)$

若谱峰宽度大于预设的线谱宽度门限Th_W，转到步骤S8；

若谱峰宽度小于设定的线谱宽度门限Th_W，执行步骤S7；

S7、将螺旋桨匀速旋转与螺旋桨变速旋转的时间段进行区分：

对G(x，y)沿频率号x_M所代表的线谱按时间从前向后进行谱线跟踪，Y1取最先跟踪到线谱的行数，Y2取最后跟踪到线谱的行数，跟踪的实现方法为：对第y行及其两侧相邻的预设Y_W行中搜索G(x，y)在x_M±x_P范围内的非零值，其中x_P为允许频率的抖动频点数，若非零值的个数Y_N与Y_W+1的比值大于0.5，则认为第y行跟踪到线谱，否则认为第y行未跟踪到线谱；

若Y1≥0，Y2＝Y，变速时间区间识别过程结束；

若Y1≥0，Y2＝Y_R≥Y1，将G(x，y)的前Y_R行的频谱值置零，即G(x，y∈[0，Y_R])＝0，并置标志位W1＝Y_R，并返回步骤S5；

若Y2≤Y，0＜Y1＝Y_R≤Y2，将G(x，y)的后Y-Y_R行的频谱值置零，即G(x，y∈[Y_R+1，Y])＝0，并置标志位W2＝Y_R，并返回步骤S6；

S6、根据W1和W2的值计算螺旋桨加速或减速旋转时间段对应的起始采样点数K₁和结束采样点数K₂分别为：

$\left\{\begin{array}{c}{K}_1=F_S{t}_{D}W1+1\\ K_2=F_S(W2\·t_D+T_W)\end{array}---\left(8\right).\right.$

优选地，所述T_W取秒级，步进t_D取T_W/100，Th_W取频率分辨率1/T_W的3至6倍。

(三)有益效果

本发明基于船舶水下辐射的螺旋桨空化噪声受螺旋桨转速调制的特点，以及螺旋桨匀速旋转和变速旋转时调制包络信号的特点与区别，运用图像处理的方法对经恒虚警处理的调制包络谱瀑布图上代表其叶频分量的谱线的识别与判断，将螺旋桨匀速旋转的时间段与变速旋转的时间段予以区分，从而获得螺旋桨变速旋转的起止时刻。本发明提出的一种自动获取船舶螺旋桨变速旋转时间区间方法解决了无人值守的被动水声探测系统无法自动获取螺旋桨加速或减速旋转的时间段，进而无法获取螺旋桨完全处于变速旋转过程中辐射噪声数据的问题，并且未增加额外的硬件设备，从而为被动水声探测系统进一步正确提取螺旋桨旋转速度变化特征提供了保障。

附图说明

图1是本发明自动获取船舶螺旋桨变速旋转时间区间方法的总体实现流程图；

图2是本发明获取经频域CFAR处理的灰度调制包络谱图的实现流程；

图3是本发明自动识别螺旋桨变旋转时间段的信号处理实现流程；

图4是在两种试验工况下，螺旋桨旋转速度随时间变化曲线：(a)为试验数据1，(b)为试验数据2；

图5是在两种试验工况下，试验数据经0角度Radon变换后的输出：(a)为试验数据1，(b)为试验数据2。

具体实施方式

为使本发明的目的、内容、和优点更加清楚，下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。

本发明实施的前提条件是：被动水声探测系统已经检测到船舶目标，并能够对目标进行连续跟踪，且接收的数据中仅包含单个船舶的辐射噪声信号。同时，被动水声探测系统已获得在未出现船舶时，环境噪声信号的数据。

本发明的基本原理以船舶的加速过程为例进行叙述，减速过程与此类似。水面商用船舶在加速过程中，船舶螺旋桨会从一个较低的转速提高至设计航速。在转速提高前或在转速提高的过程中，螺旋桨的旋转速度会超过空化起始转速，从而产生受到螺旋桨转速调制的螺旋桨空化噪声，并成为被动水声探测系统的主要信号源。并且，螺旋桨在达到设计转速后，将会保持匀速旋转。匀速旋转螺旋桨产生空化噪声幅度调制包络信号在频域具有良好的能量聚集性，在调制包络谱的瀑布图上，调制包络信号的各次谐波分量表现为与频率轴相垂直的亮线。同时，由于螺旋桨空化噪声调制包络的叶频分量通常具有最大能量，因此可以运用图像处理方法对调制谱瀑布图上代表叶频分量的与频率轴相垂直的谱线进行判断和识别，将螺旋桨匀速旋转和变速旋转的时间段进行区分，再将匀速旋转的时间段剔除，从而获得螺旋桨变速旋转的时间区间。

如图1所示，本发明的方法包括以下几个步骤：

S1、对采样后的观测信号r(n)(r(n)为一个序列，n为正整数，总观测时间为T)进行带通滤波(用于去噪)、绝对值包络检波，并消除直流分量，以获取包络信号。绝对值检波的过程可描述为：

$E_0\left(l\right)=\frac{1}{M}\underset{j=1}{\overset{M}{\mathrm{\∑}}}\left|r(\mathrm{Mn}+j)\right|---\left(1\right)$

其中M是绝对值检波的平均点数，满足：

M≤int(F_s/f_h)                                  (2)

式中F_S为采样频率，f_h是带通滤波器的最高频率，int表示取整函数(以四舍五入的方式取整)。在绝对值检波后，采样率降为F_S/M。

消除直流分量的过程可描述为：

$E\left(l\right)=E_0\left(l\right)-\frac{1}{L}\underset{l=1}{\overset{L}{\mathrm{\∑}}}{E}_0\left(l\right)---\left(3\right)$

式中L是E₀(l)的点数。

S2、以宽度为T_W的滑动窗(T_W取秒级)对E(l)进行截取，而后进行2^K(K为预设整数)点快速傅里叶变换(Fast Fourier Transform，FFT)和频域恒虚警率(Constant False Alarm Rate，CFAR)处理，且经CFAR处理后超过阈值的点直接取幅度谱的原始谱值(即保留原值)，低于阈值的点置零。再对输出结果进行灰度转换，形成一行调制包络谱分析结果。频域CFAR处理的阈值Th由如下公式确定：

Th＝m+dσ                                    (4)

式中d是由恒虚警概率决定的检测因子，m为无船舶目标时步骤S1的包络解调处理后频谱幅度值的均值，σ为无船舶目标时步骤S1的包络解调处理后频谱幅度值的方差。

步骤S3、以步进时间t_D将滑动窗进行滑动，步进t_D可取T_W/100，重复步骤S2，获得多行调制包络谱分析结果，并形成瀑布图的形式，从而获得接收信号的经CFAR处理的灰度调制包络谱瀑布图，记为G(x，y)，其中x＝X₁，X₁+1，…，X_H和y＝1，2…Y分别代表了时间和频率轴上的变量：

$\left\{\begin{array}{c}{X}_1=\mathrm{int}(2^{K-1}M{F}_l/F_s)\\ X_H=\mathrm{int}(2^{K-1}M{F}_h/F_s)\\ H=X_H-X_1+1\\ Y=\mathrm{int}((T-T_W)/t_D)\end{array}---\left(5\right)\right.$

式中(F_l，F_h)为包络解调处理所分析的频率范围；

步骤S4、初始化标志位：Y1＝0，Y2＝Y，W1＝0，W2＝Y；

步骤S5、对G(x，y)进行0角度的Radon变换，并将Radon变换的横坐标s转换成频率号x，转换公式为：

$x=\mathrm{int}\left(\frac{2^{K-1}M(F_h+F_l)}{F_S}\right)+s---\left(6\right)$

从而可以得到随频率变化的0角度Radon变换输出R(x)；

步骤S6、对R(x)的峰值进行分析。

若R(x)＝0，说明当前尚未接收到目标产生的受调制的螺旋桨空化噪声，变速时间区间识别过程结束。

当R(x)≠0时，确定R(x)最大值的位置，即峰值位置的频率号x_M。并计算该峰值两侧相临的连续超过该峰值一半的频率点数N_P，以获得谱峰宽度：

$W_R=\frac{N_P{F}_S}{2^{K}M}---\left(7\right)$

若谱峰宽度大于预设的线谱宽度门限Th_W，Th_W取3至6倍的频率分辨率1/T_W，则该谱峰不是线谱，表明当前的G(x，y)所表征的时间范围内完全处于加速或减速的状态，并转到步骤S8。

若谱峰宽度小于设定的线谱宽度门限Th_W，说明形成调制包络谱图的时间内包含螺旋桨匀速旋转的过程，并执行步骤S7；

步骤S7、将螺旋桨匀速旋转与螺旋桨变速旋转的时间段进行区分。对G(x，y)沿频率号x_M所代表的线谱按时间从前向后进行谱线跟踪，Y1取最先跟踪到线谱的行数，Y2取最后跟踪到线谱的行数。跟踪的实现方法为：对第y行及其两侧相邻的Y_W行(为预设值)中搜索G(x，y)在x_M±x_P范围内的非零值，其中x_P为允许频率的抖动频点数，可以取x_P≤2，若非零值的个数Y_N与Y_W+1的比值大于0.5，则认为第y行跟踪到线谱，否则认为第y行未跟踪到线谱。

1)若Y1≥0，Y2＝Y，则表明当前形成调制包络谱图的时间内螺旋桨始终处于匀速旋转状态，变速时间区间识别过程结束。

2)若Y1≥0，Y2＝Y_R≥Y1，则表明形成调制包络谱图的时间内螺旋桨首先处于匀速旋转状态。将G(x，y)的前Y_R行的频谱值置零，即G(x，y∈[0，Y_R])＝0，并置标志位W1＝Y_R，并直接返回步骤S5。

3)若Y2≤Y，0＜Y1＝Y_R≤Y2，则表明形成调制包络谱图的时间内螺旋桨最终处于匀速旋转状态。将G(x，y)的后Y-Y_R行的频谱值置零，即G(x，y∈[Y_R+1，Y])＝0，并置标志位W2＝Y_R，并直接返回步骤S6；

步骤S8、根据W1和W2的值计算螺旋桨加速或减速旋转时间段对应的起始采样点数K₁和结束采样点数K₂分别为：

$\left\{\begin{array}{c}{K}_1=F_S{t}_{D}W1+1\\ K_2=F_S(W2\·t_D+T_W)\end{array}---\left(8\right)\right.$

下面以在符合国际标准的大型空泡循环水洞进行的加速旋转螺旋桨噪声测试试验获得的试验数据为实例，对本发明的具体实现方式进行说明，试验用目标螺旋桨是为某集装箱船设计的五叶缩比模型桨，缩尺比28，通过整船模型形成周向非定常的尾流场，整船模型长约7.2m。取两种不同工况下的试验数据进行分析，两组试验工况的参数如表1所示。试验中螺旋桨转速变化范围根据水洞试验相似准则的要求设计，转速随时间的变化曲线如附图4所示。

表1两组试验工况参数

在此基础上，实现步骤归纳如下：

(1)在不同工况下的观测时间内，均取滑动窗长度T_W＝2.261s，步进时间t_D＝0.0226，绝对值检波的平均点数M＝30，调制包络信号的频率分析范围F_l＝80，F_h＝180Hz，进行8192点的快速傅里叶变换，根据步骤S1至步骤S3，获得两种试验工况下，试验数据经频域CFAR处理的灰度调制包络谱图G(x，y)，如附图4所示；

(2)初始化标志位：

工况1：Y1＝0，Y2＝342，W1＝0，W2＝342，Y＝342

工况2：Y1＝0，Y2＝165，W1＝0，W2＝165，Y＝165

(3)对G(x，y)进行0角度的Radon变换，得到随频率变化的0角度Radon变换输出R(x)，两组试验工况下R（x)；

(4)对R(x)的峰值进行分析。两组不同试验工况下，R(x)输出如附图5所示。

工况1：R(x)≠0，计算得峰值位置的频率序号x_M＝278，N_P＝2，因此谱峰宽度W_R＝≤Th_W＝6/T_W，直接进行步骤S7；

工况2：R(x)≠0，计算得峰值位置的频率序号x_M＝257，N_P＝10，因此谱峰宽度W_R＝≥Th_W＝6/T_W，直接进行步骤S8；

(5)对于工况1，将螺旋桨匀速旋转与螺旋桨变速旋转的时间段进行区分。进行谱线跟踪得到：Y2≤Y＝342，0＜Y1＝Y_R＝199≤Y2。将G(x，y)的后Y-Y_R＝143行的频谱值置零，即G(x，y∈[200，342])＝0，并置标志位W2＝199，重复步骤S5；

(6)对于工况1，再次对R(x)的峰值进行分析，有R(x)≠0，计算得N_P＝13，谱峰宽度W_R＝≥Th_W＝6/T_W，直接进行步骤S8，且使W2＝199；

(7)计算螺旋桨加速或减速旋转时间段对应的起止采样点数。

工况1：K₁＝1，K₂＝734608；

工况2：K₁＝1，K₂＝651266；

需要说明的是，对于工况1，尽管在最开始的一小段时间内，螺旋桨处于匀速旋转的过程，但由于此时螺旋桨转速尚未达到空化起始转速，因而不产生螺旋桨空化噪声，即该段时间内的数据在G(x，y)和R(x)上均无体现，因此运用本发明的步骤不会识别出开始时刻螺旋桨匀速旋转过程。但由于在这些时段内无螺旋桨空化噪声的调制包络信号，在进行分数阶傅里叶变换时不会产生调频率不匹配的问题，因而对于从调制包络信号中提取转速变化特征没有影响，即对于没有产生空化噪声的螺旋桨匀速旋转时间区间无需剔除。

根据变速旋转时间段对应的起止采样点数，可求出工况1下螺旋桨达到最终转速时的时刻分别为6.7294s，并根据R(x)峰值位置的频率序号x_M可求出工况1下螺旋桨达到的最终转速为24.8Hz，再将计算得到的最终转速及其起始时刻用虚线分别标示于螺旋桨转速-时间曲线上，如附图4所示，以更加直观地检验对变速时间段估计的正确性。

通过以上步骤，在两种不同的工况下，自动识别了螺旋桨加速旋转的时间段，并由附图4可以看出，对螺旋桨变速旋转时间段的估计基本准确，说明本发明所提出的自动获取船舶螺旋桨变速旋转时间区间方法解决了被动水声监测系统无法自动获取螺旋桨变速旋转时间段的问题，进而为被动水声探测系统对螺旋桨旋转速度变化特征的准确提取提供了保障。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。

Claims (2)

1.一种自动获取船舶螺旋桨变速旋转时间区间方法，其特征在于，包括以下步骤：

S1、对采样后的观测信号r(n)依次进行带通滤波、绝对值包络检波，并消除直流分量，得到包络信号，其中，n为正整数，所述绝对值包络检波的过程为：

$E_0\left(l\right)=\frac{1}{M}\underset{j=1}{\overset{M}{\mathrm{\∑}}}\left|r(\mathrm{Mn}+j)\right|---\left(1\right)$

其中M是绝对值检波的平均点数，满足：

M≤int(F_s/f_h)                                   (2)

式(2)中F_S为采样频率，f_h是带通滤波器的最高频率，int表示取整，在绝对值检波后，采样率降为F_S/M；

消除直流分量的过程为：

$E\left(l\right)=E_0\left(l\right)-\frac{1}{L}\underset{l=1}{\overset{L}{\mathrm{\∑}}}{E}_0\left(l\right)---\left(3\right)$

式中L是E₀(l)的点数；

S2、以宽度为T_W的滑动窗对E(l)进行截取，然后依次进行2^K点快速傅里叶变换和频域恒虚警率处理，K为正整数，经频域恒虚警率处理后超过阈值的点保留，低于阈值的点置零，再对输出结果进行灰度转换，形成一行调制包络谱分析结果，频域恒虚警率处理的阈值Th由公式(4)确定：

Th＝m+dσ                                         (4)

式(4)中d是由恒虚警概率决定的检测因子，m为无船舶目标时经步骤S1处理后频谱幅度值的均值，σ为无船舶目标时经步骤S1处理后频谱幅度值的方差；

S3、以步进时间t_D将滑动窗进行滑动，重复步骤S2，获得多行调制包络谱分析结果，并形成瀑布图的形式，从而获得包络谱瀑布图，记为G(x，y)，其中x＝X₁，X₁+1，…，X_H和y＝1，2…Y分别代表时间和频率轴上的变量：

$\left\{\begin{array}{c}{X}_1=\mathrm{int}(2^{K-1}M{F}_l/F_s)\\ X_H=\mathrm{int}(2^{K-1}M{F}_h/F_s)\\ H=X_H-X_1+1\\ Y=\mathrm{int}((T-T_W)/t_D)\end{array}---\left(5\right)\right.$

式(5)中(F_l，F_h)为步骤S1的处理过程所分析的频率范围，T为观测信号r(n)的总观测时间；

S4、初始化标志位：Y1＝0，Y2＝Y，W1＝0，W2＝Y；

S5、对G(x，y)进行0角度的Radon变换，并将Radon变换的横坐标s转换成频率号x，转换公式为：

$x=\mathrm{int}\left(\frac{2^{K-1}M(F_h+F_l)}{F_S}\right)+s---\left(6\right)$

从而得到随频率变化的0角度Radon变换输出R(x)；

S6、对R(x)的峰值进行分析：

若R(x)＝0，变速时间区间识别过程结束；

当R(x)≠0时，确定R(x)最大值的位置，即峰值位置的频率号x_M，并计算该峰值位置两侧相临的连续超过该峰值一半的频率点数N_P，以获得谱峰宽度：

$W_R=\frac{N_P{F}_S}{2^{K}M}---\left(7\right)$

若谱峰宽度大于预设的线谱宽度门限Th_W，转到步骤S8；

若谱峰宽度小于设定的线谱宽度门限Th_W，执行步骤S7；

S7、将螺旋桨匀速旋转与螺旋桨变速旋转的时间段进行区分：

对G(x，y)沿频率号x_M所代表的线谱按时间从前向后进行谱线跟踪，Y1取最先跟踪到线谱的行数，Y2取最后跟踪到线谱的行数，跟踪的实现方法为：对第y行及其两侧相邻的预设Y_W行中搜索G(x，y)在x_M±x_P范围内的非零值，其中x_P为允许频率的抖动频点数，若非零值的个数Y_N与Y_W+1的比值大于0.5，则认为第y行跟踪到线谱，否则认为第y行未跟踪到线谱；

若Y1≥0，Y2＝Y，变速时间区间识别过程结束；

若Y1≥0，Y2＝Y_R≥Y1，将G(x，y)的前Y_R行的频谱值置零，即G(x，y∈[0，Y_R])＝0，并置标志位W1＝Y_R，并返回步骤S5；

若Y2≤Y，0＜Y1＝Y_R≤Y2，将G(x，y)的后Y-Y_R行的频谱值置零，即G(x，y∈[Y_R+1，Y])＝0，并置标志位W2＝Y_R，并返回步骤S6；

S8、根据W1和W2的值计算螺旋桨加速或减速旋转时间段对应的起始采样点数K₁和结束采样点数K₂分别为：

$\left\{\begin{array}{c}{K}_1=F_S{t}_{D}W1+1\\ K_2=F_S(W2\·t_D+T_W)\end{array}---\left(8\right).\right.$

2.如权利要求1所述的方法，其特征在于，所述T_W取秒级，步进t_D取T_W/100，Th_W取频率分辨率1/T_W的3至6倍。

Images