2021FME博客大赛 —— SAR在土地利用变化检测中的应用

作者：王会

1 摘要

广州多云多雨的气候，仅凭光学影像难以满足全域多时序观测。合成孔径雷达(SAR)作为一种主动式微波传感器，其具有全天候、全天时成像，可穿透云雾等特点，能弥补光学遥感在常年多云多雨地区数据获取问题，非常适用广州这样的地区。

基于FME实现基于哨兵1号的包括数据处理、滤波、变化检测、统计分析的变化检测全过程。经核查，提取的变化图斑均为有效变化图斑，表明该方法行之有效，进一步证实了SAR影像在变化检测领域的可操作性，本方案可应用在土地执法、灾害预警等工作中。

2 总体思路

合成孔径雷达(SAR)作为一种主动式微波传感器，由于具有全天候、全天时成像的特点，尤其适用于南方这样多云多雨地区。通常，在南方很难获得完整的高分辨率遥感影像，但是对于3米分辨率的光学影像，可以在每个月拼出一张。对于光学遥感影像来说，做变化检测，3米分辨率是无法保证精度的，但是SAR影像具有穿透云雾的特点，且观测不分昼夜，因此基于多时相的SAR影像，可以快速得到全域范围变化范围，获取地表变化情况，此时在利用3米分辨率的影像进行人工判读，即可得到图斑变化类型。

本方案主要分为两部分来实现，首先是利用FME对哨兵1号SAR数据进行变化检测，提取变化图斑，然后是利用FME对提取的变化图斑和两期光学影像进行套合切图，生成变化监测专题图以及统计结果。

 

                                                                                                                               图 1 总体实现流程

 

 

                                                                                                                                图 2 总体实现FME模板

 

3 方案实现

3.1 数据准备

哨兵1号的数据可以来EARTHDATA网站（https://search.earthdata.nasa.gov/search）下载。这里需要注意的是，哨兵1号数据分为Level-0、Level-1、Level-2三级产品，这里我们使用Level-1级产品，在Level-1中，又分为SLC和GRD产品，GRD产品包含有经过多视处理、采用WGS84椭球投影至地距的聚焦数据，因此其数据本身是包含了GCP参数的，如果是下载的SLC产品，则需要先用SNAP等专业软件进行配准和地理编码，才能在FME里使用。

本次工作我们选用2020年12月21日和2021年02月19日两期的数据，下载下来后，数据结构如下图所示。

 

                                                                                                                                      图 3 数据结构

这里有两种方式打开数据，第一种是直接在FME里选择SLC读模块，然后选择manifest.safe文件；第二种是直接打开measurement文件夹下的tiff即可，其实，GRD产品包含了VH和VV两种极化方式数据，均在measurement以tiff保存。

 

                                                                                                                               图 4 SAR数据示意图

这里我们直接把两期的VV极化数据拖进FME里，用GeoTiff读模块打开。

3.2 数据预处理

前面讲了我们拿到的GRD产品是包含GCP信息的，所以数据预处理这一步最重要的就是对数据进行GCP纠正。

 

                                                                                                                               图 5 数据预处理FME模板

1）FME Hub里提供了RasterGCPApplier转换器，可以自动将GCP信息提取出来并对影像进行纠正。

2）由于前后两期的分辨率其实是不一致的，在FME Feature Information Window中可以看到这两期分辨率，分别是9.456619555294263e-5，9.464697119229038e-5，因此我们把分辨率高的那一景降采样到低的那一景，便于后面做波段运算。

3）裁剪。整景影像参与运算计算量很大，这里先进行裁剪，由于我的数据是广州市区县行政区的，所以我先用Dissolve转换器先合一下，然后在对SAR影像进行裁剪。

3.3 变化检测

这里最重要的就是如何提取变化图斑了。由于SAR是利用回波的相位和强度信息重建出被测区域的散射强度，因此若该区域没有发生变化，那前后两期的强度比值应该在1附件（不完全为1，有成像误差）；相反，若发生了变化，后一时相与前一时相的比值应远小于1或远大于1。

 

                                                                                                                              图 6 变化检测FME模板

1）在SAR图像中，斑点噪声是其很重要的误差来源，通常可以用中值滤波消除独立的斑点噪声，并能能够对整幅影像进行平滑，FME提供了RasterConvolver转换器，可以对影像进行卷积。

 

                                                                                                                            图 7 利用FME进行中值滤波

2）利用RasterExpressionEvaluator，将两幅影像做比值，得到变化比率。

 

                                                                                                                                     图 8 波段运算

3）再做一次中值滤波平滑数据，然后再次利用RasterExpressionEvaluator转换器，将[0.5,2)以为的栅格设置为0。此时，非0栅格即为变化区域。

