做沿海洪水淹没评估、海平面上升情景推演、风暴潮风险制图这类工作第一道坎几乎都卡在高程数据上。很多新手一开始去下载普通陆地DEM拿到沿海区域一画就露馅陆地部分看着还行一到滩涂、河口、潮间带要么是空的要么高程乱跳。直到我接触了 CoNED——沿海国家高程数据库Coastal National Elevation Database才算真正把“陆地高程”和“近岸水深”这半张地图拼完整。这篇文章就聊聊我在 GEEGoogle Earth Engine里实际使用 CoNED 做沿海分析的经验、脚本以及那些不会写在官方文档里的坑。1. 为什么沿海分析一上来就卡在高程数据上1.1 陆地有高程浅海没深度普通DEM的“半边图”问题普通DEM产品通常只覆盖陆地有些全球影像产品即便画到了海岸线近岸浅海部分要么是空洞要么是粗插值结果。但沿海洪水淹没分析的物理过程恰恰发生在海陆交界风暴潮把水面从海里抬起来漫过滩涂沿着低洼地形向陆地推进。这个过程需要知道水面下的坡度和陆地上升的坡度缺了任何一边模拟出的淹没边界都是“半拉子工程”。CoNED这类产品做的事情简单说就是“海陆一张图”把来自航空激光雷达、船载测深、多波束声呐、卫星测高等多渠道的地形数据统一对齐到一个连续的栅格上。所以你在GEE里加载它时能同时看到沙滩、盐沼、潮沟、水下河道和近岸浅滩而不是一条生硬的“陆止于此、海始于斯”切割线。1.2 同一个高程数值可能是“两把不同的尺子”量的沿海DEM对垂直基准特别敏感这里需要多解释几句。陆地高程通常基于某个高程基准面比如常用的EGM2008大地水准面模型来计算“海拔”而海图水深通常基于当地理论最低潮面来计算“深度”。两者本身就不是同一把尺子。如果你直接拿一个海图水深栅格和陆地DEM拼接交界处可能出现几十厘米甚至一两米的高程跳变。在我处理过的一个沿海项目中两套数据在潮间带接缝处跳变了约0.35m这个误差放在内陆基本无感放在沿海淹没模拟里直接导致某条设计堤防的“是否漫顶”结论反转。所以判断一份沿海高程数据是否可用第一件事不是看分辨率而是看它的垂直基准面说明。1.3 什么场景下非用CoNED不可需要这类海陆一体高程数据的场景很明确海岸带洪水风险评估、海平面上升情景制图、风暴潮增水模拟、近岸生态栖息地规划、港口和跨海通道选址预研、海底电缆路由初选等等。这类任务大多要求对“水淹到哪里、淹多深、哪些资产暴露”做出定量估计。普通全球DEM在陆地部分足够应付宏观尺度但到了近岸工程尺度比如一段几公里长的海湾堤防评估1m、3m、10m分辨率的高精度CoNED产品就是刚需。GEE的好处在于你不用把十几GB的栅格下载到本地直接在云端处理还能批量叠加人口、土地利用、路网等数据做统计分析。2. 在GEE里找到CoNED数据并把它变成可计算图层2.1 别靠记忆填Asset ID要在数据目录里确认很多教程会直接甩一个数据集ID让你加载。实际经验是GEE的数据目录经常调整同一条数据在不同时期可能被分到不同子目录Asset ID也有变化。第一次上手建议打开GEE代码编辑器右侧的Search面板输入“CoNED”或“coastal elevation”先把可用结果刷一遍再决定用哪一份。我看一份CoNED候选影像时会重点检查四件事分辨率常见的有1m、5m、10m分辨率越高覆盖范围越小。垂直基准面说明数据说明里会写清楚是否已经做过基准面统一。波段构成主波段通常是elevation附带的hillshade等波段用于可视化验证。覆盖范围用Map.addLayer先叠加到地图上看一眼确认确实覆盖目标海岸段。加载示例代码如下ID仅作演示请按你在目录中查到的实际ID填写// 加载一份沿海高程影像示意ID实际以目录查询为准 var dem ee.Image(A/CoNED/coastal_dem_v1) .select(elevation) .rename(dem); // 添加可视化图层 Map.addLayer(dem, {min: -5, max: 10, palette: [#0d0887, #46039f, #f0f921]}, CoNED elevation);这里用负到正的范围是因为沿海DEM同时包含浅海负值区和陆地正值区色带两端要留够空间否则近岸地形细节会被压平。2.2 先打印元数据别急着做计算加载之后我几乎一定会做一步“傻操作”print(dem)把影像的波段、投影、分辨率、时间戳全部打印到Console里看一眼。别觉得多余这一步能拦住大部分后面要返工的问题。比如我就遇到过一份影像的描述字段写着“海陆无缝衔接”实际打印后发现它的近岸水体像元压根是0值而不是NoData导致后续淹没模拟在水域部分算出一片“假陆地”。这类信息只有看原始元数据才清楚靠目视抽查很容易漏。// 打印影像元数据 print(CoNED metadata:, dem); // 查看投影与分辨率 var proj dem.select(dem).projection(); print(投影信息:, proj); var scale dem.select(dem).projection().nominalScale(); print(名义分辨率(米):, scale);2.3 切分成多幅时的合并策略有些区域的CoNED数据不是单张影像而是按图幅切分的ImageCollection。这时不要手动一张张加载再拼接直接在集合上做mosaic即可var collection ee.ImageCollection(A/CoNED/tiles_v1) .filterBounds(roi); var mosaic collection.mosaic() .rename(dem) .float(); Map.addLayer(mosaic, {min: -2, max: 8, palette: [#0d0887, #f0f921, #fca636]}, CoNED mosaic);注意mosaic对重叠区域的默认处理是“后加载的覆盖先加载的”正常情况没问题但若发现重叠带出现阶梯状差异就要回头检查各图幅是不是使用了不同基准面或不同年份的测深数据。3. 跑通第一个淹没模拟常数潮位到多情景对比3.1 静态淹没的本质是一场“比高”计算有了CoNED高程第一个能立刻落地的应用是静态淹没模拟。它的原理不复杂给定一个水面高度凡是DEM高程低于这个水面的陆地区域就是潜在淹没区再用水面高度减去地面高程就是淹没水深。为什么不直接说“低于海平面的地方就是淹没区”因为沿海很多防御设施和人工填海地块本身就低于平均海面却被堤防保护着。静态淹没模拟反映的是“头顶没有堤防”时的最坏情形常被用来做警示性评估而不是精确到具体某一次风暴事件的实际内涝。理解这个前提后面的结果才不会用偏。3.2 基础版本用一个常数潮位算淹没范围最朴素的写法是构造一个常数水面与DEM做减法// 设定一个水位值单位米相对于DEM的垂直基准 var tideLevel 1.8; // 计算水深正数表示被水淹没 var depth ee.Image.constant(tideLevel) .subtract(dem) .rename(depth); // 淹没掩膜水深大于0且只保留近岸陆地区域 var floodMask depth.gt(0) .and(dem.gt(-1)); // 避免把深水区也算成“淹没陆地” var flood depth.updateMask(floodMask); Map.addLayer(flood, {min: 0, max: 3, palette: [#1d91c0, #7fcdbb, #ffffcc]}, Flood depth 1.8m);这里dem.gt(-1)这条很容易被忽略作用却很大。如果不加这个限制一片水深50m的近岸海域也会被划进“淹没区”面积统计会完全失真。我习惯把下限设为潮间带范围具体值根据当地潮差调整。3.3 把模拟包装成函数一键生成多情景对比沿海淹没从来不是只算一个水位就结束的。通常要对比平均高潮位、极端风暴潮位、海平面上升极端潮位等不同情景。写成函数最省事function simulateFlood(level, name) { var depth ee.Image.constant(level) .subtract(dem) .max(0) .updateMask(ee.Image.constant(level).subtract(dem).gt(0) .and(dem.gt(-1))); return depth.set(scenario, name).rename(depth_ name); } var s1 simulateFlood(1.2, avg_high); var s2 simulateFlood(2.5, storm_surge); var s3 simulateFlood(3.4, slr_storm); Map.addLayer(s1, {palette: [#c6dbef]}, 平均高潮); Map.addLayer(s2, {palette: [#6baed6]}, 风暴潮); Map.addLayer(s3, {palette: [#08519c]}, SLR风暴潮);多情景并排看的时候最直观的产物就是“边界不重叠度”哪些区域在1.2m时还是干的到3.4m时成了重灾区这些区域往往是基础设施规划最该回避的高风险地带。3.4 用面积统计量化“淹了多少”只出图不够领导汇报要数字。GEE里统计淹没面积的标准做法是乘像元面积var areaImage flood.multiply(ee.Image.pixelArea()); var stats areaImage.reduceRegion({ reducer: ee.Reducer.sum(), geometry: roi, scale: 10, maxPixels: 1e9 }); print(淹没面积(平方米):, stats.get(depth));这里必须用ee.Image.pixelArea()而不是自己用分辨率估算因为影像在不同纬度、不同投影下像元实际地面尺寸并不完全一致。4. 剖面线与坡度验证数据质量也预判淹没边界4.1 为什么我强烈建议先画一条横穿海岸的剖面线加载一份CoNED数据后我做的第一件“验证性工作”永远是画一条横穿海岸的剖面线而不是直接开算。原因很简单任何淹没模拟都建立在“高程数字可信”的前提下而数字是否可信一条线就能暴露七七八八。剖面线还能帮你理解地貌控制机制。比如某段海岸如果前50m坡度很缓后面开始急剧抬升那2m风暴潮位可能淹不到太远反之如果是一条平缓潮沟延伸到内陆那淹没就会顺着沟往里钻。这些空间特征光看覆盖范围图是看不出门道的。4.2 剖面提取脚本一条线看穿高程数据在GEE里用ui.Chart.image.profile可以直接对一条线采样并绘出剖面图// 创建一条横穿海岸的示意线请替换为你的实际坐标 var line ee.FeatureCollection( ee.Feature(ee.Geometry.LineString([[-75.5, 38.2], [-75.3, 38.0]]))); // 绘制高程剖面 var profile ui.Chart.image.profile({ image: dem, region: line, scale: 10, reducer: ee.Reducer.first() }); print(海岸剖面图:, profile);这条命令会输出一张横轴为采样距离、纵轴为高程的折线图。拿到图以后重点看三处海陆交界处有没有明显台阶可能提示基准面冲突。滩涂部分是否出现锯齿状跳动可能提示测深与激光雷达拼接异常。陆地近岸段有没有条带状洼地这些位置往往是淹没模拟里的“咽喉”。4.3 坡度数据找出淹没边界变化最剧烈的地方淹没边界对高程误差的敏感程度和地形坡度直接相关。在坡度很缓的平原高程差0.5m淹没边界可能横向移动几百米在陡峭的海崖差1m也几乎不动。因此我会配合剖面线把坡度渲染出来看var slope ee.Terrain.slope(dem); Map.addLayer(slope, {min: 0, max: 5, palette: [#ffffff, #fee08b, #d53e4f]}, Slope 0-5deg);坡度图层特别适合找“隐性地势低点”那些绝对高程并不低、但因为周围全是隆起而被围住的小洼地。风暴潮一旦漫过外围的垭口这些洼地就会积水且退不出去属于容易被常规淹没图忽略的次生风险区。5. 把淹没结果转成暴露评估人、路、用地叠加5.1 人口暴露量统计重灾区覆盖多少人有了淹没水深图层接下来通常要问“淹到了谁”GEE里最常见的做法是叠加人口栅格数据统计不同淹没等级下受影响人口数量。操作上我一般先把水深分等级例如0-0.5m、0.5-1.5m、1.5m以上再按等级分别统计人口// 按水深分级 var class1 flood.gt(0).and(flood.lte(0.5)).selfMask(); var class2 flood.gt(0.5).and(flood.lte(1.5)).selfMask(); var class3 flood.gt(1.5).selfMask(); // 假设这里有一个全局人口栅格影像请替换为实际Asset var population ee.Image(A/GlobalPop/v1).rename(pop); function sumPop(mask, name) { var maskedPop population.mask(mask); var sum maskedPop.reduceRegion({ reducer: ee.Reducer.sum(), geometry: roi, scale: 100, maxPixels: 1e9 }); return ui.Label(name 受影响人口: sum.get(pop)); } print(sumPop(class1, 0-0.5m)); print(sumPop(class2, 0.5-1.5m)); print(sumPop(class3, 1.5m));这里要注意人口栅格的分辨率往往比CoNED低很多统计结果只适合作为宏观排序不适合直接套进工程级风险计算。5.2 基础设施暴露道路与建设用地“一选就中”基础设施暴露评估的思路完全类似把路网、管线、学校、医院等矢量数据放在淹没层下做一个空间连接。比如查某段海岸在极端情景下有多少米的道路被淹可以直接用reduceRegions// 构造一个示意路网要素集实际请加载你的道路数据 var roads ee.FeatureCollection(A/Roads/Example); var roadsExposure roads.map(function(f) { var floodedLength flood.gt(0) .reduceRegions({ collection: f, reducer: ee.Reducer.mean(), scale: 10 }); return floodedLength; }); print(道路淹没栅格均值:, roadsExposure);这个统计输出的是每个道路要素在淹没层上的平均水深如果只想挑选“完全断通”的路段可以再加一个“最大水深是否超过阈值”的过滤条件。阈值通常取0.3m——多数小型车辆在这之上已经无法安全通行。5.3 暴露结果如何回馈到高程分析暴露评估做完经常会出现一道“回头题”某个关键路段或规划地块恰好落在淹没边界附近而边界处的高程数据又只有1-2m的精度这时就需要回到CoNED数据本身去查证边界位置。这也是我强调“在工作流里留一条剖面图出口”的原因。遇到关键决策点画一条精细剖面亲眼确认该位置的微地形比盲目相信一次全局模拟更稳妥。6. 实际项目里踩过的几个坑和我的处理SOP6.1 坑一水面像元的NaN和0值处理不当CoNED里的水域像元在不同区域版本里表现完全不同有的写成NoData有的写成0有的写成负值水深。如果把0值当成真实地形水位0.5m时所有0值的小水体都会变成“被淹区域”图上全是斑点。我的处理习惯是加载后先用直方图检查像元值分布如果看到0值处有明显峰值再对照海岸线矢量确认这些0值是否位于水体内部。确认后在模拟脚本里通过掩膜把所有非目标区域剔除而不是直接选算了事。6.2 坑二基准面差异会隔着一个区域“发作”不同分幅的CoNED数据可能采用不同的垂直基准。特别是跨省级、跨机构成果合并时某一段海岸整体偏高或偏低的现象并不罕见。遇到这种情况我会先做一个“接缝诊断”分别渲染相邻两幅影像设置相同的色带肉眼观察接缝两侧是否有系统性色差。如果色差存在再分别取两侧的统计值计算偏移量然后先偏移校准再合并。这个偏移量怎么来最靠谱的是找沿岸已有的潮位站高程点或精密GNSS水准点做控制而不是凭空猜。6.3 坑三分幅拼接后出现“接缝辫子”即使基准面一致不同图幅的分辨率和获取年份也可能不同。mosaic之后接缝处会出现星星点点的异常——一条辫子状的纹理。这种区域如果在淹没边界附近会导致模拟结果出现沿着接缝的虚假“狭长淹没带”。检测方法很简单做一次坡向提取接缝处往往会出现高度规则的长条状异常高坡。筛掉它们的方法是在后续计算里用邻域统计比如均值滤波做一次平滑或者干脆在关键研究区内只使用单一图幅避免跨幅拼接。6.4 坑四数据版本更新导致结果无法复现沿海高程数据会不定期更新同一条海岸去年和今年的版本可能差不少。项目评审阶段如果要求复现半年前的模拟结果而Asset已经更新数字就会对不上。所以我现在的习惯是所有重要模拟都会导出一份中间GeoTIFF存档同时在项目说明里写明数据Asset ID、加载日期和模拟脚本版本。将来有人问“这个数怎么来的”我能直接拿档案出来对质而不是再跑一遍某个可能已经变了的新数据。6.5 我现在用于沿海高程分析的固定流程这套流程看起来不复杂但它是在反复返工中打磨出来的1. 搜索并确认CoNED版本检查分辨率、波段、基准面说明。 2. 打印元数据做直方图检查确认0值和NoData分布。 3. 画一条横穿海岸的剖面线目视检查海陆接缝质量。 4. 计算坡度标记潜在高低洼区和接缝异常。 5. 做多水位淹没模拟导出淹没深度图层。 6. 叠加人口、路网、用地数据统计暴露。 7. 对关键边界位置用剖面图做人工复核。 8. 导出全部中间结果和元数据存档。这套流程里最花时间的其实不是GEE代码而是第2到第4步的数据体检。一台手术做得好不好七分靠术前检查一份沿海高程分析靠不靠谱至少一半靠数据体检。最后说点个人体会用CoNED做沿海分析最深的体会是算法模型决定上限数据质量决定底限。GEE里的淹没模拟脚本满地都是但你换一份基准面没对齐的数据进去所有输出可能全都存在系统性偏差——这种偏差不是靠更复杂的模型能救回来的它会跟着每个结果走最后传导到决策依据上。如果你刚开始接触这类数据我建议从一个很小的范围下手挑一段几公里长的海岸加载CoNED画一条剖面线做一个水位模拟再手动把模拟边界和卫星影像上的湿地、滩涂边界叠在一起看。做这一步只需要一两个小时但对理解沿海高程数据的脾气比读十篇文档都管用。等这一步顺畅了再往人口暴露、多情景对比这些方向扩展会顺手很多。