17.10. 图像预处理的应用

本部分介绍遥感中常用的图像预处理应用。 目标是从原始的、未经处理的产品(几何学和辐射测量学)生成可在地理信息系统中使用的图像。

17.10.1. 光学图像的辐射处理

本节回顾了纠正多光谱卫星图像辐射测量的常用方法,主要受到Ose 等人 的启发。

[OSE 16]. Atmospheric effects and a number of other mechanisms alter theradiometry of images: topographic effects, atmospheric scattering,directional effects (related to solar angle, observation angle, or thenature of the surface) .本章中我们不提供可用的纠正方法。

传感器接收的光能量是与大气层和地球表面多次相互作用的结果。在大气中,扩散和吸收机制减弱辐射。 能量然后被地球表面吸收、传输和/或反射(图5.21):

  • reflected energy, surface radiance, or Top Of Canopy (TOC) radiance crosses the atmosphere a second time;

  • then the sensor records an apparent radiance or Top Of Atmosphere (TOA) radiance.

_images/image205.png

图5.21.太阳辐射、大气和地球表面之间的相互作用(来自 [OSE 16])

在遥感中,将原始像素转换为TOA或TOC反射率(这是反射能量的一部分,相对于射入能量) 通常用于处理时间序列上的辐射指数分析:

  • Raw pixel to radiance: the relationship between the numeric value recorded by the sensor and the radiance is a linear function, for each channel. Parameters of this set of functions derive from the sensor calibration curve, and are generally available in the metadata of the product.

  • Radiance to TOA reflectance: the TOA reflectance is the part of energy reflected to the sensor by the Earth's surface and the atmosphere. It thus has continuous values in the range [0, 1] that enables the relative comparison between two or more images. The relationship between the TOA reflectance and the TOA radiance is a linear function for each channel. The parameters derive from the Earth– Sun distance, solar illumination and elevation angle. They are generally available in the metadata of the product.

  • TOC reflectance: the TOC reflectance depends of the reflectance and the transmittance of the atmosphere. Thus, the calculation is difficult due to two main difficulties: first, the atmosphere conditions are not homogeneous. Second, they are not always available inside the ROI. There are two main approaches to calculating TOC reflectance:

  • radiative transfer modeling (6S method): this method requires meteorological information (pressure, temperature, quantities of the different gases, composition of the aerosols, etc.), optical properties (distance of visibility, degree of transparency, etc.), and surface properties (homogeneous or non-homogeneous reflectance, directional effects, land use type, etc.) plus the parameters for computing the radiance;

  • correction using invariant targets: this set of methods rely only on available information, which is inside the image, to deduct atmosphere properties. Among the most popular, the Dark Object Subtraction method is based on targets that are radiometrically invariant (e.g. hill shades) to compute correction parameters.

在下面的练习中,我们使用 OpticalCalibration 应用程序计算Spot 6图像的反射率。

步骤说明

浏览文件

浏览文件系统以到达要进行放射性修正的图像位置。

在Spot 6产品目录树中,图像的名称如下: IMG_... TIF(用于TIFF图像)或IMG _... JP2(用于jpeg2000图像)。 其元数据命名如下:DIM _..包含用于计算OTA反射率的参数的HTML。

图像和元数据文件必须在Dimap兼容性表格中具有名称:

  • DIM_SPOT6_MS_ ..._1.XML for metadata file

  • IMG_SPOT6_MS_ ..._1_RxCy.TIF for image

如果文件的命名不同,OTB 函数可能无法正确访问metrix文件。 在这种情况下,用户必须指定所有参数(日期、角度、收益、偏差等)手工

注:OTB 支持以下传感器的光学校准:QuickBird、Ikonos、WorldView(1、2和3)、 Formosat、Spot 4/5/6/7、Pleiades。对于这些传感器,所有参数都在产品的元数据中自动读取。

显示帮助

执行 otbcli_OpticalCalibration (不带参数)以显示应用程序描述:

otbcli_OpticalCalibration

然后打印应用程序描述:

This is the OpticalCalibration application, version
5.2.1
Perform optical calibration TOA/TOC (Top Of
Atmosphere/Top Of Canopy). Supported sensors: QuickBird,
Ikonos, WorldView2, Formosat, Spot5, Pleiades, Spot6.

For other sensors the application also allows providing
calibration parameters manually.
Complete documentation: http://www.orfeo-toolbox.org/Applications/OpticalCalibration.html
Parameters:
 -progress <boolean>
Report progress
MISSING -in <string> Input
(mandatory)
MISSING -out <string> [pixel]
Output
[pixel=uint8/uint16/int16/uint32/int32/float/double]
(default value is float) (mandatory)
 -ram <int32>
Available RAM (Mb) (optional, off by default, default
value is 128)
 -level <string>
Calibration Level [toa/toatoim/toc] (mandatory, default
value is toa)
 -milli <boolean>
Convert to milli reflectance (optional, off by default)
 -clamp <boolean> Clamp
of reflectivity values between [0, 100] (optional, on
by default)
 -acqui.minute <int32>
Minute (mandatory, default value is 0)
 -acqui.hour <int32> Hour
(mandatory, default value is 12)
 -acqui.day <int32> Day
(mandatory, default value is 1)
 -acqui.month <int32> Month
(mandatory, default value is 1)
 -acqui.year <int32> Year
(mandatory, default value is 2000)
 -acqui.fluxnormcoeff <float> Flux
Normalization (optional, off by default)
 -acqui.sun.elev <float> Sun
elevation angle (°) (mandatory, default value is 90)
 -acqui.sun.azim <float> Sun
azimuth angle (°) (mandatory, default value is 0)
 -acqui.view.elev <float>
