姜晓鹏,高志强,吴晓青,王跃启,宁吉才
(1.中国科学院海岸带环境过程与生态修复重点实验室,山东烟台264003;2.中国科学院大学,北京100049;3.中国科学院烟台海岸带研究所山东省海岸带环境过程重点实验室,山东烟台264003)
自2007年夏季,黄海绿潮已连续暴发14年,成为黄海海域一种常态化的海洋生态灾害[1–2],对青岛、海阳、日照等沿海城市的滨海旅游、海水养殖、生态环境等造成了严重的危害[3]。已有研究确认了黄海绿潮的原因种为浒苔(Ulva prolifera)[4–5],通过野外调查、遥感以及数值模拟等手段证明了黄海绿潮最初起源于苏北浅滩海域[6–7],然后向北进入黄海,边漂移边生长,并最终在山东半岛东部沿海聚集堆积[8–10]。黄海绿潮具有长时间跨区域输运、异地成灾的特点,准确掌握其漂移速率、方向与轨迹等动态信息,能为绿潮灾害预警与防控提供重要参考。
以往绿潮漂移速度的研究数据多源于卫星遥感。如衣立等[11]利用2009年6月24−25日的MODIS的绿潮影像,提取了黄海3个区域的绿潮漂移速率,分别为0.2 m/s、0.4 m/s和0.12 m/s;夏深圳[12]基于MODIS与GOCI数据,估算了2013−2015年黄海绿潮漂移速率的变化范围,为0.01 ~0.98 m/s,并指出绿潮漂移速度存在明显的时变特征。陈晓英等[13]利用GF-4卫星影像,基于最大相关系数法追踪绿潮,分析了9时与15时绿潮漂移速度的不同变化。从MODIS的1 d,到GOCI的1 h,再到GF-4卫星的20 s,虽然所用影像的时间分辨率不断提高,但卫星遥感极易受到天气及云层的干扰,如2013−2015年5−7月的MODIS数据中,可用于绿潮速度研究的影像只占20.7%[12]。GF-4卫星在南黄海区域的重访周期并不恒定,不能满足绿潮灾害应急监测的要求,且MODIS与GOCI的空间分辨率分别为250m与500m,低分辨率的影像中存在着大量的混合像元[14−15],从而影响绿潮提取及速度估算的精度。
无人机遥感具有时空分辨率高、不受云层干扰等优势,能根据任务需求灵活开展航测作业,在陆地灾害评估、测绘、作物估产等方面发挥了重要的作用,是地面调查与卫星遥感的重要补充[16]。近年来,随着小型旋翼无人机技术的发展,无人机应用也逐渐扩展到海洋环境监测[17]、绿潮初始生物量估算[18]等海洋研究领域。无人机与科考船的联合作业,可发挥船载设备能同时获取海上气象、水文等数据的优势。
本文利用科考船作为无人机起降平台,首次获取了双时相的绿潮正射影像,基于可见光波段的漂浮藻类指数(RGB-FAI)[19]提取了无人机可见光影像中的绿潮,计算出高精度的绿潮漂移速度,并与GOCI影像反演的速度结果进行对比。同时结合船载气象仪的风速数据与已有研究的潮流数据,分析了绿潮漂移对风与潮流的响应。本文以期扩展无人机在绿潮监测方面的应用,为绿潮灾害的应急预警与防灾减灾提供高精度的技术支撑。
研究区位于南黄海海域(图1)。南黄海为山东半岛成山角与朝鲜半岛长山串连线南部的黄海区域,面积约为30.9×104km2,平均水深约46m,其西部沿岸海域为规则半日潮。南黄海受东亚季风影响,秋、冬季盛行北风或西北风,平均风速达6~8m/s,海流多为偏南向流。春季,风向多变,风力减弱;夏季,东亚夏季风与副热带高压带来平均风速为4~6m/s的偏南风,海流为偏北向流,海表面温度达到一年中的最大值[13,20]。2007年后,南黄海在每年的春夏季都会出现大规模的绿潮,其通常开始于3月末至4月的苏北浅滩海域,约5月中下旬进入暴发期,自7月则逐渐衰退,直至最终消亡[21]。
利用大疆“悟1”无人机获取3个站点(图1)的海上绿潮影像。该四旋翼无人机体积小,为垂直起降方式,可在小面积场地起降,其惯性测量单元模块以及最大5级的抗风能力,能保证无人机在海上复杂气象条件下的飞行姿态稳定。该无人机采用GLONASS+GPS双系统定位,搭载4K可见光高清镜头,获取的航拍图像具有坐标信息,在小范围的大比例尺测图中已得到了应用与精度验证[22]。无人机对绿潮进行航测作业时,以科考船甲板为起降平台,飞行高度尽可能使单幅图片能容纳所拍摄的绿潮斑块宽度,并保证图像的航向重叠度大于60%,以满足后期图像拼接的要求。对漂移的同一绿潮斑块分别执行2次航拍(图2),生成其在2个时相的数字正射影像。
图1 航拍站点、验潮站及研究区示意图Fig.1 Map of aerialphotography sites,tide station and study area
图2 船载无人机获取绿潮正射影像示意图Fig.2 Green-tidalorthophoto taken by ship-borne unmanned aerialvehicle (UAV)
采用RGB-FAI与阈值分割法提取无人机正射影像中的绿潮。RGB-FAI是一种基于基线校正思想,用于从可见光影像中提取绿潮的植被指数,其提取精度高,并能适用于薄雾气象条件下对零散、小面积绿潮斑块的提取与监测[19]。RGB-FAI的算术公式如下:
式中,R表示反射率或像元值;λ表示中心波长,下标red、green、blue分别表示红光波段、绿光波段、蓝光波段,其中,λred=700 nm, λgreen= 546.1 nm,λb=lu4e35.8 nm[23]。
由于绿潮不断漂移,为减小速度估算的误差,在整个绿潮斑块选取一定数量的子斑块,通过计算子斑块速度的平均值来表示该站点绿潮的漂移速度。首先,子斑块的选取原则为:(1)尽可能均匀分布在整个绿潮范围内;(2)子斑块形状清晰、明显,容易识别;(3)形态特征在两幅影像中基本不变。其次,为保证两个时相影像中相同子斑块的代表点位一致,利用A rcmap软件的Mean Center工具生成这些子斑块的地理中心点。然后,量测各地理中心点(子斑块)的位移,通过子斑块所在图像的拍摄时间获得其漂移时长,最终可计算出各个子斑块的漂移速率与方向,并求出它们的平均值。整个绿潮漂移速率Vave为这些子斑块漂移速率的算术平均值,而漂移方向Dirave取子斑块的矢量平均值,两者的计算公式如下:
式中,i与n分别代表子斑块及其数量;S、t和Dir分别表示漂移的距离、时间和方向。
GOCI(Geostationary Ocean Color Imager)影像从韩国海洋卫星中心网站(http://kosc.kiost.ac.kr/)下载,影像包含从412~865 nm的可见光及近红外光共8个波段,每天可提供从北京时间08:30−15:30的8景逐小时影像,其时间分辨率(1 h)和空间分辨率(500m)优势可用于黄、渤海及东海海域的海洋环境监测[24]。
GOCI影像首先利用GDPS(GOCIData Processing System)进行辐射校正,然后基于GLT(Geographic Lookup Table)进行几何校正,选取并裁切出2幅影像中同时未被云层遮挡的海域。参考已有的绿潮提取算法研究,采用NDVI指数[25−26],并基于Jenks的自然断裂点分级法和目视解译,提取海水中的绿潮。其中,NDVI=(RNIR−RRED)/(RNIR+RRED),RNIR、RRED分别对应于GOCI影像的第7近红外波段(波长为745 nm)与第5红光波段(波长为660 nm)。
2.5.1 数据来源
本文的风速数据分为航次站点的风速数据与南黄海区域的风场数据。航次站点的风速数据来源于中国科学院海洋研究所海洋大数据中心(http://msdc.qdio.ac.cn/)在2019年6月的透明海洋共享航次,由“科学三号”的“XZC6-1型船用自动气象仪”所获取的1min平均真风速、真风向数据。海域的风场数据来源于欧洲中期天气预报中心(ECWMF)发布的再分析资料(ERA 5),从中提取出南黄海海域在6月19日当天1 000 hPa大气压下,海面10m高处的逐小时风场数据,空间分辨率为0.25 °×0.25 °。
潮流数据引用了已发表的成果,为南黄海6月表层的涨落潮流数据。此外,潮汐表用于查询站点作业时的潮位,其来源于中国海事服务网(https://www.cnss.com.cn/tide/),考虑到站点S1、S2所在的纬度(36°N)与青岛港验潮站(36°04′N)相近,站点S3的纬度(35°N)与岚山港验潮站(35°05′N)相近(图1),因此分别选取了青岛港2019年6月19日与岚山港2019年6月16日的潮汐表的潮位数据。
2.5.2 站点有效风速数据的选取与处理
由船载气象仪记录的风速数据可知,风速在作业时段变化较大。若航拍所记录的绿潮漂移时间段为T0~T1,则作用于该时段绿潮漂移的有效风速时段,并不等同于船载气象仪风速资料中的T0~T1。这是由于以下两点:第一,各个站点中,科考船与航拍绿潮的距离有数百米至1 km,风通过这段距离,需要一定时间(设为△t,单位:m in);第二,绿潮漂移具有惯性,对风等外力驱动作用有一定的滞后时间,本文将该响应时间设为2m in。通过量测各站点船载气象仪与绿潮中心在风向上的投影位移,结合风速,可估算出△t。当科考船位于绿潮的上风向时,将航拍记录的时段整体前移△t,选取的有效风速时段为T0−△t−2~T1−△t;下 风 向 时,则 整 体 后 移△t,为T0+△t−2~T1+△t。则站点S1、S2和S3的△t的值分别为−1m in、5m in、5m in。各站点选取的有效风时段如图3 中的红色竖线间所示。
图3 3个站点作业期间的1m in平均风速矢量Fig.3 The 1-minuteaverage wind data of threesites during operation
最后,用选取时间段内的1 min平均真风速和风向数据来计算风速的算术平均值V风速和风向的矢量平均值Dir风向,表示作用于绿潮漂移的有效风。两者的计算公式如下:
式中,Vi代表第i时刻的1m in平均真风速;Diri代表第i时刻的1min平均真风向。
用RGB-FAI分别提取3个站点中的绿潮。为检验该指数提取绿潮的精度,从S1站点的绿潮分布图中,截取宽为317 m、长为662 m的矩形区作为样例(图4 a),计算其RGB-FAI(图4 b),并设定阈值提取出绿潮,同时在图像区域内生成999个随机验证点(图4 c),利用无人机影像的高空间分辨率(0.12 m)、易目视识别的特点,对该矩形区域的绿潮提取结果进行精度评价,基于误差矩阵表获得生产者精度与用户精度(表1),并计算提取结果的总精度与kappa系数。
表1 绿潮提取误差矩阵表Table 1 Error matrix of green-tide extraction
图4 绿潮提取步骤Fig.4 Steps of extracting green-tide
计算出绿潮提取的总精度为97.9%,kappa系数为0.95。kappa系数为一致性检验指标,值域为−1~1,当系数在0.81 ~1时,认为分类结果与真实状况几乎完全一致[27],因此,基于以上验证结果,认为RGBFAI能够适用于无人机可见光影像的海上绿潮提取。
基于RGB-FAI,对3个站点的双时相无人机正射影像进行绿潮提取,并在3个站点各选出6~8个子斑块,量测出子斑块地理中心点的位移、时间、方向(图5),计算出各子斑块的漂移速度。利用公式(2)、公式(3),依次计算出各站点的绿潮漂移速率和方向,并按照2.5 节中确定的有效风速时间段,计算各站点的平均风速与平均风向,结果如表2 所示。
表2 绿潮斑块漂移速度与风速的统计Tab le 2 Statistics of green-tide drifting velocity and the corresponding wind velocity
图5 3个站点的绿潮斑块漂移图Fig.5 Thegreen-tide drifting map of three sites
在站点S1,绿潮漂移速率为0.257 m/s,其北向速度分量较大,符合绿潮北漂的整体趋势。站点S2,绿潮漂移速率达到0.439 m/s,其东向速度分量较大。在站点S3,绿潮漂移速率为0.256 m/s,方向为东南方向,有南向的速度分量。结果表明,绿潮在从苏北浅滩海域到山东近岸的北向漂移过程中,甚至会有暂时向南漂移的时段(S3);在一天中的不同时刻,邻近海域的绿潮漂移速率和方向也会明显不同(S1、S2),无人机遥感可以精确地揭示出绿潮的瞬时漂移状态。
一般而言,海面漂浮物在长时段(大于12 h)的漂移轨迹受风与局部洋流的共同影响,但短时间内的漂移则会明显受潮流的驱动[28]。本文计算绿潮速度所记录的漂移时长为2~29m in,远小于潮汐涨落1次所需的约12 h,因此,绿潮的驱动需要结合潮汐流进行分析。此外,风是绿潮漂移的重要驱动因素。本文利用船载气象仪记录的风速数据、验潮站记录的潮汐表以及南黄海6月表层的涨落潮流场,将绿潮漂移与同一时刻的风速(表2)、潮流进行一致性分析。
南黄海区域涨潮与落潮时的潮流方向并非简单的朝向正西、正东。于晓杰[29]基于FVCOM(Finite Volume Coastand Ocean Model)海洋数值模式计算出了南黄海6月M2分潮的涨落潮流场,即1个潮周期内4个主要时刻(高潮时、落潮中间时、低潮时、涨潮中间时)的表层流场分布。本文的航拍作业同为6月的南黄海海域,查询S1、S2、S3站点同纬度近岸验潮站的潮汐表,根据各站点的作业时间查询出航拍时段的潮汐状态(图6)。据此,在南黄海6月表层“低潮时流场分布图”中可得到S1站点的潮流方向(图7 b),在“落潮中间时流场分布图”中可得到S2、S3站点航拍时的潮流方向(图7 d)。
图6 3个站点航拍作业时的潮位Fig.6 Tidal height during aerial photography of three sites
在站点S1,绿潮漂移方向位于风向右侧约57°,与潮流同为东北向(图7 a,图7 b),风在绿潮漂移方向上的速率分量为1.7 m/s,约为绿潮速率的6.7 倍。在站点S2,绿潮方向偏于风向右侧约62°,与潮流方向基本一致,夹角约为20°,为偏东向(图7 c,图7 d),风在绿潮漂移方向上的速率分量约为1.5 m/s,为绿潮速率的3.5 倍。在站点S3,绿潮漂移方向在风向右侧约1°,与潮流方向夹角约为24°,同为东南向(图7 d,图7 e),风在绿潮漂移方向上的速率分量为1.7 m/s,约为绿潮速率的6.8 倍。可以看出,3个站点的绿潮漂移与潮流方向基本一致,位于风向右侧1°~62°。风在绿潮漂移方向上的速率分量为绿潮速率的3.5 ~6.8 倍,非恒定系数。
图7 绿潮、风、潮流的移动矢量图Fig.7 Motion vector of green-tide,wind and tidal current
GOCI影像的时间分辨率为1 h,同样可应用于绿潮漂移速度的估算,因此,可与同时段下无人机遥感的速度做对比分析。
S1、S3站点航拍时段的海面绿潮,在GOCI影像上被云层覆盖,本文选取了与S2站点同一时段,即19日09:30与10:30的两景GOCI影像反演绿潮速度。由于S2站点被云遮挡,在另两处无云海域提取绿潮。利用Arcmap软件计算两个区域绿潮在两景影像的中心点,量测出中心点在09:30−10:30的移动距离和方向,分别算出两个区域绿潮的漂移速度(图8),最后取两者的均值作为GOCI影像反演的绿潮漂移速度。经计算,绿潮漂移速率均值为0.547 m/s,方向均值为北偏东58.4°。
图8 基于GOCI影像的绿潮逐时漂移图Fig.8 Theone-hour drifting map of green-tide from GOCI images
与S2站点的无人机遥感结果相比(表3),两种方式反演的绿潮漂移方向基本一致,同为东北向,漂移速率上相差0.108 m/s。需要说明的是,由当天10时的风场数据(图8 a)可知,S2站点与GOCI反演海域在风速与风向的相对偏差分别为0.58 m/s和3%,因此两个海域绿潮漂移速度的实际值可能也存在差异。由于无人机影像的时空分辨率远高于卫星影像,若能获取同时段、同位置绿潮漂移的两种影像,则无人机可以为卫星反演大范围海域的绿潮漂移速度提供高精度校准。
表3 基于UAV和GOCI影像的绿潮漂移速度提取结果的对比Table 3 Com parison of green-tide drifting velocity derived from UAV and GOCI images
本文通过无人机遥感与船载气象仪获取的绿潮漂移速度与风速数据具有较高精度,但在以下方面存在不确定性。首先,对于绿潮漂移的驱动分析,用于分析的数据样本只有3个,对绿潮漂移的定量分析还需要更多样本数据的支持。第二,在绿潮对风驱动作用的滞后时间上,根据主观经验设置为2 min,而这需要相关专业后期理论与实验的探索。第三,本文引用了基于海洋数值模式的潮流数据,能用于绿潮漂移方向的驱动分析,未来研究可考虑利用船载的声学多普勒流速剖面仪(ADCP),获取更精准的海水表层流场数据,为绿潮漂移的定量分析提供更精准的数据支撑。
受风、浪、流的影响,海上绿潮斑块的形态在超过6 h时变化较大[13],因此绿潮漂移速度的研究需要高时间分辨率的遥感影像。常用的卫星影像中MODIS、GOCI与GF-4能满足时间分辨率的要求,但MODIS和GOCI影像的空间分辨率较低,GF-4卫星的成像范围并不固定在绿潮发生区,且卫星影像易受云层影响。
而基于船载起降的无人机,可依靠科考船实时追踪绿潮,具备作业灵活、获取影像的时空分辨率高、不受云层干扰等优势,可快速、连续、更为精细地监测绿潮漂移。这些优势,一是可用于对重点区域的绿潮应急监测与灾害预警,如海滨景区、重要设施(核电站进水口等)的附近海域;其次,可校正卫星遥感提取结果,实现更大范围海域上绿潮漂移的精确监测。最后,可为数值模型预测绿潮漂移等研究提供高精度的验证数据与方法。同时,科考船搭载的温盐深仪、ADCP、自动气象仪等设备,可以测量海水的温度、盐度、叶绿素浓度、海流速度、流向以及大气温度、风、气压等诸多水文与气象参数,这些数据与船载无人机获取的绿潮面积、形态、漂移速度等高精度结果相结合,可以为绿潮生长、漂移相关的研究提供更加充分的数据支持。
本文基于船载无人机获取了黄海绿潮双时相的数字正射影像,开展了绿潮漂移速度的估算,并与GOCI卫星影像反演的绿潮速度进行对比,同时分析了风、潮流对绿潮漂移的驱动。结论如下:
(1)RGB-FAI可有效提取可见光影像中的绿潮信息(kappa系数为0.95),可辅助无人机高精度的估算绿潮漂移速度;3个站点的绿潮漂移速率为0.26 ~0.44 m/s,漂移方向在一天中会发生明显变化,甚至有南向分量的漂移。
(2)绿潮短时间内的漂移受风与潮流的共同影响。绿潮漂移方向与M2分潮的潮流方向基本一致,位于风向的右侧1°~62°。
(3)船载无人机对绿潮漂移状态的监测,可校正卫星遥感的结果,能为绿潮灾害的应急监测、预警及相关研究提供重要的数据与技术支撑。
致谢:感谢中国科学院海洋研究所海洋大数据中心(http://msdc.qdio.ac.cn/)提供数据支持。