4）先将0设为为NoData，在利用RasterToPolygonCoercer和Dissolver转换器将第（3）步的结果转成矢量，到这里，变化图斑就提取出来了，接下来基于对应两期的光学影像进行切图并生成报告。

3.4 切图

对于SAR来说，他的最重要的作用就是可以在大面积遥感影像中快速帮你找到变化线索，然后在利用中低分辨率的影像人工进行判读，这个时候，如果可以把每个图斑对于的前后两期的光学影像给裁出来，在进行判读，可以提升判读效率。

实现思路则是首先提取每个图斑的最小外接矩形，然后做一个50米的缓冲区，然后利用Clipper进行裁剪，紧接着用MapnikRasterizer进行制图，最后利用FeatureWriter分类输出就完成了。

 

                                                                                                                                图 9 切图效果示意图

 

 

                                                                                                                               图 10 切图FME模板

1）先用Reprojector将数据转成平面坐标系（因为FME的默认参数单位是Ground Units，和坐标系单位保持一致），在用Counter转换器给每个图斑编号，方便后面外接矩形和图斑的对应。然后用BoundsExtractor提取出每个图斑的外接矩形左下和右上角的坐标，利用2DBoxReplacer构建一个向外拓展50米的外接矩形，在利用Clipper进行裁剪。

 

                                                                                                                                图 11 外接矩形构建

2）利用MapnikRasterizer进行制图，设置原始图斑样式为Line，红色，Group By设置为_count，这样MapnikRasterizer就会将裁剪的影像和变化图斑进行一一对应，并分别输出。

3）在利用FeatureWriter进行输出，同样的方法，对第二期光学影像进行处理，输出到另外一个文件夹即可。这样，每个变化图斑我们可快速前后期的影像，在人工判别的时候可大幅提升效率。

    

                                                                                                                                      图 12 输出结果

 

3.5 统计

目的是为了统计各行政区的变化图斑数量和面积。

 

                                                                                                                            图 13 统计FME模板

1）用将行政区界限读进来，然后利用Clipper进行裁剪。

2）然后用AreaCalculator转换器重新计算面积（这里主要要转成平面坐标系），在把_area构建成一个List，Group By设置为行政区名称，在用ListSummer和ListElementCounter分别统计各区变化图斑面积和个数。

4 广州市SAR变化检测结果

本次使用2020年12月21日和2021年02月19日两期的哨兵1号数据对广州市进行变化检测，共检测到15885个变化图斑，总面积达26.95平方公里。

 

                                                                                                                                                    图 14

                                                                                                                                           表 1 各区统计结果

行政区	南沙	天河	越秀	荔湾	番禺	增城	花都	黄埔	白云	从化	海珠
变化面积 （平方公里）	9.45	0.34	0.03	0.31	2.65	4.19	3.01	2.52	2.22	1.88	0.37
图斑个数	3555	251	26	205	1355	2984	2243	1335	1658	2033	240

 

从图13可以看到在水体、港口、林地是存在大量变化的，这些变化可能是因为船离港进港、树木砍伐等原因造成的，并不一定是用地性质的变化，这个到时候要根据具体需求在叠合一些GIS数据将其剔除。

同时，我们也结合对应时期的光学影像进行了核对，发现均为有效变化。

第一时相	第二时相

虽然最终提取的图斑在1万多，其实很多图斑可以进一步合并，并且1万多的图斑，加上我们切好的图，核对的时候还是比较方便的。另外就是在没有质量高的光学影像时，该方法可以作为变化检测的一个有效补充手段。

5 不足

1）实际上对于提取出的变化图斑，有很多细碎图斑其实是反映了同一块用地的变化，因此可以将其合并，这样可以大幅减少提取的图斑数量，也便于后期的人工判别。

2）本次变化检测主要是采用比值算子做运算，可以更好的抑制相干噪声，但是专业软件里，还有对数比算子、均值比算子，而这些算子也很容易在FME里实现，后续也会进行实现，对比一下检测精度。

6 结语

已经在工作中使用FME有2年了，在这2年里，用FME做了很多事情，爬虫、大数据分析、格式转换还有一些基于FME的系统开发，这次也是因为一些工作尝试使用FME来进行SAR的处理。在之前一直都是用一些专业软件来处理SAR数据，但是每次都要把过程重复一遍，很难形成一个标准化的工作流程。其实，这也是FME的优势所在，最开始我是抱着试试看的态度在FME里处理SAR，发现在FME里可以做卷积，可以直接读取SAFE格式的SAR数据，可以GCP纠正，经过这次实验，FME包括FME Hub在遥感方面其实还是提供了很多算子的。

FME的优势绝不仅仅在于你知道的领域，可以说在你做其他领域工作时，不放也来试一下FME，一定会给你惊喜！在这里，我还要感谢FME中国技术交流群里的各位技术达人，也正是与他们一次一次的交流，让我能顺利完成这项研究，谢谢FME！