海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法技术

一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，属于海浪监测技术领域。为解决现有的海浪信息反演技术用采的技术手段单一导致反演的海浪参数较少，获得的海浪信息的精度较低的问题而提出的。本发明专利技术将海浪谱分析法和雷达图像几何阴影统计法结合起来，可以在目标海域波浪浮标缺失或故障的情形下反演出有效波高、波峰峰向和波峰周期。利用试验场实测数据反演得到海浪信息，并与同海域内投放的波浪浮标数据进行实验误差分析，结果表明海浪谱分析法和雷达图像几何阴影统计法相结合的海浪信息反演技术达到预期精度。本发明专利技术的海浪信息反演技术用采的技术手段合理且反演的海浪参数较多，从而获得的海浪信息的精度较高。获得的海浪信息的精度较高。获得的海浪信息的精度较高。

【技术实现步骤摘要】
海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法


[0001]本专利技术涉及一种利用海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，属于海浪监测


技术介绍

[0002]目前，基于X波段航海雷达图像进行海浪信息反演的主流技术有两类，分别为基于海浪谱分析法的海浪信息反演技术和基于雷达图像几何阴影图像统计法的海浪信息反演技术。海浪谱分析法反演海浪信息的精度更高但应用条件受限，需要外部设备如浮标观测的实测数据来进行拟合；几何阴影图像统计法可以无需外部设备直接对雷达图像进行处理来反演有效波高。两种技术各有侧重，相辅相成实现试验场全海域覆盖的海浪观测。
[0003]基于海浪谱分析法的海浪信息反演技术首先得到发展。利用遥感技术对海浪进行观测始于二十世纪六十年代。直到1985年，Young和Rosenthal等人首次引入三维傅里叶变换的方法应用于雷达海杂波图像序列，在三维的角度上分析了海浪方向和海表面流的特征，同时考虑了海浪信号在时间上和空间上的变化特征，利用海浪谱分析的方法实现对海浪参数的反演。后来学者Nieto等人将应用在合成孔径雷达中反演有效波高的方法引入到X波段雷达中来，该方法基本思想是有效波高和雷达图像信噪比开方之间线性相关，通过计算得到的信噪比与浮标数据实测值进行线性拟合，得到经验参数。至此基于海浪谱分析法的海浪信息反演技术逐渐成熟，并成为了X波段雷达反演海浪参数的主要理论，后来的学者开始研究如何继续提高算法的精度。2002年，Rune针对色散关系带通滤波器进行改进，提出了一种基于加权最小二乘的滤波器，结果表明该滤波器的应用有效提高了有效波高的反演精度。王淑娟、王剑等人先后提出了最小二乘拟合和分段拟合的方法标定有效波高的线性关系参数，从而提高了海浪参数反演的精度。针对降雨和同频干扰噪声的影响，唐艳红改进了噪声去除及滤波算法，准确提取了海浪频谱信息，进一步提高了反演精度。李英通过野值剔除和稀疏化方法对实验数据进行处理，并提出了基于PSO的分段拟合方法，提高了波高反演拟合的精度。基于X波段雷达反演海浪参数多在短脉冲下采集的数据进行处理，Ludeno等人研究了对中脉冲下采集的数据进行反演的可行性，结果表明两种情况下的结果一致。Al
‑
Habashneh等人对笛卡尔坐标系下的雷达图像进行了研究。Wendy等人提出了一种基于滤波和插值的新方法，通过X波段雷达估计浅水区域的海浪参数。针对可能受到其他传播和雷达的电磁干扰的问题，Pavel等人提出了一种滤除伪波能量算法，结果表明去除干扰特征显著改善了海浪反演精度。
[0004]另一类海浪信息反演技术是基于雷达图像几何阴影的图像统计法。1990年，L.B.Wetzel首次引入粗糙面几何光学原理，为了将海面雷达图像分割为阴影区和非阴影区，对低掠射角入射下的海面后向散射模型进行研究，选取灰度阈值并计算阴影面积比例来反演有效波高。2001年，J.R Buckley针对阴影分割阈值与径向距离之间关系的问题，提出随着雷达径向距离的变化，也即掠射角在1
°
～10
°
内变化时，阴影分割阈值的变化不大，
也就是说应该在1
°
～10
°
的低掠射角条件下分析雷达图像。2014年，Rune针对实际应用中雷达图像阴影难以辨别的问题，利用边缘检测思想对雷达图像的阴影区和非阴影区进行分割，并计算得到阴影比例，通过Smith函数拟合在无需数据校准的条件下反演得到有效波高。至此基于雷达图像几何阴影的图像统计法基本成型，后来的学者在阴影分割的问题上不断深入研究，使得该技术发展成熟。卢志忠等人针对传统雷达图像固定分块法的准确性问题，提出了一种基于差分边缘检测思想的自适应分块方法，提高了有效波高反演精度。而后卢志忠等人又提出了基于风夹角特征的反演方法，结果验证了该方法的有效性。基于原始雷达图像和未校准波仰角图像之间的相关性，Ludeno在Rune方法基础上提出一种新的波高无校准反演方法。卫延波在利用波陡计算有效波高时考虑了海表面流速的影响，提高了有效波高的反演精度。姚智超等人针对利用差分算子的方法估计阴影分割阈值时，其准确度和稳定性不高的问题，采用Prewitt算子进行改进，结果表明该算子相较于差分算子准确度更高。
[0005]由于欧美国家在海浪信息监测领域发展较早，现在已经开发出商业化应用的基于X波段航海雷达的测浪系统，其中包括挪威Miros公司开发的WAVEX系统和德国GKSS公司开发的WaMoS系统，实现了对海浪和海流信息进行观测。
[0006]目前国内还未出现成熟商用的基于X波段航海雷达的测浪系统，大多进口了国外设备系统，我国需要开发具有独立知识产权的海浪监测设备，亟需将两种反演技术结合，在以试验场为依托的大规模试验应用方面进行验证。
[0007]现有的海浪信息反演技术用采的技术手段单一导致反演的海浪参数较少，获得的海浪信息的精度较低。