Viewing elevation angle (°) (mandatory, default value
is 90)
 -acqui.view.azim <float>

Viewing azimuth angle (°) (mandatory, default value is
0)
 -acqui.gainbias <string> Gains
| biases (optional, off by default)
 -acqui.solarilluminations <string> Solar
illuminations (optional, off by default)
 -atmo.aerosol <string>
Aerosol Model
[noaersol/continental/maritime/urban/desertic]
(mandatory, default value is noaersol)
 -atmo.oz <float> Ozone
Amount (optional, on by default, default value is 0)
 -atmo.wa <float> Water
Vapor Amount (optional, on by default, default value is
2.5)
 -atmo.pressure <float>
Atmospheric Pressure (optional, on by default, default
value is 1030)
 -atmo.opt <float>
Aerosol Optical Thickness (optional, on by default,
default value is 0.2)
 -atmo.aeronet <string>
Aeronet File (optional, off by default)
 -atmo.rsr <string>
Relative Spectral Response File (optional, off by
default)
 -atmo.radius <int32>
Window radius (adjacency effects) (optional, off by
default, default value is 2)
 -atmo.pixsize <float> Pixel
size (in km) (optional, on by default, default value is
1)
 -inxml <string> Load
otb application from xml file (optional, off by
default)
Examples:
otbcli_OpticalCalibration -in QB_1_ortho.tif -level toa
-out OpticalCalibration.tif

设置参数

如果输入图像符合dimap格式,则自动设置“acqui”参数组。参数组“atmo”仅用于TOC反射率。

其他参数:

  • “level” refers to the type of correction (the possible values are “toa”, “toc” or “toatoim”;

  • the value of a reflectance is a floating decimal number, less than 1. The “milli” parameter is used to multiply this value by 1,000 so that the values can be encoded on 16 bits, enabling a saving in disk space (half as many bits as floating).

在本练习中,仅设置以下参数:

  • “in”: filename of the input image “ROI_2_XS_Ortho_GSD2016”;

  • “level” to “toa”;

  • “out”: filename for the output TOA reflectance image.

执行

运行算法。

分析

在原始图像(TOA修正之前)和TOA反射率图像上计算NDVI。

比较植被上的两个值:在QGIS中打开两个输出图像并使用“值工具”(来自“查看”>“面板”)。 如果不可用,请使用插件管理器搜索并安装)。

表5.15.使用“OpticalCalibration”应用程序进行光学校准

17.10.2. SAR辐射测量处理

合成口径雷达(SAR)是一种用于通过利用雷达的运动来模拟更大的天线以更高的分辨率来数字化图像的雷达。 图5.22描述了移动平台观察目标P三次但采集角度不同的采集原理。 SAR处理结合了这些获取,允许您提高长度L的原始天线的分辨率, 以获得长度M的等效合成天线,该天线对应于方位角方向上的平台路径。 回声目标P将位于射程维度中的距离h处。

SAR原理是一种主动微波遥感方法,允许用户白天或晚上捕获图像,在操作应用中具有明显的优势。 它被广泛用于从陆地到海洋应用的多种环境中。例如,我们可以列出农业监测(通常与光学观察相结合)、 洪水区检测、漏油(例如海上非法排放)或冰监测。

还有更先进的方法可以从SAR采集原理中派生出来; 我们可以提到干涉测量法,它使您能够生成数字地形图或测量小位移或旋光测量法, 使您能够描述诸如植被类型等物体的特征。

_images/image206_xwa.jpeg

图5.22. SAR采集原理。

移动平台上的雷达多次看到目标P(距离h),提供与长度M的合成天线相对应的分辨率(来源;维基百科)

SAR图像的辐射处理能够生成与反向散射信号直接相关的像素值。 事实上,反向散射信号的行为取决于表面特征, 尤其是凹凸度、土壤湿度以及几何形状(偶极、角反射器)。

对于光学图像,辐射校准将数字转换为物理量。 一般来说,1级地球观测产品不包括这些更正。 因此,应用这种校准来获取反向散射信号的物理值和表面属性非常重要。 它还允许您在来自潜在不同传感器的不同图像之间执行协调。

SAR校准主要包括计算以下量:

  • Beta Nought ( °): or radar brightness, which corresponds to the reflectivity per unit area, expressed in slant range.

  • Sigma Nought ( °): this corresponds to the radar backscatter coefficient and is directly linked to the observed surface properties. This represents the Equivalent Radar surface of the pixel considered with a flat surface hypothesis. This unitless value is generally expressed in decibels (dB).

  • Gamma_Nought ( °): radar backscatterer coefficient normalized by the incidence angle.

为了计算这些操作,SAR图像产品通常提供元数据,允许您将像素值转换为这些物理量。

这些量可以使用OTB 应用程序 SARCalibration 计算。 该应用程序的文档可在CookBook [OTB 17 c]的“SAR处理”部分中找到。 支持 TerraSAR-XSentinel-1Radarsat-2 传感器。

17.10.3. 遥感图像的几何处理

矫形矫正

卫星图像的正射纠正包括将其呈现为在最低点(即垂直)处获取的图像。 此操作的目的是消除图像透视(或磁贴,如图5.23所示)和浮雕(图5.24)的影响。 正射矫正产生的图像能够测量距离、角度和表面。

_images/image207_xiy.jpeg

图5.23.倾斜位移是物体在倾斜照片上的图像位置与其在真正垂直照片上的理论位置的移动。 这是由于曝光时照片平面相对于基准平面倾斜的结果。这种效果在上图中得到了证明。 尽管A-B和C-D之间的距离在基准平面上相同, 但它们在照片平面上的相应表示则不同(即a带c-d之间的距离不相等)(来源:OSSSIM)。

_images/image208_xkz.jpeg

