Reading Geospatial Raster files with GDAL(使用GDAL读取地理空间栅格文件)

Reading Geospatial Raster files with GDAL(使用GDAL读取地理空间栅格文件)

参考教程及hgt,tif等图片资源下载链接:
https://download.csdn.net/download/weixin_43940163/11188632

一.先决条件

1.Digital Elevation Model:

​ 数字高程模型（Digital Elevation Model)，简称DEM，是通过有限的地形高程数据实现对地面地形的数字化模拟（即地形表面形态的数字化表达），它是用一组有序数值阵列形式表示地面高程的一种实体地面模型，是数字地形模型(Digital Terrain Model，简称DTM)的一个分支，其它各种地形特征值均可由此派生。

参考网址:https://baike.baidu.com/item/数字高程模型

2.GeoTiff:

​ 通俗来讲，GeoTIFF 是 TIFF 格式的一种扩展，是支持地理空间参考的 TIFF. GeoTIFF 已经是 TIFF 规范的一部分。

​ 标签图像文件格式（Tag Image File Format，简写为TIFF）是一种灵活的位图格式，主要用来存储包括照片和艺术图在内的图像。它最初由Aldus公司与微软公司一起为PostScript打印开发。TIFF与JPEG和PNG一起成为流行的高位彩色图像格式。

参考网址:https://baike.baidu.com/item/TIFF/2106?fr=aladdin

3.C++ std::pair:

​ std::pair主要的作用是将两个数据组合成一个数据，两个数据可以是同一类型或者不同类型。例如std::pair<int,float> 或者 std：：pair<double,double>等。pair实质上是一个结构体，其主要的两个成员变量是first和second，这两个变量可以直接使用。初始化一个pair可以使用构造函数

4.Lerp

​ Lerp,就是返回两个值之间的插值，一般有三个参数。第一个参数为初始值，第二个参数为最终值，插值为0~1d的一个浮点数值，为0时为初始值，1时为最终值，为0到1之间的数值时返回一个混合数值。若第三个参数为Time.daltaTime,这返回随时间推移的一个值。

二.代码调试

程序的输入为一张影像图和一张DEM图。在读入图像时，要想例程中一样，在imread()函数中加入flag标签，表明我们要使用GDAL取读取影响。

两个数据读取函数如下：

cv::Mat image = cv::imread(argv[1], cv::IMREAD_LOAD_GDAL | cv::IMREAD_COLOR );

cv::Mat dem = cv::imread(argv[2], cv::IMREAD_LOAD_GDAL | cv::IMREAD_ANYDEPTH );

其中argv[1]为gdal_input.jpg,argv[2]为N37W123.tif(修改过)

转换函数:

cv::Point2d world2dem( const cv::Point2d&, const cv::Size&);

该函数完成了由影像图的像素坐标到世界坐标的转换，利用的是双线性插值算法。

cv::Point2d pixel2world( const int&, const int&, const cv::Size& );

该函数完成了由世界坐标向DEM中像素坐标系的转换。

其余的函数就是完成了一个颜色映射表，依据高程信息，分级填充输出像素的BGR值。

调试这个代码时，有以下两点要注意：**

1、DEM数据N37W123.hgt，不能用imread()函数直接读取。需要你使用ArcGIS或ENVI等软件转换为*.tif格式数据才能读取。

2、 if( dem.type() != CV_16SC1 ){throw std::runtime_error(“DEM image type must be CV_16SC1”); }

一直报错,我在我的代码里注释掉了

下面是我的源码和图片源(包括hgt,tif,输入的jpg等)预览和链接

ps:以上是linux下的目录结构,windows用户直接运行程序即可,不过记得输入图片路径

三.代码详解

/*
 * gdal_image.cpp -- Load GIS data into OpenCV Containers using the Geospatial Data Abstraction Library
*/

// OpenCV Headers
#include "opencv2/core.hpp"
#include "opencv2/imgproc.hpp"
#include "opencv2/highgui.hpp"

// C++ Standard Libraries
#include <cmath>
#include <iostream>
#include <stdexcept>
#include <vector>

using namespace std;

/****函数声明****/
// define the corner points 定义角点
//    Note that GDAL library can natively determine this GDAL库可以自动判断出以下的值
cv::Point2d tl( -122.441017, 37.815664 );
cv::Point2d tr( -122.370919, 37.815311 );
cv::Point2d bl( -122.441533, 37.747167 );
cv::Point2d br( -122.3715,   37.746814 );

// determine dem corners
cv::Point2d dem_bl( -122.0, 38);
cv::Point2d dem_tr( -123.0, 37);

// range of the heat map colors
std::vector<std::pair<cv::Vec3b,double> > color_range;


// List of all function prototypes
cv::Point2d lerp( const cv::Point2d&, const cv::Point2d&, const double& );

cv::Vec3b get_dem_color( const double& );

cv::Point2d world2dem( const cv::Point2d&, const cv::Size&);

cv::Point2d pixel2world( const int&, const int&, const cv::Size& );

void add_color( cv::Vec3b& pix, const uchar& b