本文档详细介绍了如何使用GDAL API进行栅格数据的打开、读取、元数据获取、波段操作、创建和复制等关键步骤，包括数据集信息的获取、栅格数据类型处理和文件格式支持情况。
**
栅格API教程
**
1. 打开文件
在打开GDAL支持的栅格数据存储之前，必须注册驱动程序。每个支持的格式都有一个驱动程序。通常这是通过 GDALAllRegister() 函数，该函数尝试注册所有已知的驱动程序，包括使用 GDALDriverManager::AutoLoadDrivers() . 注册，申请应调用独立 GDALOpen() 函数打开数据集，传递数据集的名称和所需的访问权限（GA_ReadOnly只读或GA_Update）。
#include "gdal_priv.h"
#include "cpl_conv.h" // for CPLMalloc()
int main()
{
GDALDataset *poDataset;
GDALAllRegister();
poDataset = (GDALDataset *) GDALOpen( pszFilename, GA_ReadOnly );
if( poDataset == NULL )
{
...;
}
注意如果 GDALOpen() 返回NULL表示打开失败，并且已经通过 CPLError() . 请注意pszFilename实际上不需要是物理文件的名称（尽管它通常是）。它的解释依赖于驱动程序，它可能是一个URL，一个文件名，在结尾添加了额外的参数来控制open或几乎任何东西。请不要将GDAL文件选择对话框限制为仅选择物理文件。
测试代码：
//1.打开文件
gdaldataset* podataset;
gdalallregister(); //注册驱动
podataset = (gdaldataset*)gdalopen("d:/class/chendu201055.25.tif", ga_readonly);//打开数据集，打开方式ga_readonly或ga_update
//判断是否为空
if (podataset == null)
{
cout << "无" << endl;
}
else
{
cout << "succ" << endl;
}
2. 获取数据集（Dataset）信息
GDALDataset 包含栅格标注栏列表，所有标注栏都属于同一区域，并且具有相同的分辨率。它还具有元数据、坐标系、地理参照转换、栅格大小和各种其他信息。在没有任何旋转或剪切的“北向上”图像的特殊但常见的情况下，使用地理参考变换采用以下形式：
adfGeoTransform[0] /* 左上角x */
adfGeoTransform[1] /* 东西方向一个像素对应的距离 */
adfGeoTransform[2] /* 旋转，0代表上面为北方 */
adfGeoTransform[3] /* 左上角y*/
adfGeoTransform[4] /* 旋转，0代表上面为北方*/
adfGeoTransform[5] /*南北方向一个像素对应的距离 */
如果要打印有关数据集的一些常规信息，可以执行以下操作：
double adfGeoTransform[6];
printf( "Driver: %s/%s\n",
poDataset->GetDriver()->GetDescription(),
poDataset->GetDriver()->GetMetadataItem( GDAL_DMD_LONGNAME ) );
printf( "Size is %dx%dx%d\n",
poDataset->GetRasterXSize(), poDataset->GetRasterYSize(),
poDataset->GetRasterCount() );
if( poDataset->GetProjectionRef() != NULL )
printf( "Projection is `%s'\n", poDataset->GetProjectionRef() );
if( poDataset->GetGeoTransform( adfGeoTransform ) == CE_None )
{
printf( "Origin = (%.6f,%.6f)\n",
adfGeoTransform[0], adfGeoTransform[3] );
printf( "Pixel Size = (%.6f,%.6f)\n",
adfGeoTransform[1], adfGeoTransform[5] );
}
测试代码：
//2、获取Dataset信息.输出
double addGeoTransform[6];
printf("Driver:%s/%s\n", poDataset->GetDriver()->GetDescription(),
poDataset->GetDriver()->GetMetadataItem(GDAL_DMD_LONGNAME));
//输出图像大小和波段个数
printf("大小:%dx%dx%d\n", poDataset->GetRasterXSize(), poDataset->GetRasterYSize(),
poDataset->GetRasterCount());
//输出图片的投影信息
if (poDataset->GetProjectionRef() != NULL)
{
printf("投影信息`%s'\n", poDataset->GetProjectionRef());
}
//输出图像的坐标和分辨率信息
if (poDataset->GetGeoTransform(addGeoTransform) == CE_None)
{
printf("原始坐标：%.6f,%.6f\n", addGeoTransform[0], addGeoTransform[3]);
printf("像素大小： %.6f,%.6f\n", addGeoTransform[1], addGeoTransform[5]);
}
3. 获取栅格波段
通过GDAL对栅格数据的访问是一次一个波段完成的。此外，还提供了元数据、块大小、颜色表和各种其他信息，这些信息都是按波段提供的。以下代码获取 GDALRasterBand 数据集中的对象（编号1到 GDALRasterBand::GetRasterCount() )并显示一些关于它的信息。
GDALRasterBand *poBand;
int nBlockXSize, nBlockYSize;
int bGotMin, bGotMax;
double adfMinMax[2];
poBand = poDataset->GetRasterBand( 1 );
poBand->GetBlockSize( &nBlockXSize, &nBlockYSize );
printf( "Block=%dx%d Type=%s, ColorInterp=%s\n",
nBlockXSize, nBlockYSize,
GDALGetDataTypeName(poBand->GetRasterDataType()),
GDALGetColorInterpretationName(
poBand->GetColorInterpretation()) );
adfMinMax[0] = poBand->GetMinimum( &bGotMin );
adfMinMax[1] = poBand->GetMaximum( &bGotMax );
if( ! (bGotMin && bGotMax) )
GDALComputeRasterMinMax((GDALRasterBandH)poBand, TRUE, adfMinMax);
printf( "Min=%.3fd, Max=%.3f\n", adfMinMax[0], adfMinMax[1] );
if( poBand->GetOverviewCount() > 0 )
printf( "Band has %d overviews.\n", poBand->GetOverviewCount() );
if( poBand->GetColorTable() != NULL )
printf( "Band has a color table with %d entries.\n",
poBand->GetColorTable()->GetColorEntryCount() );
测试代码：
//3、获取一个光栅波段 包含元数据、块大小、颜色表
//GDALRasterBand对象 读取波段
// GDALRasterBand* poBand;
//获取第一个波段
GDALRasterBand* poBand = poDataset->GetRasterBand(1);
//获取该波段的最大值最小值，获取失败，进行统计
int nBlockXSize, nBlockYSize;
int bGotMin, bGotMax;
double adfMinMax[2];
poBand = poDataset->GetRasterBand(1);
poBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
printf("Block=%dx%d Type=%s, ColorInterp=%s\n",
nBlockXSize, nBlockYSize,
GDALGetDataTypeName(poBand->GetRasterDataType()),
GDALGetColorInterpretationName(
poBand->GetColorInterpretation()));
adfMinMax[0] = poBand->GetMinimum(&bGotMin);
adfMinMax[1] = poBand->GetMaximum(&bGotMax);
if (!(bGotMin && bGotMax))
GDALComputeRasterMinMax((GDALRasterBandH)poBand, TRUE, adfMinMax);
printf("Min=%.3fd, Max=%.3f\n", adfMinMax[0], adfMinMax[1]);
if (poBand->GetOverviewCount() > 0)
printf("Band has %d overviews.\n", poBand->GetOverviewCount());
if (poBand->GetColorTable() != NULL)
printf("Band has a color table with %d entries.\n",
poBand->GetColorTable()->GetColorEntryCount());
4. 读取栅格数据
有几种读取栅格数据的方法，但最常见的是通过 GDALRasterBand::RasterIO() 方法。此方法将自动处理数据类型转换、上/下采样和窗口。下面的代码将把第一个扫描行的数据读入一个类似大小的缓冲区，作为操作的一部分将其转换为浮点。
float *pafScanline;
int nXSize = poBand->GetXSize();
pafScanline = (float *) CPLMalloc(sizeof(float)*nXSize);
poBand->RasterIO( GF_Read, 0, 0, nXSize, 1,
pafScanline, nXSize, 1, GDT_Float32,
0, 0 );
当pafScanline缓冲区不再使用时，应使用CPLFree（）释放它。
请注意，基于eRWFlag的设置（GF-read或GF-write），相同的RasterIO（）调用用于读取或写入。nXOff、nYOff、nXSize、nYSize参数描述要读取（或写入）的磁盘上栅格数据的窗口。它不必落在平铺边界上，尽管如果这样做，访问可能会更高效。
pData是数据读入或写入的内存缓冲区。它的真正类型必须是作为eBufType传递的任何类型，例如GDT_Float32或GDT_Byte。RasterIO（）调用将负责在缓冲区的数据类型和频带的数据类型之间进行转换。请注意，当将浮点数据转换为整数RasterIO（）时，向下舍入；当将源值转换到输出的合法范围之外时，将使用最接近的合法值。这意味着，例如，读取到GDT_字节缓冲区中的16位数据将映射所有大于255到255的值，数据不会缩放！
nBufXSize和nBufYSize值描述缓冲区的大小。当以全分辨率加载数据时，这将与窗口大小相同。但是，要加载分辨率降低的概述，可以将其设置为小于磁盘上的窗口。在这种情况下，如果概览合适，RasterIO（）将利用概览来更有效地执行IO。
nPixelSpace和nLineSpace通常为零，表示应使用默认值。然而，它们可用于控制对存储器数据缓冲器的访问，例如允许读入包含其它像素交错数据的缓冲器。
5.关闭数据集
//关闭文件
GDALClose((GDALDatasetH)poDataset);
6.创建文件
如果格式驱动程序支持创建，则可以创建GDAL支持格式的新文件。使用CreateCopy（）和Create（）创建文件有两种通用技术。CreateCopy方法涉及调用格式驱动程序上的CreateCopy（）方法，并传入应复制的源数据集。Create方法涉及在驱动程序上调用Create（）方法，然后使用单独的调用显式地写入所有元数据和栅格数据。所有支持创建新文件的驱动程序都支持CreateCopy（）方法，但只有少数驱动程序支持Create（）方法。
要确定特定格式是否支持Create或CreateCopy，可以检查格式驱动程序对象上的DCAP_Create和DCAP_CreateCopy元数据。确保 GDALAllRegister() 在打电话之前 GDALDriverManager::GetDriverByName() . 在本例中，我们获取一个驱动程序，并确定它是否支持Create（）和/或CreateCopy（）。
判断如下：
#include "cpl_string.h"
...
const char *pszFormat = "GTiff";
GDALDriver *poDriver;
char **papszMetadata;
poDriver = GetGDALDriverManager()->GetDriverByName(pszFormat);
if( poDriver == NULL )
exit( 1 );
papszMetadata = poDriver->GetMetadata();
if( CSLFetchBoolean( papszMetadata, GDAL_DCAP_CREATE, FALSE ) )
printf( "Driver %s supports Create() method.\n", pszFormat );
if( CSLFetchBoolean( papszMetadata, GDAL_DCAP_CREATECOPY, FALSE ) )
printf( "Driver %s supports CreateCopy() method.\n", pszFormat );
注意，许多驱动程序是只读的，不支持Create（）或CreateCopy（）。
7.使用CreateCopy（）
这个 GDALDriver::CreateCopy() 方法的使用相当简单，因为大多数信息都是从源数据集中收集的。但是，它包括用于传递特定于格式的创建选项的选项，以及用于在进行长数据集复制时向用户报告进度的选项。从一个名为pszSrcFilename的文件到一个名为pszDstFilename的新文件的简单复制，使用先前获取驱动程序的格式上的默认选项，可能如下所示：
GDALDataset *poSrcDS =
(GDALDataset *) GDALOpen( pszSrcFilename, GA_ReadOnly );
GDALDataset *poDstDS;
poDstDS = poDriver->CreateCopy( pszDstFilename, poSrcDS, FALSE,
NULL, NULL, NULL );
/* Once we're done, close properly the dataset */
if( poDstDS != NULL )
GDALClose( (GDALDatasetH) poDstDS );
GDALClose( (GDALDatasetH) poSrcDS );
注意，CreateCopy（）方法返回一个可写的数据集，必须正确关闭它才能完成数据集的写入并将其刷新到磁盘。在Python的情况下，当“dst_ds”超出范围时会自动发生这种情况。在CreateCopy（）调用中的目标文件名之后，bStrict选项使用的FALSE（或0）值指示，即使无法创建与输入数据集完全匹配的目标数据集，CreateCopy（）调用也应继续进行，而不会出现致命错误。这可能是因为输出格式不支持输入数据集的像素数据类型，或者是因为目标无法支持例如写入地理参考。
更复杂的情况可能涉及传递创建选项和使用预定义的进度监视器，如下所示：
#include "cpl_string.h"
...
char **papszOptions = NULL;
papszOptions = CSLSetNameValue( papszOptions, "TILED", "YES" );
papszOptions = CSLSetNameValue( papszOptions, "COMPRESS", "PACKBITS" );
poDstDS = poDriver->CreateCopy( pszDstFilename, poSrcDS, FALSE,
papszOptions, GDALTermProgress, NULL );
/* Once we're done, close properly the dataset */
if( poDstDS != NULL )
GDALClose( (GDALDatasetH) poDstDS );
CSLDestroy( papszOptions );
测试：
GDALDataset* poSrcDS =(GDALDataset*)GDALOpen("D:/class/chendu201055.25.tif", GA_ReadOnly);
GDALDataset* poDstDS;
poDstDS = poDriver->CreateCopy("D:\SJ", poSrcDS, FALSE, NULL, NULL, NULL);
Once we're done, close properly the dataset
if (poDstDS != NULL)
GDALClose((GDALDatasetH)poDstDS);
GDALClose((GDALDatasetH)poSrcDS)*/;
char** papszOptions = NULL;
papszOptions = CSLSetNameValue(papszOptions, "TILED", "YES");
papszOptions = CSLSetNameValue(papszOptions, "COMPRESS", "PACKBITS");
poDstDS = poDriver->CreateCopy("D:\SJ\jjsagd.tif", poSrcDS, FALSE, papszOptions, GDALTermProgress, NULL);
if (poDstDS != NULL) {
cout << "success" << endl;
printf("大小:%dx%dx%d\n", poDstDS->GetRasterXSize(), poDstDS->GetRasterYSize(),poDstDS->GetRasterCount());
printf("Driver:%s/%s\n", poDstDS->GetDriver()->GetDescription(),poDstDS->GetDriver()->GetMetadataItem(GDAL_DMD_LONGNAME));
GDALClose((GDALDatasetH)poDstDS);
}
GDALClose((GDALDatasetH)poSrcDS);
8.使用Create（）
对于不只是将现有文件导出到新文件的情况，通常需要使用 GDALDriver::Create() 方法（尽管通过使用虚拟文件或内存文件可以使用一些有趣的选项）。Create（）方法采用与CreateCopy（）类似的选项列表，但必须显式提供图像大小、带区数和带区类型。
GDALDataset *poDstDS;
char **papszOptions = NULL;
poDstDS = poDriver->Create( pszDstFilename, 512, 512, 1, GDT_Byte,
papszOptions );
成功创建数据集后，必须将所有适当的元数据和栅格数据写入文件。
double adfGeoTransform[6] = { 444720, 30, 0, 3751320, 0, -30 };
OGRSpatialReference oSRS;
char *pszSRS_WKT = NULL;
GDALRasterBand *poBand;
GByte abyRaster[512*512];
poDstDS->SetGeoTransform( adfGeoTransform );
oSRS.SetUTM( 11, TRUE );
oSRS.SetWellKnownGeogCS( "NAD27" );
oSRS.exportToWkt( &pszSRS_WKT );
poDstDS->SetProjection( pszSRS_WKT );
CPLFree( pszSRS_WKT );
poBand = poDstDS->GetRasterBand(1);
poBand->RasterIO( GF_Write, 0, 0, 512, 512,
abyRaster, 512, 512, GDT_Byte, 0, 0 );
/* Once we're done, close properly the dataset */
GDALClose( (GDALDatasetH) poDstDS );
测试代码：
GDALDataset* poDstDS;
char** papszOptions = NULL;
poDstDS = poDriver->Create("D:/SJ/asgdas.tif", 512, 512, 1, GDT_Byte,papszOptions);
double adfGeoTransform[6] = { 444720, 30, 0, 3751320, 0, -30 };
OGRSpatialReference oSRS;
char* pszSRS_WKT = NULL;
// GDALRasterBand* poBand;
GByte abyRaster[512 * 512];
poDstDS->SetGeoTransform(adfGeoTransform);
oSRS.SetUTM(11, TRUE);
oSRS.SetWellKnownGeogCS("NAD27");
oSRS.exportToWkt(&pszSRS_WKT);
poDstDS->SetProjection(pszSRS_WKT);
CPLFree(pszSRS_WKT);
poBand = poDstDS->GetRasterBand(1);
poBand->RasterIO(GF_Write, 0, 0, 512, 512,
abyRaster, 512, 512, GDT_Byte, 0, 0);
/* Once we're done, close properly the dataset */
GDALClose((GDALDatasetH)poDstDS);
}
注意：按自己需求选取代码，有些代码是测试功能的，需要自我检索
扫一扫
313
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
