一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法
技术领域
本发明属于海洋监测技术领域,具体涉及一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法。
背景技术
在诸多海洋环境因素中起主导作用的是海浪,海浪是发生在海洋表面的一种波动现象,是海洋监测要素中最重要且最复杂的一种自然现象[1]。海浪监测是人类认识、研究、开发利用海洋的重要任务之一,现代航海、造船、军事、海洋工程、灾害监测预警等领域,都需要获取高精度的实时海浪资料,用来描述海浪特征的量,如波高、波向、周期、相速、波长等,统称为海浪要素[2]。波向是海浪的重要参数之一,波向测量研究对于发展海洋工程和海洋科学来说都是十分重要的[3]。见参考文献[1-3](王言英,徐继文,赫亮.海浪能量谱密度函数计算方法[J].哈尔滨工程大学学报,2011,32(10):1283-1289.齐勇,闫星魁,郑姗姗,等.海浪监测技术与设备概述[J].气象水文海洋仪器,2015,33(3):113-117.Wang L,ChengY,Hong L,et al.Anovel algorithm for ocean wave direction inversion fromX-band radar images based on optical flow method[J].Acta Oceanologica Sinica,2018,37(3):88-93)。
目前国内外进行波向数据观测主要有两类方法:第一类采用遥感技术,其中有代表性的是利用X波段航海雷达海浪回波反演计算波向等信息[4]。第二类采用接触式传感器,目前的技术主流是利用波浪浮标测量波向,这是一种自动、定点、定时的测量系统,具有更高的准确性和置信度[5]。见参考文献[4-5](Guozhen,Qingyou,GUAN,et al.Acurvelet-based method to determine wave directions from nautical X-band radar images[J].Acta Oceanologica Sinica,2018,37(1):11-19.刘建国,周辉云,王晓亮,etal.SZF2-1A型多参数波浪试验浮标波浪要素海上比测试验与初步分析[J].海洋技术学报,2016,35(3):74-77)。
国内外公开文献上能够查到的利用波浪浮标进行波向测量的仪器原理上主要分成如下两类[6][7]:第一类是利用浮标体内的三轴加速度传感器和电子罗盘,测得浮标体的方位角和水平两个垂直方向的瞬时位移,进而计算出主波向,国产的SBF型波浪浮标即属于该类型;第二类是利用浮标内的波高倾斜一体化传感器和方位传感器,分别测得浮标体的瞬时方位角、纵倾角和横倾角信息,进而计算出主波向,国产的SZF型波浪浮标即属于该类型。见参考文献[6-7](唐原广,康倩.波浪浮标测波方法比较[J].现代电子技术,2014,v.37;No.422(15):121-122.毛祖松.我国近海波浪浮标的历史、现状与发展[J].海洋技术学报,2007,26(2):23-27)。
针对第一类工作原理的浮标,目前能够从公开文献上查到的测量波向方法如下:以国产的SBF型波浪浮标为例,该浮标是由山东省科学院海洋仪器仪表研究所研制完成的,测量波向方法是采用浮标体内的三轴加速度传感器和电子罗盘相结合的方式,在一个波浪周期内的跨零点处,计算浮标的水平矢量合成加速度方向作为单个波浪周期的参考方向。以加速度x轴正轴和y轴正轴建立相对正交坐标系,计算出水平矢量合成加速度与x轴的夹角即为相对波向,再结合电子罗盘计算得到的方位角就可得到波浪的地理波向[8][9]。见参考文献[8-9](赵杰,惠力,初士博,等.基于三轴加速度的波浪测量技术研究[J].海洋技术学报,2015,34(5):66-70.刘国栋.波浪浮标数据处理方法研究[D].天津大学,2012:49-50)。
国内针对第二类工作原理的浮标,目前能够从公开文献上查到的测量波向方法如下:以国产的SZF型波浪浮标为例,该浮标是由中国海洋大学研制生产的,测量波向方法与国外文献公布的方法类似,首先根据上跨零点来确定单个波浪周期,再确定每一个单个波浪周期的下跨零点,利用浮标测得的这一时刻瞬时方位角、纵倾角和横倾角信息,计算得到该单个波浪周期的参考波向;然后把各个波浪周期测得的单个波浪波向进行归并统计,将360°分成16个角度区间,统计16个角度区间中每个区间内单个波浪波向出现的个数,进而得到波向出现率,最后,选择波向出现率最大的角度区间的中值作为主波向10][11][12][13][14]。见参考文献[10-14](唐原广,王金平.SZF型波浪浮标系统[J].海洋技术学报,2008,27(2):31-33.王娟.波浪浮标传感器综合测试系统[D].中国海洋大学,2010:28-30.王金平.SZF型波浪浮标数据采集系统及数据处理软件设计[D].中国海洋大学,2008:14-15.张莹.SZF型波浪浮标数据采集、处理与无线传输系统[D].中国海洋大学,2010:8-9.贺成柱.海浪谱估计方法研究及SZF型波浪浮标数据分析软件设计[D].中国海洋大学,2007:42-43)。
国外针对第二类工作原理的浮标,目前能够从公开文献上查到的测量波向方法与上述国内相关文献的方法基本一致,最早查到的文献是由美国国家数据浮标中心公布了一种采用横倾、纵倾传感器进行波浪波向测量的研究工作,但未公布具体的算法[15]。测量波向的方法由Kenneth Steele,DWWang等人首次公布,该方法首先测得浮标的方位角、横倾和纵倾,然后根据这些参数计算浮标体的北向和东向斜度,最后得到波向[16][17][18]。通过比对,国内针对第二类原理浮标的波向测量方法与Kenneth Steele方法类似。见参考文献[15-18](Steele K,Lau J,Hsu Y H.Theory and application of calibrationtechniques for an NDBC directional wave measurements buoy[J].IEEE Journal ofOceanic Engineering,1985,10(4):382-396.Steele K E.Pitch-roll buoy mean wavedirections from heave acceleration,bow magnetism,and starboard magnetism[J].Ocean Engineering,2003,30(17):2179-2199.Steele K E,Wang D W,M.D.Earle,etal.Buoy pitch and roll computed using three angular rate sensors[J].CoastalEngineering,1998,35(1–2):123-139.Steele K E,Earle M D.Directional ocean wavespectra using buoy azimuth,pitch,and roll derived from magnetic fieldcomponents[J].IEEE Journal ofOceanic Engineering,1991,16(4):427-433)。
发明内容
本发明的目的在于提供一种可靠性高、准确性高的基于纵横倾角测量野值修正的波浪浮标主波向计算方法。
一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,具体包括以下步骤:
步骤1、野值修正:读取浮标测量记录的一组连续原始数据,包含瞬时波高、方位角、纵倾角、横倾角等,利用筛选和替补方法对浮标测得的原始纵倾角和横倾角数据进行野值处理,提高数据的准确度;
步骤2、数据处理:根据浮标测得的波高信息,利用上跨零点将野值修正后的该组数据切分成若干个单个波浪周期,并确定单个波浪周期的总个数;
步骤3、计算单个波浪波向:对每一个切分好的单个波浪周期,首先选取参与计算单个波浪波向的采样点区间,然后利用瞬时方位角以及筛选和替补方法处理后的纵倾角、横倾角信息,计算出该区间内每个采样点的瞬时波向,再通过分类和矢量平均方法对这些采样点的瞬时波向数据进行处理,从而计算得到单个波浪波向。重复该步骤,直到得到该组数据中所有的单个波浪波向;
步骤4、统计上述单个波浪波向测量结果在16个角度区间的计数比例:先将360°划分成16个角度区间,统计16个角度区间中每个区间内单个波浪波向出现的个数,进而得到计数比例;
步骤5、计算主波向:选择计数比例值最大的角度区间的中值作为主波向。
所述一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,步骤1包括以下步骤:
步骤1.1、读取波浪浮标在海上进行海浪监测时测得的一组连续实验数据(包含瞬时波高、方位角、纵倾角、横倾角等),然后先利用筛选方法选出波浪浮标测量纵倾角和横倾角产生的野值,筛选方法为:
遍历一组实验数据的所有纵倾角和横倾角数据值,确定超出-20°至20°范围的数据为野值,进行替补;
步骤1.2、经过步骤1.1筛选得到野值
对野值数据进行替补,替补方法为:
针对一组数据,数据序列为X
1,...,X
2048,从第一个野值开始,针对任一野值X
oi(i=1,...,N
o),选择此野值前M个数据值X
oi-M,...,X
oi-1,其权重为
再选择此野值后M个数据值X
oi+1,...,X
oi+M,其权重为
若选择的X
oi+1,...,X
oi+M中有超出-20°至20°范围的数据值,则选择超出范围数据值后一个数据值,保证选取的数据值不超出-20°至20°的范围;
则该野值的替补值的计算公式为:
通常M=2,则该野值的替补值的计算公式为:
利用上述替补方法处理完纵倾角和横倾角数据中所有的野值,得到处理后的纵倾角和横倾角信息。
所述一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,步骤2包括以下步骤:
步骤2.1、根据浮标测得的波高信息,采用上跨零点法,把波面上升与零点的交点,即瞬时波高数据由小于0变至大于0的点作为起点,当波浪运行到零点以下,即瞬时波高出现了小于0的情况后,将其后再次波浪上升与零点的相交点作为终点,则起点与终点之间为一个单个波浪周期;
步骤2.2、重复步骤2.1的方法,从第一个上跨零点开始,将该组实验数据切分成若干个单个波浪周期,并统计出该组数据中单个波浪周期的个数为N。
所述一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,步骤3包括以下步骤:
步骤3.1、根据步骤2切分好的单个波浪周期,在一个波浪周期内,选择波高数据大于0m的范围,作为参与计算单个波浪波向的采样点区间;
步骤3.2、对步骤3.1确定的采样点区间内的所有采样点,重复下列操作直至获得所有采样点的瞬时波向:利用瞬时方位角以及筛选和替补方法处理后的纵倾角、横倾角信息,先计算出采样点瞬时东向斜度和北向斜度,再计算出该采样点的瞬时波向;
式中:A为方位角,P为纵倾角,R为横倾角;
式中:A为方位角,P为纵倾角,R为横倾角;
采样点瞬时波向θp的计算公式为:
式中:
为浮标体采样点瞬时东向斜度,
为浮标体采样点瞬时北向斜度;
步骤3.3,对步骤3.2计算得到的采样点瞬时波向在空间域上进行分类,分类方法为:
步骤3.31、设置固定夹角为θc的扇形作为初试分类区,统计采样点瞬时波向值落在0~θc区间的个数U,θc取90°~180°;
步骤3.32、分类区从0°开始,顺时针旋转并保持夹角仍为θc,每旋转1°统计一次该分类区中采样点瞬时波向的个数值U并记录类中的采样点瞬时波向值,旋转360°后,得U1,...,U360;
步骤3.33、从U1,...,U360选出最大值UM,则此最大值所对应分类区中的采样点瞬时波向即为分类后的采样点瞬时波向θp1,...,θpa(a=1,2,3...UM);
步骤3.4、利用矢量叠加法处理步骤3.3获得的分类后的采样点瞬时波向,可得到瞬时波向北向分量ΔθN和瞬时波向东向分量ΔθE,进而计算出单个波浪波向θw;
瞬时波向北向分量ΔθN的计算公式为:
式中:θpi为分类后的第i个采样点瞬时波向;
瞬时波向东向分量ΔθE的计算公式为:
式中:θpi为分类后的第i个采样点瞬时波向;
单个波浪波向θw的计算公式为:
式中:ΔθE为瞬时波向东向分量,ΔθN为瞬时波向北向分量;
步骤3.5,重复步骤3.1、步骤3.2、步骤3.3和步骤3.4,直到得到该组数据中所有的单个波浪波向θw1…θwN,式中N为切分的单个波浪周期个数。
所述一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,步骤4包括以下步骤:
步骤4.1、以22.5°为间隔将360°均匀分成16个角度区间,以地理正北方位为零度,按照顺时针方向来划分;
步骤4.2,统计步骤3得到的所有单个波浪波向θw1…θwN在16个角度区间内分别出现的个数n1…n16;
步骤4.3,计算波浪波向在16个角度区间的计数比例;
计数比例pi的计算公式为:
式中:ni为单个波浪波向在第i个角度区间出现的个数,N为波浪周期的总个数。
所述一种基于纵横倾角测量野值修正的波浪浮标主波向计算方法,步骤5包括以下步骤:
步骤5.1、根据步骤4得到的计数比例,得到计数比例最大值K;
步骤5.2、将计数比例最大值K所对应的角度区间的中值作为主波向θ。
本发明的有益效果在于:
本发明通过传感器测得的载体原始方位、横倾和纵倾数据,对Kenneth Steele和唐原广团队等提出的方法进行了改进,实现了一种基于纵横倾角测量野值剔除的波浪浮标主波向计算方法。由于实际海洋中的海浪具有强随机性,而波浪浮标在海上工作的时候,处于一个十分复杂的海洋环境,会受到除波浪作用力以外的各种外力因素的影响,所以波浪浮标测得的纵倾角和横倾角信息有时会出现野值。国内波浪浮标对于海上工作出现的野值未曾进行深入研究,因此,国内波浪浮标测得的波向与真实波向相比具有较大的误差。本发明利用筛选和替补方法处理测得的纵倾角和横倾角的野值,并通过选取单个波浪过程中的多点波向数据进行融合处理的方法,有效提高了浮标测量波向的准确度。采用国产的SZF型波浪浮标在海上进行了现场比测试验,结果表明本发明计算的结果与现场实测波向结果的吻合度更高,改进方法优于浮标自身计算输出结果,即唐原广团队提出的波向测量方法,提高了主波向测量的可靠性和准确性。
附图说明
图1为浮标实测的原始纵倾角数据;
图2为浮标实测的原始横倾角数据;
图3为处理后的纵倾角数据;
图4为处理后的横倾角数据;
图5为上跨零点划分波浪周期方法;
图6为某一个波浪周期的实测波高数据;
图7为计算的某一单个波浪周期中全部采样点的瞬时波向;
图8为计数比例柱状图;
图9为本发明流程图。
具体实施方式
下面将结合附图对本发明作进一步描述。
如附图9所示,为本发明流程图;具体可以分为以下几步,第一步为野值修正,第二步为数据处理,第三步为计算单个波浪波向,第四步为统计单个波浪波向测量结果在16个角度区间的计数比例,第五步为计算主波向。
本发明实施例所用的波浪浮标为国产的SZF型波浪浮标,部署在平潭海域进行海上测波实验。工作方位为1小时定时测量方式,浮标在每天的整时进行波浪测量工作,采样间隔为0.5s,采样点为2048,浮标每1小时会连续采集1024秒的波浪数据,即一组实验数据。
上述SZF型波浪浮标的主要技术参数如表一所示:
表一 SZF型波浪浮标的技术参数
结合附图1~9,本发明具体实施步骤为:
第一步为野值修正,包括以下步骤:
步骤1.1、读取波浪浮标在海上进行海洋监测测得的实验数据,本例选取的一组连续实验数据观测时间为2010年10月6日00:59,附图1和附图2为浮标实测的原始纵倾角和横倾角数据,先利用筛选方法选出波浪浮标测量纵倾角和横倾角产生的野值,筛选方法如下:
遍历一组实验数据的所有纵倾角和横倾角数据值,确定超出-20°至20°范围的数据为野值,需要进行替补。
根据筛选方法,获得野值数据,野值数据的统计结果如表二所示:
表二 野值数据的统计结果
步骤1.2、经过步骤1.1筛选得到野值为
纵倾角的野值个数N
o为106,纵倾角的野值个数N
o为93,而野值数据不能简单的删去,否则数据值缺失仍然会降低波向测量的准确性,因此需要对其进行替补,替补方法如下:
针对一组数据,数据序列为X
1,...,X
2048,从第一个野值开始,针对任一野值X
oi(i=1,...,N
o),选择此野值前M个数据值X
oi-M,...,X
oi-1,其权重为
再选择此野值后M个数据值X
oi+1,...,X
oi+M,其权重为
若选择的X
oi+1,...,X
oi+M中有超出-20°至20°范围的数据值,则选择超出范围数据值后一个数据值,保证选取的数据值不超出-20°至20°的范围。
则该野值的替补值的计算公式为:
通常M=2,则该野值的替补值的计算公式为:
利用上述替补方法处理完纵倾角和横倾角数据中所有的野值,得到处理后的纵倾角和横倾角信息。附图3和附图4为处理后的纵倾角和横倾角数据。
第二步为数据处理,包括以下步骤:
步骤2.1、根据浮标测得的波高信息,采用上跨零点法,把波面上升与零点的交点,即瞬时波高数据由小于0变至大于0的点作为起点,当波浪运行到零点以下,即瞬时波高出现了小于0的情况后,将其后再次波浪上升与零点的相交点作为终点,则起点与终点之间为一个单个波浪周期。附图5为上跨零点划分波浪周期方法,图示中标注的圆圈所示的点为上跨零点,双箭头所示的区间为上跨零点法确定的单个波浪周期。
步骤2.2、利用上跨零点法,从第一个上跨零点开始,可将该组实验数据切分成138个单个波浪周期,即统计该组数据中单个波浪周期的个数为138。
第三步为计算单个波浪波向,包括以下步骤:
步骤3.1、根据步骤2切分好的单个波浪周期,在一个波浪周期内,选择波高数据大于0m的范围,作为参与计算单个波浪波向的采样点区间。附图6为某一个波浪周期的实测波高数据,图示中标注了参与计算单个波浪波向的采样点区间。
步骤3.2、对步骤3.1确定的采样点区间内的所有采样点,重复下列操作直至获得所有采样点的瞬时波向:利用瞬时方位角以及筛选和替补方法处理后的纵倾角、横倾角信息,先计算出采样点瞬时东向斜度和北向斜度,再计算出该采样点的瞬时波向。
式中:A为方位角,P为纵倾角,R为横倾角;
式中:A为方位角,P为纵倾角,R为横倾角;
采样点波向θp的计算公式为:
式中:
为浮标体采样点瞬时东向斜度,
为浮标体采样点瞬时北向斜度;
步骤3.3,对步骤3.2计算得到的采样点瞬时波向在空间域上进行分类,分类方法如下:
1、设置固定夹角为θc的扇形作为初试分类区,统计采样点瞬时波向值落在0~θc区间的个数U,通常θc取90°~180°;
2、分类区从0°开始,顺时针旋转并保持夹角仍为θc,每旋转1°统计一次该分类区中采样点瞬时波向的个数值U并记录类中的采样点瞬时波向值,旋转360°后,得U1,...,U360;
3、从U1,...,U360选出最大值UM,则此最大值所对应分类区中的采样点瞬时波向为分类后的采样点瞬时波向θp1,...,θpa(a=1,2,3...UM)。
本例中θc取90°,根据分类方法本例中得到某一采样点波向个数最大值UM为6,则分类后的各采样点瞬时波向为75°,89°,96°,105°,111°,123°。附图7为计算的某一单个波浪周期中全部采样点的瞬时波向,分类后的6个采样点波向落在60~150°之间。
步骤3.4、利用矢量叠加法处理步骤3.3获得的分类后的采样点瞬时波向,可得到瞬时波向北向分量ΔθN=-0.9901,瞬时波向东向分量ΔθE=5.6985,进而计算出单个波浪波向θw=100°。
瞬时波向北向分量ΔθN的计算公式为:
式中:θpi为分类后的第i个采样点瞬时波向;
瞬时波向东向分量ΔθE的计算公式为:
式中:θpi为分类后的第i个采样点瞬时波向;
单个波浪波向θw的计算公式为:
式中:ΔθE为瞬时波向东向分量,ΔθN为瞬时波向北向分量;
步骤3.5、重复步骤3.1、步骤3.2、步骤3.3和步骤3.4,直到得到该组数据中所有的单个波浪波向。
第四步为统计波浪波向在16个角度区间的计数比例。包括以下步骤:
步骤4.1、以22.5°为间隔将360°均匀分成16个角度区间,以地理正北方位为零度,按照顺时针方向来划分。角度区间划分标准如表三所示。
表三 角度区间划分标准
步骤4.2、统计步骤3得到的所有单个波浪波向θw1,...,θw138在16个角度区间内分别出现的个数为4,5,6,11,25,46,18,10,3,1,0,0,2,0,3,4。
步骤4.3、计算波浪波向在16个角度区间的计数比例,附图4为计数比例柱状图。
计数比例pi的计算公式为:
式中:ni为单个波浪波向在第i个角度区间出现的个数,N为波浪周期的总个数,本例中为138;
第五步为计算主波向。包括以下步骤:
步骤5.1、根据步骤4得到的计数比例,得到计数比例最大值K=33.3%;
步骤5.2、计数比例最大值所对应的角度区间为第6角度区间,该角度区间的中值为112.5°,则主波向为112.5°;
在本海区的长期观测结果表明,利用wavex雷达测波仪测量的波峰峰向与具有多年观测经验的观测员人工观测结果一致性非常好,因此利用wavex实测的波峰峰向数据和本算法计算的浮标自测结果进行了比测,以检验算法性能。
1、在该时刻,wavex实测结果为105°。
2、传统的浮标测量波向方法如下:首先根据上跨零点来确定单个波浪周期,再确定每一个单个波浪周期的下跨零点,利用浮标测得的这一时刻瞬时方位角、纵倾角和横倾角信息,计算得到该单个波浪周期的参考波向;然后把各个波浪周期测得的单个波浪波向进行归并统计,将360°分成16个角度区间,统计16个角度区间中每个区间内单个波浪波向出现的个数,进而得到波向出现率,最后,选择波向出现率最大的角度区间的中值作为主波向。该算法计算的结果为67.5°。
误差分析表明,本发明计算结果与wavex实测结果的绝对误差E1=7.5°,浮标传统算法测量结果与wavex实测结果的绝对误差E2=37.5°。本发明计算结果与wavex实测结果基本吻合,表明本改进方法有效,且优于传统算法。
本发明计算结果与wavex实测结果的绝对误差E1计算公式为:
E1=|θ -θwavex|
式中:θ为本发明计算结果,θwavex为wavex实测结果;
浮标传统算法测量结果与wavex实测结果的绝对误差E2计算公式为:
E2=|θbuoy -θwavex|
式中:θbuoy为浮标传统算法测量结果,θwavex为wavex实测结果;
为避免偶然性,选取不同时段的多组实验数据,按照上述步骤进行实验验证。本用例选取了40组不同时段的波浪浮标实验数据进行实验验证,本发明计算结果、浮标传统算法测量结果和wavex实测结果对比如表四所示。
表四 本发明计算结果,浮标传统算法测量结果和wavex实测结果对比
实验结果表明,与wavex实测结果对比,本发明计算结果的精度为12°,浮标传统算法计算结果的精度为41°。
本发明所提出的基于纵横倾角测量野值修正的波浪浮标主波向计算方法,不仅在实际测量中有效,与Wavex实测结果吻合度更高,且在精度上更加优于浮标传统算法测量的主波向结果,提高了波浪浮标测量主波向的稳定性和准确率,可在海洋观测设备中进行广泛的推广和应用。