技术实现思路

[0008]本专利技术要解决的技术问题为：现有的海浪信息反演技术用采的技术手段单一导致反演的海浪参数较少，获得的海浪信息的精度较低。
[0009]本专利技术为解决上述技术问题所采用的技术方案为：基于海浪谱分析法和雷达图像几何阴影统计法相结合的海浪反演技术的基本流程为：
[0010]步骤1、获取海杂波数据，海杂波数据为X波段航海雷达发射信号入射海浪时发生后向散射的回波数据，并对数据进行预处理，预处理包括数字下变频和坐标变换；
[0011]步骤2、选取雷达图像序列；
[0012]步骤3、利用三维傅里叶变换的方法变换得到三维波数频率能量谱；
[0013]步骤4、然后通过基于色散关系的带通滤波器滤除同频干扰、背景噪声信号；
[0014]步骤5、对三维波数频率能量谱的海浪角频率ω进行积分得到二维图像谱I(k
x
,k
y
)，利用调制传递函数对不同波数条件下的图像谱进行修正，修正后可以直接计算海浪能量，进而计算得到海浪信噪比；
[0015]步骤6、利用谱变换方法得到波峰峰向和波峰周期。
[0016]步骤7、根据步骤1采集的原始导航雷达图像信息判断是否降雨；
[0017]步骤8、剔除降雨时的雷达图像，在不含降雨信息的图像中选取扇形数据分析区域，进行海浪波高反演预处理；
[0018]步骤9、选取合适的雷达图像分析区域并对雷达数据进行方位向稀疏化预处理；
[0019]步骤10、随后对预处理后的雷达图像分析区域进行基于Canny算子的边缘检测，得到边缘梯度图像并建立灰度统计直方图，在阈值处理中，使用Otsu方法的最优全局阈值处理寻优得到阴影分割阈值，进而得到阴影分割图像；
[0020]步骤11、在方位向和径向上进行分块，并计算各分块的阴影比例；
[0021]步骤12、随后基于史密斯函数对曲线进行拟合来估计均方根波陡；
[0022]步骤13、最终利用有效波高与均方根本文档来自技高网...

【技术保护点】

