说明:本文章转载自:http://zhan.renren.com/chinalee?gid=3602888497997597705&checked=true#nogo

GDAL栅格图像操作

GDAL是一个操作各种栅格和矢量(由ogr这个库实现)地理数据格式的开源库。包括读取、写入、转换、处理各种栅格和矢量数据格式(有些特定的格式对一些操作如写入等不支持)。

即使不是进行地理遥感方面的应用研究,GDAL也是一个非常有用的库,因为它可以支持大量我们常见的图像数据,比如jpg,gif之类的。完整的格式清单可以到此链接查看http://www.gdal.org/formats_list.html。而且已经有包括GoogleEarth在内的很多软件都是在使用GDAL作为后台库。

本文就以VC为开发平台介绍GDAL对栅格数据的操作方法。

Include目录是开发中需要的头文件,lib中是所需要的lib文件,在VC8中应当将其存放目录添加到目录列表中,选择菜单的“工具-选项-项目和解决方案-VC++目录”,分别在“包含文件”和“库文件”中将此两个目录添加进去。

在项目的属性页中,选择“配置属性-链接器-输入”,在“附加依赖项”中添加gdal_i-vc8.lib和gdal_id-vc8.lib两个使用GDAL中需要的静态库文件,或者在程序中添加以下两行代码也可以。

#pragma comment(lib, "gdal_i-vc8.lib")

#pragma comment(lib, "gdal_id-vc8.lib")

Bin目录下的动态链接库文件应当放置于程序能够访问的位置,比如windows\system32中。

此外,在程序中需要引入的头文件是gdal_priv.h。

现在开始用C++来对图像文件进行操作。

在打开文件之前需要首先注册所需要的驱动程序,一般来说我们可以默认注册所有支持的格式驱动,所使用的函数是GDALAllRegister()。然后就是打开文件操作。

这里要说一个数据集的概念,也就是所谓的Dataset。在GDAL中可以说数据的核心就是Dataset,简单来说可以将Dataset就理解为图像文件,比如说一个jpeg格式的文件就是一个数据集,当然其他一些文件格式可能在一个数据集中包含多于一个文件,比如可能除了图像数据文件外还可能会有一些附加信息文件等。

在数据集下最重要组成部分就是所谓的波段band,波段可多可少,比如一个RGB真彩色的图像就有3个波段,分别代表红色绿色和蓝色波段,如果是灰度图,那可能就只有一个波段,而很多遥感图像可能就会多于3个波段。

除了波段外,数据集中还含有图像相关的坐标系投影信息,元数据信息等数据。

文件的打开使用的是GDALOpen ( const char *  pszFilename, GDALAccess  eAccess  ),pszFilename是文件路径,eAccess 是访问权限,可以是GA_ReadOnly只读,也可以是GA_Update来对文件进行修改。比如我们以只读模式打开一个tif文件:

GDALDataset *poDataset; //数据集对象指针

GDALAllRegister();//注册驱动

poDataset = (GDALDataset *) GDALOpen( "c:\\terra335h_EV_250_Aggr500_RefSB_b0.tif", GA_ReadOnly );

if( poDataset != NULL /*检查是否正常打开文件*/)

{

//do something

}

delete poDataset; //释放资源

在确认poDataset不是NULL的情况下就可以对图像数据集进行操作了。

首先来看一下有关这个图像的基本信息,比如长宽和波段数。

cout<<"RasterXSize:"<<poDataset->GetRasterXSize()<<endl;//x方向长度

cout<<"RasterYSize:"<<poDataset->GetRasterYSize()<<endl;//y方向长度

cout<<"RasterCount:"<<poDataset->GetRasterCount()<<endl;//波段数量

我的这个文件输出结果如下:

RasterXSize:5684

RasterYSize:3655

RasterCount:1

这是个5684×3655的灰度图。

有关的投影信息可以由GetProjectionRef函数获得。

cout<<poDataset->GetProjectionRef()<<endl;

输出一个字符串,是WKT格式的坐标投影信息,形如:

GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.2572235630016,AUT

HORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["deg

ree",0.0174532925199433],AUTHORITY["EPSG","4326"]]

有关仿射变换和GCP的信息可以通过GetGeoTransform,GetGCPCount,GetGCPProjection,GetGCPs获得。比如取得仿射变换的信息如下:

double adfGeoTransform[6];//仿射参数

poDataset->GetGeoTransform(adfGeoTransform);//获取参数

然后就可以使用以下公式来计算仿射坐标了:

Xgeo = adfGeoTransform[0] + Xpixel* adfGeoTransform[1] + Yline* adfGeoTransform[2]

Ygeo = adfGeoTransform[3] + Xpixel* adfGeoTransform[4] + Yline* adfGeoTransform[5]

这里的坐标系是以左上角为(0,0)点,向右向下为正向。

Xgeo表示转换后的地理x坐标,Ygeo是其相应y坐标,Xpixel是图像中像素的列坐标,Yline则是行坐标。而且Xpixel和Yline的(0,0)坐标是图像左上角像素的左上角,因此左上角像素的中心坐标是Xpixel=0.5,Yline=0.5。

接下来,我们就开始进入到波段处理。波段的获取使用GetRasterBand函数,参数是波段序号,从1开始。比如还是上面那个图,只有一个波段,那要得到这个波段的操作如下:

GDALRasterBand *poBand;

poBand = poDataset->GetRasterBand( 1 );

关于波段的信息也很多,我们捡几个主要的看看。

首先也是尺寸,在一个Dataset里所有波段的尺寸都一样,也就是上面用Dataset获取的数值,就不多数了。

然后是数据类型。不同的图像格式的存储的数据类型是不一样的,比如我的这个tif用的是UInt16,也就是无符号的短整型,其他的用的可能就是byte型或其他类型的了。获取波段数据的数据类型信息可以如下进行:

cout<<"Data Type:"<<poBand->GetRasterDataType()<<endl;

这里返回的信息是一个数据类型的编号,对应枚举类型GDALDataType中的一个值,在此tif的数值是2,对应GDT_UInt16,如果希望返回可读性比较强的信息,可以如下

cout<<"Data Type:"<<GDALGetDataTypeName(poBand->GetRasterDataType())<<endl;

返回的是数据类型名的字符串:Data Type:UInt16

关于波段色彩类型的情况类似数据类型,也是可以由GetColorInterpretation获取关于类型分类的GDALColorInterp枚举值,同时由GDALGetColorInterpretationName来返回名称字符串:

cout<<"Color interpretation :"<<poBand->GetColorInterpretation()<<endl;

cout<<"Color interpretation:"<<GDALGetColorInterpretationName(poBand->GetColorInterpretation())<<endl;

输出:

Color interpretation :1

Color interpretation:Gray

表明这个波段是个灰度图波段。如果是真彩色图可能得到的结果就是Red,Blue或Green了。

GetMaximum和GetMinimum分别返回的是波段数据值可能的最大值和最小值,比如UInt16类型的波段最大值是65535,最小值就是0,知道这两个值就可以将图像转换成0~255的数值用于图像的显示操作了。如:

cout<<"Max :"<<poBand->GetMaximum()<<endl;

cout<<"Min :"<<poBand->GetMinimum()<<endl;

终于要说到具体数据的读写了。关于波段数据的续写核心函数就是RasterIO。这个函数可以将图像的某一个子块读入或写入。

CPLErr GDALRasterBand::RasterIO ( GDALRWFlag  eRWFlag, int  nXOff, int  nYOff, int  nXSize, int nYSize, void *  pData, int  nBufXSize, int  nBufYSize, GDALDataType  eBufType, int  nPixelSpace, int  nLineSpace  )

eRWFlag:读写标志位,GF_Read读,GF_Write写。

nXOff,nYOff:读写数据块的左上角坐标

nXSize,nYSize:读写数据块的xy方向像素尺寸

pData:读写数据块缓冲区指针

nBufXSize,nBufYSize:缓冲区的像素宽度和高度

eBufType:数据类型,如GDT_UInt16等

nPixelSpace:pData中一个像素的字节大小,值为0的时候就默认是sizeof(eBufType)

nLineSpace :pData中一行数据的字节大小,值为0的时候就默认是sizeof(eBufType)*nBufXSize

比如说我们要读取第一行数据:

short *pafScanline;

int   nXSize = poBand->GetXSize();

pafScanline = ( short *) CPLMalloc(sizeof( short)*nXSize); //分配缓冲区空间

poBand->RasterIO( GF_Read, 0, 0, nXSize, 1,

pafScanline, nXSize, 1, GDT_UInt16, 0, 0 );

CPLFree(pafScanline); //释放缓冲区空间

获取波段的像素数据后就可以使用GDI或GL的像素操作函数绘制图像了。另外,如果希望一次性读出整个波段的所有数据,根据GDAL文档的说法使用ReadBlock函数效率会更好一些。

至于要对有关数据和信息进行写入,一般都有读取方法对应的Set*或Write*函数可以使用,在此不再赘述。

【C/C++】超大遥感影像读取和存储 GDAL相关推荐

  1. 覆盖5大任务,30+特色模型,高性能、全流程开发套件PaddleRS助力遥感影像智能解译化繁为简...

    近年来,随着卫星技术的发展和深度学习的火热,基于深度学习的遥感影像智能解译得到了前所未有的关注,并已成功应用于建筑物变化检测.SAR影像船舶检测.道路提取.多光谱影像分类等任务中.高精度.高速度.自动 ...

  2. 8影像计算ndvi landsat_使用GDAL读取遥感影像的信息

    读取影像数据集的元数据 GDAL已经提供了足够方便的函数,可以读取影像的一些元数据信息, 从而方便对数据进行处理.GDAL一般是以字典的形式对元数据进行组织的, 但是对于不同的栅格数据类型,元数据的类 ...

  3. python 读取geotiff_科学网—利用python GDAL库读写geotiff格式的遥感影像方法 - 张伟的博文...

    (1)利用python GDAL库读写geotiff格式的遥感影像方法,具有很好的参考价值,不错! from osgeo import gdal import numpy as np def read ...

  4. python/gdal处理遥感影像(读取、投影转换、裁剪、建立图像金字塔等)

    python/gdal处理遥感影像(读取.投影转换.裁剪.建立图像金字塔等) gdal库简单介绍 python使用gdal 一.安装python环境 二.安装gdal库 三.使用gdal处理遥感影像 ...

  5. 监督分类:SVM即支持向量机实现遥感影像监督分类(更新:添加机器学习模型存储、大影像划框拼接)

    前面已经有一个版本了,但是影像太大内存顶不住,而且训练和预测没有分离,后面批量用这个不可能每次每张影像都训练了再预测,这次正好有需求,我就最后把这个整理一下,算是终版吧,以后也不会再花时间整这个了 这 ...

  6. python读取tiff影像_科学网—利用python GDAL库读写geotiff格式的遥感影像方法 - 张伟的博文...

    (1)利用python GDAL库读写geotiff格式的遥感影像方法,具有很好的参考价值,不错! from osgeo import gdal import numpy as np def read ...

  7. 3_读取遥感影像信息

    读取遥感影像信息 本文中的所有方法可以读取遥感影像信息,但是仅限于概括地获取元数据信息,无法具体到对像元进行处理. 打开已有的GeoTiff GDAL不能通过with打开影像,会报错 GDAL打开影像 ...

  8. COG云原生优化遥感影像,瓦片切分的应用实践

    摘要:云上遥感影像文件Cloud optimized GeoTIFF(COG)格式的详细介绍,大量数据上云面临的挑战,并分享了获得云原生影像最佳性能的实践经验. 本文分享自华为云社区<COG云原 ...

  9. TM遥感影像波段/通道bands

    遥感影像波段band 一.遥感影像波段 1.原理 2.举例说明 二.TM影像各波段简介 1.TM影像概述 2.各波段影像特征 3.波段组合 4.类型提取 5.光谱差异 三.遥感图像--多波段数据存储的 ...

最新文章

  1. k8s使用kube-router网络插件并监控流量状态
  2. poj-1384 Piggy-Bank
  3. 1-1 什么是微信小程序
  4. Progressive Web App(PWA)
  5. python3的float数精度_python浮点数精度问题
  6. linux gcc 简单使用记录01
  7. 翻译:程序员数据结构基础:选择正确的数据结构
  8. 浅谈Spring的AOP实现-代理机制
  9. 循环 计算数值的整数次幂
  10. python机器学习案例系列教程——集成学习(Bagging、Boosting、随机森林RF、AdaBoost、GBDT、xgboost)
  11. SpringCloud组件:Eureka服务注册是采用主机名还是IP地址?
  12. Direct3D 11 Devices 之 Using Direct3D 11 feature data to supplement Direct3D feature levels
  13. 最新狂雨小说CmsV1.5.2漂亮的小说网站源码
  14. [人工智能-深度学习-29]:卷积神经网络CNN - 全连接网络与卷积网络结构的互为等效与性能比较
  15. 35岁鞋不合脚的问题
  16. maven项目,pom.xml文件变成小虫子(蜘蛛)解决办法
  17. oppor15android10怎么降级,OPPOR15系统降级教程
  18. 【Matlab】数据插值
  19. CSS line-height概念与举例
  20. Linux Performance Analysis and Tools(Linux性能分析和工具)

热门文章

  1. mac word2016 去除页眉下面的横线
  2. 《PHP挖宝》1—再论框架
  3. 索尼计算机bios正确设置,索尼笔记本电脑怎么进入Bios,小编教你如何四步完成
  4. 移植AT91Bootstrap1.15
  5. Silverlight 2教程(四):Chiron.exe:Silverlight 2打包和动态语言部署工具
  6. 粒子能量、量子波动方程、狄拉克方程、量子态【量子力学基础知识学习笔记_3】
  7. matlab张正友摄像机标定算法应用,张正友摄像机标定的研究(MATLAB+OpenCV)
  8. matlab程序按哪里运行,脱离matlab运行可执行程序的步骤
  9. 红米4手机(其它小米应该一样)adb 调试(usb ,tcp)
  10. 泰拉瑞亚指令代码大全 无限钱无敌作弊码一览