图5.24.浮雕位移是物体在特定基准上方的高度引起的图像位置的移动。 对于垂直或近垂直摄影,移动从最低点呈放射状发生。这种效果如上图所示。 尽管A-B和C-D之间的距离在基准上相同

平面,它们在照片平面上的相应表示不相同(即A-B和C-D之间的距离不相等)(来源:OSSSIM)

通常,封装非正射纠正产品的图像文件不具有地理参考元数据。元数据通常在另一个文件中提供, 专门用于封装。例如,Spot 6 Raster Sensor图像的标题不包含任何地理参考信息; 所有与采集相关的元数据都包含在ML文档(dimap格式)中。 该图像无法直接在地理信息系统中使用,因为它没有附加坐标参考系。 只能以像素坐标显示(图5.25(a))。正射矫正包括对图像应用几何变换。 经过正射纠正后,一个基准(测地系统)就会附加到图像上。 包含正射纠正图像(例如,GeoTivf)的文件就具有描述空间信息的元数据。 地理信息系统中常用的元数据包括:

  • spatial reference system;

  • upper left pixel (image origin) coordinates in spatial reference system units;

  • pixel size in both X and Y dimensions in spatial reference system units.

该空间信息可以计算每个像素的地图坐标。 这允许将正射纠正图像叠加到地图或任何其他正射纠正图像(图5.25(b))。

图5.25. a)多光谱栅格传感器图像(斑点6); b)正射矫正后的图像

为了对卫星图像进行正射纠正,通常使用以下信息:

  • 图像和传感器模型;

  • DEM .

传感器模型描述了所有传感器采集参数,因此允许在地球表面上精确定位像素。 物理传感器模型源自轨道信息(获取时传感器的位置和方位), 但还有其他模型包括大气影响或接收站进行的处理(例如,在地图投影中对图像进行重新采样时)。 通常,模型的实施是特定于特定传感器的(例如,扫帚、推扫帚、针孔、帧摄像头、合成孔径雷达等)。 也有用数学近似(如多项式、几何级数等)表示的传感器模型, 它们取代了模拟物理变换的方程(“替换模型”)。另一类传感器模型包括提供图像的变换网格。 这些模型总是从采集的物理参数中得出, 但它们的目标是提高执行速度而损害精度(变换是简单的内插,而不是复杂的、严格的和非线性的建模)。 数字地形模型是适于由图像处理工具(即,通常是栅格)使用的陆地区域的地形的表示。 数字地形模型需要能够在三维空间中定位地球表面上的一个点,这是使用传感器模型所必需的。

在下面的练习中,我们对Spot 6 Raster传感器图像执行正射纠正(“ROI_3_P_Sensor_GSD2016”)。 将使用产品的传感器模型(如元数据中定义的)以及数字降雨模型。 接下来,OTB 应用程序将在命令线上被调用。

首先,我们使用 DownloadSRTMTiles 应用程序下载正射纠正图像所需的数字地形模型部分。 该应用程序确定原始图像及其相关传感器模型的坐标,以定位要下载的数字地形模型的相关部分。

步骤

指令

创建一个文件夹

Create a “ DEM _Tiles” directory on the file system. This directory willbe used to host the digital elevation model files that will bedownloaded by the DownloadSRTMTiles application.

在文件系统上创建“ DEM _Tiles”目录。 此目录将用于托管将由DownloadSRTMTiles应用程序下载的数字海拔模型文件。

获取EM

使用下载SRTMTiles应用程序检索数字地形模型的有用部分。作为此应用程序的输入,将输入以下参数:

“il”:要检索数字模型的原始图像(“ROI_3_P_Sensor_GSD2016”);

“mode.download.outdir”:之前创建的目录“ DEM _Tiles”, 我们想要在其中存储数字海拔模型的文件。

执行以下命令(将 &lt; &gt; 替换 DEM _Tiles 为包含数字海拔模型切片的目录的完整路径):

otbcli_DownloadSRTMTiles -mode.download.outdir << DEM _Tiles>> -il
rawInputImage.TIF

执行

运行命令行。

该应用程序显示与图像覆盖的区域对应的下载文件的名称。

2017年3月22日14:05:37:Application.logger(Info)在USGSserver上找到瓦片, 网址:http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Africa/N48W003.hgt.zip

2017年3月22日14:05:38:Application.logger(Info)在USGSserver上找到瓦片, 网址:http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Africa/N48W004.hgt.zip

2017年3月22日14:05:39:Application.logger(Info)在USGSserver上找到瓦片, 网址:http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Africa/N49W003.hgt.zip

2017年3月22日14:05:39:Application.logger(Info)在USGSserver上找到瓦片, 网址:http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/Africa/N49W004.hgt.zip

表5.16.使用“DownloadSRMTTiles”应用程序下载SRTM数字海拔模型瓦片

然后,可以使用 Orthorectification 应用程序。

步骤

指令

显示应用程序描述

执行

otbcli_Orthorectification ,没有参数来显示应用程序描述:

otbcli_OrthoRectification

然后打印应用程序描述:

This is the OrthoRectification application, version 5.2.1
This application allows ortho-rectification of optical images from
supported sensors.
Complete documentation: http://www.orfeotoolbox.org/Applications/OrthoRectification.html
Parameters:
-progress <boolean> Report progress
MISSING -io.in <string> Input
Image (mandatory)
MISSING -io.out <string> [pixel] Output
Image [pixel=uint8/uint16/int16/uint32/int32/float/double]
(default value is float) (mandatory)
-map <string> Output
Cartographic Map Projection
[utm/lambert2/lambert93/wgs/epsg] (mandatory, default value is utm)
-map.utm.zone <int32> Zone number (mandatory, default value is 31)
-map.utm.northhem <boolean> Northern Hemisphere (optional, off by
default)
-map.epsg.code <int32> EPSG Code
(mandatory, default value is 4326)
-outputs.mode <string>
Parameters estimation modes
[auto/autosize/autospacing/outputroi/orthofit] (mandatory, default value
is auto)
MISSING -outputs.ulx <float> Upper Left X (mandatory)
MISSING -outputs.uly <float> Upper Left Y (mandatory)
MISSING -outputs.sizex <int32> Size X (mandatory)
MISSING -outputs.sizey <int32> Size Y (mandatory)
MISSING -outputs.spacingx <float> Pixel Size X (mandatory)
MISSING -outputs.spacingy <float> Pixel Size Y (mandatory)
-outputs.lrx <float> Lower right X (optional, off by default)
-outputs.lry <float> Lower right Y (optional, off by default)
-outputs.ortho <string> Model ortho-image (optional, off by default)
-outputs.isotropic <boolean> Force isotropic spacing by default
(optional, on by default)
-outputs.default <float> Default pixel value (optional, on by default,
default value is 0)
-elev.dem <string>  DEM  directory (optional, off by default)
-elev.geoid <string> Geoid File (optional, off by default)
-elev.default <float> Default elevation (mandatory, default value is 0)
-interpolator <string>
Interpolation [bco/nn/linear] (mandatory, default value is bco)
-interpolator.bco.radius <int32> Radius for bicubic interpolation
(mandatory, default value is 2)
-opt.rpc <int32> RPC modeling (points per axis) (optional, off by
default, default value is 10)
-opt.ram <int32> Available
RAM (Mb) (optional, off by default, default value is 128)
-opt.gridspacing <float>
Resampling grid spacing (optional, on by default, default value is 4)
-inxml <string> Load otb application from xml file (optional, off by
default)
Examples: otbcli_OrthoRectification -io.in
QB_TOULOUSE_MUL_Extract_500_500.tif -io.out
QB_Toulouse_ortho.tif

设置参数设置以下参数:

  • “io.in”: input image filename (that is, the image that we want to orthorectify, “ROI_3_P_Sensor_GSD2016”), followed by the extended filename “?&skipcarto=true”. This extended filename is required to enable the orthorectification of OrthoReady products, that is: QuickBird, GeoEye, WV2 and WV3[] products in the OrthoReady2A[] level, that Spot-6/7 Raster Sensor. These products have a pseudo-coordinate reference system which does not enable the use of images in a geographical information system straightforwardly. The “?&skipcarto=true” extended filename ignores these unwanted metadata.

  • “io.out”: path for the output image. Since Spot 6 images are encoded on 16 unsigned bits, one can set the pixel type to uint16.

  • “map”: choose the map projection “lambert93” used for metropolitan France.

  • “elev.dem”: “ DEM _Tiles” directory, containing the digital elevation model tiles.

参数的默认值适合当前可用的大多数图像。

然而,调整其中一些可能是有用的,以减少计算时间以进行快速预览:

  • “interpolator”: as explained in section 5.3.1.2, the processing time depends a lot on the interpolator. The nearest neighborhood interpolator is the fastest: select “nn” to minimize the calculation time.

  • “opt.gridspacing”: determines the mesh size of the deformation grid, in map units. The higher this value, the more coarse the deformation will be, but also the less costly in computation time. It should be noted that it is generally not very useful to define a mesh size smaller than the size of the pixel of the image. In principle, a mesh size of the deformation grid equal to the size of a few pixels is sufficient: for a panchromatic Spot 6 image at 1.5 m, it is reasonable to define a mesh size of 15 m, corresponding to 10 pixels. Small values of deformation grid size are relevant when the digital terrain model and the image pixel spacing are small and the terrain very rough.

执行

运行该过程(将 &lt; &gt; 替换 DEM _Tiles 为包含数字立面模型切片的目录的完整路径)。

otbcli_OrthoRectification
-io.in "PATH/TO/IMG_SPOT6_..._1_R1C1.TIF?&skipcarto=true"
-io.out roiOrtho.tif -map lambert93
-elev.dem << DEM _Tiles>> -interpolator nn
-opt.gridspacing 15

可视化将正射纠正图像导入QGIS并可视化。

评估 DEM 贡献

在不指定数字地形模型的情况下再次运行应用程序:这样,在正射纠正过程中仅使用轨道参数, 而不使用任何海拔信息。

otbcli_OrthoRectification
-io.in "PATH/TO/IMG_SPOT6_..._1_R1C1.TIF?&skipcarto=true"
-io.out roiOrtho_no_dem.tif -map lambert93
-interpolator nn
-opt.gridspacing 15

将输出图像导入QGIS并可视化。

直观地比较使用和不使用数字海拔模型的正射纠正图像的平面测量。

“正矫正”申请

用地面控制点进行正射矫正

最近的传感器模型通常能够以几米的相对误差实现图像参考。然而,这种准确性并不总是可以接受的, 尤其是当图像旨在进行小规模地理定位或变化检测时。下一节介绍地面控制点的使用, 这些点是用于改进传感器模型的参考点。因此,图像以接近参考图像的准确度进行正射纠正。 本节介绍了实现这一目标的方法。

第一步包括提取不变点对。 因此,每个对由要进行正射纠正的图像上的第一不变点和参考图像上的第二点构成。

通过参考图像,可以手动识别点对。但是,这种识别可以自动执行。 计算机视觉领域越来越多的算法使得检测和识别不同数字图像之间的相似元素变得更容易。 其中包括Sift[99]和SURF[BAY 06],它们用于当前的许多应用,如数字摄影中的场景自动创建(图5.26)

_images/image211_xd8.png

