R语言将多景遥感影像拼接在一起的方法

李全胜 發表於 2024-7-11 15:03:00

R语言将多景遥感影像拼接在一起的方法

本文介绍基于R语言中的<code>raster</code>包，遍历文件夹，读取文件夹下的大量栅格遥感影像，并逐一对每一景栅格图像加以拼接、融合，使得全部栅格遥感影像拼接为完整的一景图像的方法。
 其中，本文是用R语言来进行操作的；如果希望基于Python语言实现类似的批量拼接、镶嵌操作，大家可以参考Python arcpy创建栅格、批量拼接栅格与Python ArcPy批量拼接长时间序列栅格图像这两篇文章。
 首先，来看一下本文所需实现的需求。如下图所示，现有一个文件夹，其中含有大量栅格遥感影像；这些遥感影像均为同一成像时间、不同空间范围的遥感影像。我们希望做到的，就是对这些遥感影像加以拼接，最终的结果图像就是一景将这里各个图像拼接后的大图像。
<img src="https://img2024.cnblogs.com/blog/3080295/202407/3080295-20240711150040397-2044068064.png">
 明确了需求，我们即可开始代码的撰写。本文所用到的代码如下所示。
<pre><code class="language-r">library(raster)
tif_file_name <- list.files(path = r"(E:\02_Project\01_Chlorophyll\Select\Result)", pattern = ".tif$", full.names = TRUE, ignore.case = TRUE)
tif_file_list <- list()
for (i in 1:length(tif_file_name)){
tif_file_list <- raster(tif_file_name)
}

tif_file_list$fun <- max
tif_file_list$na.rm <- TRUE
tif_mosaic <- do.call(mosaic, tif_file_list)
plot(tif_mosaic)

# tif_merge <- do.call(merge, tif_file_list)

rf <- writeRaster(tif_mosaic, filename = r"(E:\02_Project\01_Chlorophyll\Select\NewClip\LCC_SC_3.tif)", overwrite = TRUE)
</code></pre>
 首先，需要通过<code>library(raster)</code>代码，导入本文所需的R语言<code>raster</code>包；关于这一包的配置，大家可以参考基于R语言的raster包读取遥感影像。接下来，我们通过<code>list.files()</code>函数，遍历指定文件夹，从而获取当前文件夹下所包含的全部<code>.tif</code>格式的遥感影像，也就是全部待拼接的遥感影像。
 接下来，我们需要为栅格遥感影像的拼接做准备——也就是<code>for</code>循环内部的内容。此时，<code>tif_file_name</code>变量中存放的是指定文件夹下的全部栅格遥感影像的文件名称，而不是遥感影像文件自身；而接下来我们进行拼接、融合的函数，都需要保证函数参数中的遥感影像是一个栅格对象（<code>Raster* object</code>）类型的变量。因此，我们需要在这个<code>for</code>循环中，通过<code>raster()</code>函数，将每一个遥感影像的文件名（字符串类型）转为栅格对象类型。至于什么是栅格对象类型的变量，我们可以参考下图：其中<code>Formal class RasterLayer</code>即表示这一变量为栅格对象类型的。
<img src="https://img2024.cnblogs.com/blog/3080295/202407/3080295-20240711150037415-1311968084.png">
 接下来，代码分为<code>2</code>个部分。其中，<code>for</code>循环后的<code>4</code>行代码是第一部分，为栅格拼接的代码；同时为了对比栅格拼接与栅格融合的操作，这里还将栅格融合的代码也一并列出了，也就是注释掉的那一行代码。
 我们首先来看第一部分代码，这里通过<code>mosaic()</code>函数来实现栅格遥感影像的拼接。这一函数原本的参数中，只有<code>2</code>个栅格对象（<code>Raster* object</code>）类型的参数，换句话说就是原本这个函数只能同时拼接<code>2</code>个栅格遥感影像；如果我们有更多的遥感影像，就需要每一次拼接<code>2</code>个栅格图像，不断重复这一操作，直到全部的栅格遥感影像拼接完毕。这样操作无疑是比较麻烦的，因此我们需要借助<code>do.call()</code>函数来实现<code>2</code>个以上栅格的拼接工作——这个<code>do.call()</code>函数可以接受可变数量的参数，例如本文中我们需要对大量栅格遥感影像加以逐一拼接，具体有多少景遥感影像我们自己也不一定确定，且也不关心；因此就结合这一函数，将刚刚已经转为栅格对象（<code>Raster* object</code>）类型的图像所组成的列表<code>tif_file_list</code>作为参数，用<code>do.call()</code>函数来调用<code>mosaic()</code>函数，直到将<code>tif_file_list</code>列表中全部的栅格对象（<code>Raster* object</code>）类型的元素都带入到<code>mosaic()</code>函数运行后，<code>do.call()</code>函数就结束了。
 此外，由于<code>mosaic()</code>函数在运行时，除了两个栅格对象（<code>Raster* object</code>）类型的参数，还有其他的一些辅助参数，比如拼接时重叠区域该如何处理、处理时是否考虑NoData值的影响等；由于我们时通过<code>do.call()</code>函数来调用<code>mosaic()</code>函数，因此这些参数就不太好直接指定了。因此，我们可以通过<code>$</code>运算符，将<code>mosaic()</code>函数所需要的其他参数一并放入<code>tif_file_list</code>中，在后期<code>do.call()</code>函数调用<code>mosaic()</code>函数时，将同时读取这些参数，起到将参数传递到<code>mosaic()</code>函数中的功能。其中，在本文中我们需要指定<code>mosaic()</code>函数的<code>fun</code>参数与<code>na.rm</code>参数，二者分别是指拼接时重叠区域像元值的计算方法，以及计算重叠区域像元值时，是否考虑NoData值的影响；我们将这<code>2</code>个参数分别设定为<code>max</code>与<code>TRUE</code>，二者分别是指重叠区域的像元以<code>2</code>景遥感影像中的最大值像元为准，以及在计算时不考虑NoData值的影响。
 接下来，就是第二部分，即栅格融合的代码；在这里，我们通过<code>merge()</code>函数来实现遥感影像的融合。其实，这里的<code>merge()</code>函数与前述的<code>mosaic()</code>函数功能大致一样，但<code>merge()</code>函数在处理重叠区域时，默认选择位于顶层的遥感影像的像元数值，就没有<code>mosaic()</code>函数中的这么多计算方法选择了。
 最后，这里末尾的一句代码，就是将结果图像通过<code>writeRaster()</code>函数加以保存；这句代码的解释大家同样参考R语言求取大量遥感影像的平均值、标准差：raster库这篇文章即可。
 随后，运行上述代码，我们就可以获得将指定文件夹内全部栅格遥感影像加以拼接（执行代码中的第一部分）或者融合（执行代码中的第二部分）的结果了。
 至此，大功告成。 
来源：https://www.cnblogs.com/fkxxgis/p/18296251

頁: [1]

琼殿技术论坛's Archiver

R语言将多景遥感影像拼接在一起的方法