【技术特征摘要】
1.一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，所述方法的实现过程为：步骤1、获取海杂波数据，海杂波数据为X波段航海雷达发射信号入射海浪时发生后向散射的回波数据，并对数据进行预处理，预处理包括数字下变频和坐标变换；步骤2、选取雷达图像序列；步骤3、利用三维傅里叶变换的方法变换得到三维波数频率能量谱；步骤4、然后通过基于色散关系的带通滤波器滤除同频干扰、背景噪声信号；步骤5、对三维波数频率能量谱的海浪角频率ω进行积分得到二维图像谱I(k
x
,k
y
)，利用调制传递函数对不同波数条件下的图像谱进行修正，修正后可以直接计算海浪能量，进而计算得到海浪信噪比；步骤6、利用谱变换方法得到波峰峰向和波峰周期；步骤7、根据步骤1采集的原始导航雷达图像信息判断是否降雨；步骤8、剔除降雨时的雷达图像，在不含降雨信息的图像中选取扇形数据分析区域，进行海浪波高反演预处理；步骤9、选取合适的雷达图像分析区域并对雷达数据进行方位向稀疏化预处理；步骤10、随后对预处理后的雷达图像分析区域进行基于Canny算子的边缘检测，得到边缘梯度图像并建立灰度统计直方图，在阈值处理中，使用Otsu方法的最优全局阈值处理寻优得到阴影分割阈值，进而得到阴影分割图像；步骤11、在方位向和径向上进行分块，并计算各分块的阴影比例；步骤12、随后基于史密斯函数对曲线进行拟合来估计均方根波陡；步骤13、最终利用有效波高与均方根波陡之间的关系实现对有效波高的反演。2.根据权利要求1所述的一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，在雷达回波数据预处理中数字下变频具体为：采用X波段航海雷达中频信号中心频率为135MHz，带宽为25MHz；首先通过A/D转换器对中频模拟信号直接采样，根据Nyquist带通采样定理，采样频率设计为60MHz；随后由数字相干振荡器产生两路频率为135MHz的数字本振信号cos(nω
c
)和sin(nω
c
)，两路信号相位相差90
°
；然后对FIR低通滤波器进行设计，根据系统实际需求选择Hamming窗进行滤波器设计，通带截止频率为13MHz，阻带截止频率为17MHz，阶数为60阶；将采样信号与数字本振信号进行数字混频，并分别通过FIR低通滤波器得到基带同向数字信号I(n)和基带正交数字信号Q(n)进行后续处理；在雷达回波数据预处理中坐标变换具体为：使用最近邻插值的方法来进行坐标转换，利用各点的极坐标来绘制雷达图像，任一点的极坐标可以表示为(r,θ,z)，转换到笛卡尔坐标系下该点的坐标可以表示为(x,y,z)，其中z表示该点的回波强度，二维平面坐标关系如下式所示：定义矩形框中任一点的平面坐标为(x0,y0)，采用最近点插值的方法可以得到在极坐标
中与(x0,y0)点相距最近的一点(r0,θ0)，同样的方法把矩形框中的每一点都对应得到极坐标系中的一个点，可以表示为下式所示：式中round(g)——向极坐标上最近的点取整；rem(g)——取余函数；随后将极坐标(r0,θ0)的回波强度值赋值给对应的矩形框中的点(x0,y0)，完成坐标变换；雷达图像数据经坐标变换以便于后续海浪反演算法的分析。3.根据权利要求1或2所述的一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，在步骤3中，利用三维傅里叶变换的方法变换得到三维波数频率能量谱，具体为：对笛卡尔坐标系下的雷达图像序列η(x,y,t)应用三维傅里叶变换，得到三维波数频率图像谱F(k
x
,k
y
,ω)，如下式所示：式中k
x
——x方向波数(rad/m)；k
y
——y方向波数(rad/m)；假定空间上选取在大小为L
x
*L
y
的矩形区域，时间上选取时间长度为T的时间序列，来进行傅里叶变换分析，如下式所示：式中ω——海浪角频率(rad/s)；η(,x,y)t表示笛卡尔坐标系下的雷达图像序列；其离散形式如下式所示：又根据傅里叶变换的原点对称性，即F(k
x
,k
y
,ω)＝F(
‑
k
x
,
‑
k
y
,
‑
ω)于是，三维波数频率能量谱可以表示为：4.根据权利要求3所述的一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，在步骤4中，所述基于色散关系的带通滤波器为：根据海浪特性设计滤波器对噪声信号进行滤除；色散关系方程是关于角频率和波数的关系方程，海表面波的周期和波长近似满足这一关系，色散关系方程表示为：式中tanh(g)——双曲正切函数；ω——海浪角频率(rad/s)；
g——重力加速度(m/s2)；——波数向量模(m)；h——水深(m)；——相对速度向量；海浪信号符合重力波的色散关系，利用色散关系来构造滤波器，从图像谱中去除噪声信号并将海浪信号提取出来；定义基于色散关系的带通滤波器为：其中，其中，ω
p
和ω
n
为角频率ω的截止边界。5.根据权利要求4所述的一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，在步骤5具体过程为：完成基于色散关系的带通滤波器设计后，对三维波数频率能量谱I(k
x
,k
y
,ω)进行滤波处理，并对ω进行积分得到二维图像谱I(k
x
,k
y
)，积分过程表示为：I(k
x
,k
y
)＝∫
ω＞0
E(k
x
,k
y
,ω)dω对不同波数条件下的图像谱进行修正，修正后可直接计算海浪能量；得到二维海浪谱的过程表示为：E(k
x
,k
y
)＝|M(k
x
,k
y
)|2·
I(k
x
,k
y
)|M(k
x
,k
y
)|2≈k
‑
β
式中|M(k
x
,k
y
)|2——经验函数；β——经验系数；|M(k
x
,k
y
)|2表示当波数不同时，二维海浪谱与二维图像谱之间的比例关系。6.根据权利要求5所述的一种海浪谱分析法和基于Canny算子的雷达图像几何阴影统计法相结合的海浪参数反演方法，其特征在于，在步骤6中，利用谱变换方法得到波峰峰向和波峰周期，波峰周...

【专利技术属性】
技术研发人员：金涛，王搏，冯硕，卢嘉乐，魏嘉凝，周洪娟，周志权，王晨旭，
申请(专利权)人：哈尔滨工业大学威海，
类型：发明
国别省市：