图5.26.通过SIFT方法匹配两个图像的结果示例 (Fantasia ou Jeu de la poudre,Devant la ported ' entée de la ville de Méquinez,作者:Eugène Delacroix,1782年)(来源:维基百科)

SIFS和SURF都在OTB 应用程序中实现,名为 HomologousPointsDetection。 无论选择哪种算法(SFT或SURF),方法如下:

  • Extraction of the keypoints from the images: the keypoints are scale and rotation invariant. A SIFT descriptor is a vector of 128 elements, as opposed to 64 for the SURF descriptor.

  • Keypoints matching: based on the Euclidean distance between the descriptors. The goal is to form pairs of closest descriptors. The homologous points on the two images are thus matched according to their descriptors.

  • Once a number of homologous points have been obtained, it is possible to slightly modify the parameters of the sensor model in order to better fit the reference points. This step consists of generating a new sensor model that minimizes the quadratic error between pairs of points extracted in the preceding step. The OTB RefineSensorModel application performs this operation.

  • The last step consists of performing the orthorectification with the new refined sensor model.

在以下练习中,将使用额外控制点对Spot 6 RasterSenser(“ROI_3_P_Sensor_GSD2016”) 图像进行正射纠正。光学图像(“ROI_3_P_Ortho_IGN 2015”)将用作参考图像。实施的步骤如下:

1)计算关键点(SIFS)来创建对,即地面控制点( HomologousPointExtraction 应用程序);

2)使用成对的地面控制点完善传感器模型( RefineSensorModel )应用程序;

3)使用精确的传感器模型进行正射矫正(Orthorection*应用)。

在以下示例中,应用程序将从命令行使用。这种方法提供的一个有趣的可能性是创建脚本来链接上述三个步骤的能力。

为了减少每个步骤的处理时间,我们将仅处理输入的Spot 6 Raster Sensor图像的一个子集。 ExtractROI 应用程序用于提取图像的部分以及该区域特定的传感器模型。 传感器模型是由 ExtractROI 从dimap元数据中自动提取的,并将其作为输出图像, 扩展名为.geom。

步骤

指令

提取投资回报率

使用 otbcli_ExtractROI 从Spot 6 Raster Sensor图像中提取任何1000 x1000像素的区域(无论起源如何)。

使用扩展名为“?&skipcarto=true”的输入图像的扩展文件名。

输入以下命令:

otbcli_ExtractROI -startx 200 -starty 200 -sizex 1000
-sizey 1000
-in "PATH/TO/IMG_SPOT6_..._1_R1C1.TIF?&skipcarto=true"
-out roi.tif

检查生成的文件使用文件资源管理器,确保已创建以下文件:

  • the output image, in GeoTIFF format (.tif);

  • the sensor model associated with the output image, in a .geom file. This file must have the same name as the output image except the extension.

ExtractROI应用程序应该从Dimap Spot 6中提取传感器元数据,并将其合成在.geom文件中。 后者可以通过简单的文本编辑器打开。

表5.18.使用“ExtractROI”应用程序使用相关传感器模型提取图像区域

在此步骤中,我们提取了输入图像的子集,以及特定于该区域的传感器模型。 然而,此操作的主要目标是缩小图像的大小以加速过程。 如果希望在整个图像(而不仅仅是一个小区域)上继续此实践的步骤, 则有必要提取特定于整个图像的传感器模型:为此,请使用 ReadImage Info 应用程序。 该应用程序有两个强制参数, 即用于输入图像的“in”参数(输入图像)和用于存储传感器模型的“outkwl”参数(输出关键字列表)。

在进一步讨论之前,让我们对之前提取的子集图像进行正射纠正。

步骤

指令

进行正射矫正

使用 otbcli_Orthorectification 矫正使用上一节中描述的方法对子集图像执行正射矫正。

注意:您不需要使用具有扩展名“?&”的扩展文件名skipcarto=true”, 因为在使用ExtractROI应用程序提取区域时已经完成了这一操作。

otbcli_OrthoRectification -io.in roi.tif -io.out roiOrtho.tif -map
lambert93
-elev.dem << DEM _Tiles>> -interpolator nn
-opt.gridspacing 15

可视化

将生成的图像导入QGIS

表5.19.使用“Orthorectification”应用程序对提取的图像进行正射矫正

下一步(如表5.20所示)包括使用OTB 应用程序同源点来提取点对。

步骤

指令

显示帮助

在没有参数的情况下执行 otbcli_Orthorectification 以显示应用程序描述:

otbcli_HomologousPointsExtraction

然后打印应用程序描述:

This is the HomologousPointsExtraction application, version 5.2.1

Compute homologous points between images using keypoints

Complete documentation: http://www.orfeo-toolbox.org/Applications/HomologousPointsExtraction.html

Parameters:

-progress <boolean> Report progress

MISSING -in1 <string> Input Image 1 (mandatory)

-band1 <int32> Input band 1 (mandatory, default value is 1)

MISSING -in2 <string> Input Image 2 (mandatory)

-band2 <int32> Input band 2 (mandatory, default value is 1)

-algorithm <string> Keypoints detection algorithm [surf/sift]
(mandatory, default value is surf)

-threshold <float> Distance threshold for matching (mandatory, default
value is 0.6)

-backmatching <boolean> Use backmatching to filter matches. (optional,
off by default)

-mode <string> Keypoints search mode [full/geobins] (mandatory, default
value is full)

-mode.geobins.binsize <int32> Size of bin (mandatory, default value is
256)

-mode.geobins.binsizey <int32> Size of bin (y direction) (optional, off
by default)

-mode.geobins.binstep <int32> Steps between bins (mandatory, default
value is 256)

-mode.geobins.binstepy <int32> Steps between bins (y direction)
(optional, off by default)

-mode.geobins.margin <int32> Margin from image border to start/end bins
(in pixels) (mandatory, default value is 10)

-precision <float> Estimated precision of the colocalisation function
(in pixels). (mandatory, default value is 0)

-mfilter <boolean> Filter points according to geographical or sensor
based colocalisation (optional, off by default)

-2wgs84 <boolean> If enabled, points from second image will be exported
in WGS84 (optional, off by default)

-elev.dem <string>  DEM  directory (optional, off by default)

-elev.geoid <string> Geoid

File (optional, off by default)

-elev.default <float>

Default elevation (mandatory, default value is 0)

MISSING -out <string> Output file with tie points (mandatory)

-outvector <string> Output vector file with tie points (optional, off by
default)

-inxml <string> Load otb application from xml file (optional, off by
default)

Examples:

otbcli_HomologousPointsExtraction -in1 sensor_stereo_left.tif in2
sensor_stereo_right.tif -mode full -out homologous.txt

该应用程序具有许多输入参数。

首先,有两个输入图像:

  • “in1”: image 1 (the image that we want to correct geometrically);

  • “in2”: image 2 (the reference image “ROI_3_P_Ortho_IGN2015”).

对于每张图像,一个额外的参数定义要使用的通道(分别为“band 1”和“band 2”)。 然后,参数“算法”允许您在两种算法之间进行选择,用于检测不变点(SIFS或SURF)。

无论选择的算法如何,以下参数都适用:

  • “threshold”: detection threshold between 0 and 1. When close to 1, a large number of pairs of points are created. When close to 0, a small number of pairs of points are created. The association of the points is done in the descriptor space (SIFT or SURF) and the threshold is applied to the ratio between the distance to the nearest descriptor and the distance to the second nearest descriptor;

  • “backmatching”: when this Boolean is TRUE, the pairs of points are formed only when the criterion has been previously verified from image 1 to image 2, and vice versa. This minimizes the number of bad matches;

  • “mode” allows you to choose between two processing modes:

    • in “geobins” mode, images will be processes in subregions;

      − “binsize”: size of subregions;

      − “binstep”: step between two subregions;

      − “margin”: margin around each subregion;

    • In “full” mode, images are processed entirely. We do not recommend processing full images because of the amount of memory it would require. In addition, keypoint pairing can be a very long operation when there are a big number of points.

  • “precision”: when “mfilter” is true, this value is the maximum distance (in physical units) between the points of image 1 and points of image 2.

当一对点之间的距离大于该值时,拒绝一对点;

  • “2wgs84”: when true, coordinates are exported in WGS84;

  • “out” is the file name for the output tabular file containing the points.

  • “outvector” is the file name for the output vector layer. It is optional, but useful to compare visually the ground control points that have been generated.

设置参数

提取相应的点。这些点将保存在.CSV格式的表格文件中, 并保存在.shp文件(ESRI Shapefile格式)中。

  • use the SIFT algorithm;

  • use the “geobin” mode by leaving its default subparameters;

  • reject pairs of points when the distance is greater than 30 m (set “precision” to 30 and set “mfilter” to 1);

  • do not forget to set the digital elevation model in “elev.dem”.

执行

运行命令行:

执行

运行命令行:

otbcli_HomologousPointsExtraction -in1 roi_sensor.tif -in2
PATH\TO\ROI_3_P_Ortho_IGN2015\IMG_S6P_2015032636879163CP_R1C 1.tif
-2wgs84 1 -mfilter 1 -precision 25 -mode geobins -out points.txt
-outvector deplacement_points.shp -elev.dem << DEM _DIR>> -backmatching 1

分析

将参考图像导入QGIS。

导入正射矫正后图像的子集。

最后,打开之前生成的矢量图层并分析两张图像之间的位移。 可以修改渲染风格(颜色、宽度)以更好地可视化矢量图层中的点。

表5.20.使用“HomologousPointsExtraction”应用程序提取同名点

我们现在可以继续进行步骤(2)。在这一步中,使用前一步中生成的控制点来细化传感器模型。 RefineSensorModel 应用程序生成一个新的传感器模型,可以最大限度地减少平方误差。

步骤

指令

显示帮助

在没有参数的情况下执行 otbcli_Orthorectification 以显示帮助:

otbcli_RefineSensorModel

然后打印应用程序描述:

This is the RefineSensorModel application, version 5.2.1

Perform least-square fit of a sensor model to a set of tie points

Complete documentation: http://www.orfeotoolbox.org/Applications/RefineSensorModel.html

Parameters:

-progress <boolean> Report progress

MISSING -ingeom <string> Input geom file (mandatory)

MISSING -outgeom <string> Output geom file (mandatory)

MISSING -inpoints <string> Input file containing tie points (mandatory)

-outstat <string> Output file containing output precision statistics
(optional, off by default)

-outvector <string> Output vector file with residues (optional, off by
default)

-map <string> Output

Cartographic Map Projection [utm/lambert2/lambert93/wgs/epsg]
(mandatory, default value is utm)

-map.utm.zone <int32> Zone number (mandatory, default value is 31)

-map.utm.northhem <boolean> Northern Hemisphere (optional, off by
default)

-map.epsg.code <int32> EPSG Code (mandatory, default value is 4326)

-elev.dem <string>  DEM  directory (optional, off by default)

-elev.geoid <string> Geoid File

(optional, off by default)

-elev.default <float> Default elevation (mandatory, default value is 0)

-inxml <string> Load otb

application from xml file (optional, off by default)

Examples:

otbcli_RefineSensorModel -ingeom input.geom -outgeom output.geom
-inpoints points.txt -map epsg -map.epsg.code 32631

设置参数

设置以下输入参数:

  • “ingeom”: the sensor model (.geom) to refine;

  • “inpoints”: the file containing the ground control points (output file “out” generated by HomologousPointsExtraction);

  • “outgeom”: file name of the new sensor model that will be created by the application (.geom);

  • “outstats”: file (table) in which statistics on the accuracy of points will be exported. Typically, for each ground control point, the shift (in physical units), error, elevation, etc.

执行

运行命令行:

otbcli_RefineSensorModel -ingeom roi_sensor.geom -outgeom
roi_sensor_refined.geom -outstat refine_stats.csv -outvector
residus_deplacement_points.shp -map lambert93 -elev.dem
<< DEM _DIR>> -inpoints points.txt

分析

在电子表格软件中打开“outstat”输出文件(图5.27),分析地面控制点的位置误差。

_images/image212_x1e.jpeg

图5.27.统计表

将统计数据与QGIS中的是两次层进行比较。

还可以在地面控制点文件(用于“inpoints”参数)中手动删除具有较大全局误差(即错误)的点, 然后再次运行细化传感器过程。

表5.21.使用“RefineSensorModel”应用程序完善传感器模型

表5.22中描述的步骤包括执行第5.3.2.3.1节中解释的正射矫正, 但使用新的传感器模型而不是原始的传感器模型。

步骤

指令

设置参数

在控制台中,输入与第5.3.2.3.1节所述相同的参数。

按照下面的说明完成命令行。

设置输入图像

设置输入图像参数,并将新传感器型号设置为使用扩展文件名,如下所示:

otbcli_OrthoRectification
-io.in "roi.tif?&refined_sensor.geom" -io.out roiOrtho.tif
-map lambert93 -elev.dem << DEM _Tiles>> -interpolator nn
-opt.gridspacing 15

其中:

  • roi.tif is the input image subset;

  • refined_sensor.geom is the new refined sensor model generated by the RefineSensorModel application.

执行

运行命令行。

分析

在QGIS中打开输出图像,使用参考图像将物理定位的精度与其他生产图像进行比较:

  • orthorectification without DEM ;

  • orthorectification with DEM .

表5.22.使用“正射矫正”应用程序,通过精细传感器模型进行正射矫正

17.10.4. 全色锐化

对地观测卫星通常具有一个全彩色传感器和一系列空间分辨率较低的多光谱传感器。 为了利用这两种信息来源,人们开发了多种技术来产生多光谱图像, 该图像具有与全彩色图像相同的空间分辨率。这种技术被称为锐化。

在本节中,我们将介绍在OTB 中实现的两种全景细化方法:相对分量替换(RC)和Bayesian Data Fusion(BDF)。

  • RCS method (Figure 5.28): a simple way to merge the images is to consider that at the same resolution, the panchromatic image is the sum of the channels of the multispectral image. Once the images have been resampled at the same resolution, the fusion can be carried out: the idea is to apply a low-pass filter to the panchromatic channel image to give it a spectral content (in the Fourier domain) close to the multispectral image. Then, the multispectral image is normalized by the filtered panchromatic image, and multiplied with the original panchromatic image. The only parameter of this method is the low-pass filter radius.

_images/image213_x62.png

图5.28. RCS Panelishering方法(来源:Cookbook)

BDF方法基于频谱带和全彩色通道之间的统计关系。最好处理感兴趣的区域,而不是整个图像。 我们不会深入探讨该方法的技术细节;但感兴趣的读者可以获得原始出版物[FAS 08]。 该方法的一个有趣特征是可以相对于频谱信息对空间信息进行加权,从而可以根据需要调整其使用。 例如,照片解释更加关注图像的清晰度:这可以通过赋予全彩色通道中包含的信息更大的权重来实现。 请注意,在他们的文章中,作者提出了一组适应城市、林业、农业和地中海环境的参数值。

在下面的练习中,叠加全景锐化 应用程序与Spot 6图像 (“ROI_1_Bundle_Ortho_GSD2015”)的多光谱通道(XS图像)和多色通道(SPAN图像)一起使用。 ReadImage Info 应用程序用于检查图像属性,例如像素间隔和光谱带数量。

步骤

指令

获取图像信息

在XS图像上使用ReadImage Info应用程序来获取图像属性:

otbcli_ReadImageInfo -in IMG_SPOT6_MS_..._1_R1C1.TIF

该信息由应用程序显示(图5.25)。请注意图像的大小(行数、列数)、 像素间距(以地图单位)和XS图像中的通道数量。

_images/image214_xhd.jpeg

图5.29.条带数量、大小和间距(蓝色)。有关此图的彩色版本, 请参阅www.iste.co.uk/baghdadi/qgis1.zip

对PAN图像重复该操作。

比较图像的大小和分辨率(PAN图像的分辨率必须是XS图像的四倍,XS图像的分辨率只有一个通道与四倍)。

提取ROI

为了加快全景清晰操作,我们将仅处理输入图像的一个子集。

使用ExtractROI应用程序提取此PAN图像子集:

otbcli_ExtractROI -in IMG_SPOT6_P_..._1_R1C1.TIF
-sizex 1000 -sizey 1000 -out pan_roi.tif

重新采样XS图像

在此步骤中,我们以与先前创建的PANimage ROI相同的范围和像素间距对XS图像进行重新采样。

无论采用哪种方法,这一步都是全景化的先决条件

选择了事实上,XS图像的内插不包括在Panshiring应用程序中,因为它是一个独立的过程。

选择半径为两个像素的双三次字形

在四倍大的PAN图像上插值XS图像以避免伪影(请参阅第5.3.2.3.1节)。

选择以下参数:

  • “Reference input”: set the path to the PAN image ROI (extracted during the preliminary step);

  • “Image to reproject”: set the path to the XS image to be resampled;

  • “interpolation”: choose “bicubic” with a radius of two pixels.

将其他参数保留为默认值。图像已经进行正射纠正,无需输入大地水准面或数字海拔模型。

重新分配并执行在PAN图像上重新分配XS图像:

otbcli_Superimpose -inm IMG_SPOT6_MS_..._1_R1C1.TIF -inr pan_roi.tif
-out ms_roi_resamp.tif

执行RCS

使用 Panshiring 命令行应用程序。指定图像的路径:

  • XS image: set the filename of the XS image resampled on the PAN image ROI in the previous step;

  • PAN Image: set the filename of the PAN image ROI.

将“方法”参数设置为“rcs”,指定输出图像的文件名并执行算法。

otbcli_Pansharpening -inp pan_roi.tif
-inxs ms_roi_resamp.tif -out pansharp_rcs.tif -method rcs

执行

重复上一步,但将参数“方法”更改为“Bayes”:

BDS

otbcli_Pansharpening -inp pan_roi.tif
-inxs ms_roi_resamp.tif -out pansharp_bds.tif -method bayes

分析

将PAN和XS图像导入QGIS。

从RCS和BDS方法导入图像并进行视觉比较。

表5.23.使用“Pansheling”应用程序进行光学图像的Pansheling

17.10.5. 图像镶嵌

与要覆盖的广阔区域相比,遥感器的足迹通常很小:因此, 需要生成卫星图像镶嵌来将大量小足迹场景组合到单个镶嵌中。 Mosaic 应用程序满足了这一需求。可选的矢量数据可用于选择输入图像边界。 与大多数涉及图像重采样的OTB 应用程序一样,插值方法可以出现(双三次、线性或最近邻居)。

另一方面,土地利用、表面照明、大气条件和传感器引入的邻近图像之间存在颜色差异。 该应用程序包括一种颜色协调方法,可以最大限度地减少图像重叠部分中的像素差异(图5.30)。 这种差异可以在图像的原始辐射色彩空间中或在适合自然场景中人类视觉的去相关色彩空间中计算[CRE 15]。

image37image38

图5.30. a)自然色镶嵌,没有协调(132 RapidEye场景, 5 m空间分辨率); b)具有协调的自然彩色镶嵌。有关此图的彩色版本, 请参阅 www.iste.co.uk/baghdadi/qgis1.zip

在下面的练习中,我们从3个SPOT 6图像(“ROI_2_XS_Ortho_GSD2016”、 “ROI_3_XS_Ortho_GSD2015”和“ROI_4_XS_ Ortho_GSD2017”)生成输出镶嵌图案。

步骤

指令

显示帮助

不带参数地执行otbSYS_Mosaic以显示帮助:

otbcli_Mosaic

然后打印应用程序描述:

This is the Mosaic application, version 5.2.1 Perform mosaicking of
input images

Complete documentation: http://www.orfeotoolbox.org/Applications/Mosaic.html

Parameters:

-progress <boolean>

Report progress

MISSING -il <string list> Input Images (mandatory)

-vdcut <string list> Input VectorDatas for composition (optional, off by
default)

-vdstats <string list> Input

VectorDatas for statistics (optional, off by default)

-comp.feather <string> Feathering method [none/large/slim] (mandatory,
default value is none)

-comp.feather.slim.exponent <float> Transition smoothness (Unitary
exponent = linear transition) (optional, on by default, default value is
1)

-comp.feather.slim.lenght <float> Transition length (In cartographic
units) (optional, off by default)

-harmo.method <string> harmonization method [none/band/rgb] (mandatory,
default value is none)

-harmo.cost <string> harmonization cost function [rmse/musig/mu]
(mandatory, default value is rmse)

MISSING -out <string> [pixel]

Output image

[pixel=uint8/uint16/int16/uint32/int32/float/double] (default value is
float) (mandatory)

-interpolator <string>

Interpolation [nn/bco/linear] (optional, off by default,

default value is nn)

-interpolator.bco.radius <int32> Radius for bicubic interpolation
(mandatory, default value is 2)

-output.spacing <float> Pixel

Size (optional, off by default)

-tmpdir <string> Directory where to write temporary files (optional, off
by default)

-distancemap.sr <float> Distance map images sampling ratio (mandatory,
default value is 10)

-ram <int32> Available RAM (Mb) (optional, off by default, default value
is 128)

-inxml <string> Load otb application from xml file (optional, off by
default)

创建镶嵌

使用“-il”参数选择输入图像。

保留其他默认参数,输入输出文件名,并验证命令行:

otbcli_Mosaic-il PATH/TO/ROI_2…/IMG_SPOT6_..._1_R1C1.TIF
/PATH/TO/ROI_3…/IMG_SPOT6_..._1_R1C1.TIF
/PATH/TO/ROI_4…/IMG_SPOT6_..._1_R1C1.TIF -out mosaic.tif

将生成的镶嵌图片导入QGIS并可视化。

创建自然色彩和谐的镶嵌

使用“-il”参数选择输入图像。

将参数“-harmo.method”设置为“rgb”以协调镶嵌非自然色彩。 保留其他默认设置,输入输出图像的文件名,并验证命令行:

otbcli_Mosaic-il PATH/TO/ROI_2…/IMG_SPOT6_..._1_R1C1.TIF
/PATH/TO/ROI_3…/IMG_SPOT6_..._1_R1C1.TIF
/PATH/TO/ROI_4…/IMG_SPOT6_..._1_R1C1.TIF
-harmo.method rgb -out mosaic_harmo_rgb.tif

将生成的镶嵌图片导入QGIS并可视化。与之前生成的镶嵌进行比较。

表5.24.使用“Mosaic”应用程序镶嵌图像