文章目录
GDAL 读取和保存 Grd 文件
绘图软件 Golden Surfer 的网格文件(Grd)主要有三种存储格式:
- Sufer 6 Text
- Sufer 6 Binary
- Sufer 7
其中,Sufer 6 Text 格式以 ASCII 码存储,可以直接用文本编辑查看。另外两种是以二进制格式存储的,需要按一定的规则进行读取。
具体的文件格式可以参考:Grd文件格式说明_地球屋里老师的博客-CSDN博客_.grd文件
GDAL是可以直接支持读写 Golden Surfer 的 grd 文件,对应的格式字符串如下:
要实现和其他图像格式的转换也非常的简单,直接在对应的驱动下 createcopy 即可。
这里主要介绍下 Sufer 6 Binary Grd 文件的读取。
Sufer 6 Binary 文件格式
文件样式:
DSBB
nx ny
xmin xmax
ymin ymax
zmin zmax
z11 z12 ... z1nx
z21 z22 ... z2nx
⁝ ⁝ ⁝
zny1 zny2 ... znynx
数据说明:
行号 | 变量名 | 数据类型 | 说明 |
---|---|---|---|
1 | id | unsigned char[4] × \times × 4 bytes |
格式标识符,‘DSBB’ |
2 | nx | unsigend short × \times × 2 bytes |
x 方向的点数(列数) |
ny | unsigend short × \times × 2 bytes |
y 方向的点数(行数) | |
3 | xmin | double × \times × 8 bytes |
x 地理坐标的最小值 |
xmax | double × \times × 8 bytes |
x 地理坐标的最大值 | |
4 | ymin | double × \times × 8 bytes |
y 地理坐标的最小值 |
ymax | double × \times × 8 bytes |
y 地理坐标的最大值 | |
5 | zmin | double × \times × 8 bytes |
数据 z 中的最小值(不含白化点 nodata) |
zmax | double × \times × 8 bytes |
数据 z 中的最大值(不含白化点 nodata) | |
6 | z(i,j) | float × \times × 4 bytes |
对应位置的数据 z (等于nodata表示该点无效) |
⋯ \cdots ⋯ |
读取 Grd 文件
Grd 文件和 GDAL tiff 文件的区别
- 坐标原点不同
Grd 文件的坐标原点在图像的西南角,而 tiff 的坐标原点在西北角。Grd 的行的顺序是从南向北的,tiff 的行顺序则是从北向南的。因此 grd 文件 存储的数据其实是按 从南向北,由西向东 的顺序存储。 - 数据的锚点的不同
Grd 的数据锚点是在像素块的中心位置,而 tiff 的数据锚点是在像素块的左上角(西北角)。Tiff 还有 tie point 的概念,分别在图像的西北角和东南角,这两个 tie point 界定了 tiff 图像的地理范围。
根据 grd 的数据计算 tiff 的信息:
- 分辨率的计算
Δ x = x max − x min N x − 1 Δ y = y min − y max N y − 1 \begin{split} \Delta x &= \frac{x_{\max} – x_{\min}}{N_x-1}\\ \Delta y &= \frac{y_{\min} – y_{\max}}{N_y-1}\\ \end{split} ΔxΔy=Nx−1xmax−xmin=Ny−1ymin−ymax
这里 y 方向分辨率为负值,和 GDAL 中的一致。
- Tie point 的计算
{ X W N = x min − 0.5 Δ x Y W N = y max − 0.5 Δ y \begin{cases} X_{_{WN}} = x_{\min} – 0.5 \Delta x\\ Y_{_{WN}} = y_{\max} – 0.5 \Delta y \end{cases} {
XWN=xmin−0.5ΔxYWN=ymax−0.5Δy
tiff 的仿射系数 geotrans[6]
就可以确定了,即
{ g e o t r a n s [ 0 ] = X W N g e o t r a n s [ 1 ] = Δ x g e o t r a n s [ 2 ] = 0.0 g e o t r a n s [ 3 ] = Y W N g e o t r a n s [ 4 ] = 0.0 g e o t r a n s [ 5 ] = Δ y \begin{cases} geotrans[0] = X_{_{WN}}\\ geotrans[1] = \Delta x\\ geotrans[2] = 0.0\\ geotrans[3] = Y_{_{WN}}\\ geotrans[4] = 0.0\\ geotrans[5] = \Delta y\\ \end{cases} ⎩
⎨
⎧geotrans[0]=XWNgeotrans[1]=Δxgeotrans[2]=0.0geotrans[3]=YWNgeotrans[4]=0.0geotrans[5]=Δy
手动实现读取 Grd 文件
根据 [[#Sufer 6 Binary 文件格式]] 定义一个抬头结构体:
typedef unsigned char byte;
typedef unsigned char DT_8U;
typedef unsigned short DT_16U;
typedef short DT_16S;
typedef unsigned int DT_32U;
typedef int DT_32S;
typedef float DT_32F;
typedef double DT_64F;
typedef struct {
DT_8U id[4];
DT_16U sizex;
DT_16U sizey;
DT_64F xmin;
DT_64F xmax;
DT_64F ymin;
DT_64F ymax;
DT_64F zmin;
DT_64F zmax;
} DSBBHeader;
使用 fread
直接读取即可,注意写入数据的时候,需要从最后一行开始写起。
参考代码如下:
GDALDataset* grd2gdal_manual(const char* pszSrcFile, const char* pszDstFile, const char* pFormat, DT_32F* pNodata = NULL)
{
FILE* fout = NULL;
// fout = fopen(pszSrcFile, "rb");
fopen_s(&fout, pszSrcFile, "rb");
if (fout == NULL) {
printf("无法打开文件: %s\n", pszSrcFile);
return NULL;
}
// 读取文件抬头信息
DSBBHeader header;
fread(&header, sizeof(DSBBHeader), 1, fout);
// 读取数据部分
DT_32F* pdata = new DT_32F[header.sizex * header.sizey];
fread(pdata, sizeof(DT_32F), header.sizex * header.sizey, fout);
// 创建输出文件
GDALAllRegister();
GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName(pFormat);
if (pDriver == NULL) {
printf("不支持格式: %s\n", pFormat);
fclose(fout);
return NULL;
}
GDALDataset* pDstDS = pDriver->Create(pszDstFile, header.sizex, header.sizey, 1, GDT_Float32, NULL);
if (pDstDS == NULL) {
printf("创建失败!\n");
fclose(fout);
return NULL;
}
// 设置nodata
if (pNodata) {
pDstDS->GetRasterBand(1)->SetNoDataValue(*pNodata);
}
// 计算仿射参数
double adftrans[6] = {
0.0 };
adftrans[1] = (header.xmax - header.xmin) / static_cast<double>(header.sizex - 1);
adftrans[5] = (header.ymin - header.ymax) / static_cast<double>(header.sizey - 1);
adftrans[0] = header.xmin - 0.5 * adftrans[0];
adftrans[3] = header.ymax - 0.5 * adftrans[5];
pDstDS->SetGeoTransform(adftrans);
// 从最后一行开始写入数据
for (int y = header.sizey - 1; y >= 0; --y) {
pDstDS->GetRasterBand(1)->RasterIO(GF_Write, 0, y, header.sizex, 1,
pdata + header.sizex * (header.sizey - 1 - y), header.sizex, 1, GDT_Float32, 0, 0, 0);
}
fclose(fout);
return pDstDS;
}
使用GDAL直接读取Grd文件
GDAL 可以直接读取 Grd 文件,并且会自动进行仿射系数的计算,自动进行数据读取的翻转。
GDALDataset* grd2gdal(const char *pszSrcFile, const char *pszDstFile, const char *pFormat, DT_32F *pNodata = NULL)
{
GDALAllRegister();
GDALDataset* pSrcDS = (GDALDataset*)GDALOpen(pszSrcFile, GA_ReadOnly);
if (pSrcDS == NULL) {
printf("无法打开: %s\n", pszSrcFile);
return NULL;
}
GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName(pFormat);
if (pDriver == NULL) {
printf("不支持格式: %s\n", pFormat);
GDALClose(pSrcDS);
return NULL;
}
GDALDataset* pDstDS = pDriver->CreateCopy(pszDstFile, pSrcDS, false, NULL, NULL, NULL);
if (pDstDS == NULL) {
printf("创建失败!\n");
GDALClose(pDstDS);
return NULL;
}
if (pNodata) {
pDstDS->GetRasterBand(1)->SetNoDataValue(*pNodata);
}
GDALClose(pSrcDS);
return pDstDS;
}
保存为 grd 文件
手动实现保存为 grd 文件,参考代码如下:
bool gdal2grd(GDALDataset *pSrcDS, const char *pGrdFile)
{
if (pSrcDS == NULL) {
printf("数据指针为空!\n");
return false;
}
DSBBHeader header;
int sizex = pSrcDS->GetRasterXSize();
int sizey = pSrcDS->GetRasterYSize();
double adftrans[6];
if (pSrcDS->GetGeoTransform(adftrans) != CE_None) {
printf("图像无仿射系数!\n");
GDALClose(pSrcDS);
return false;
}
GDALRasterBand* pband = pSrcDS->GetRasterBand(1);
if (pband->GetRasterDataType() != GDT_Float32) {
printf("图像数据类型不是 float!\n");
GDALClose(pSrcDS);
return false;
}
int hasNodata = 0;
DT_32F nodata = pband->GetNoDataValue(&hasNodata);
header.id[0] = 'D';
header.id[1] = 'S';
header.id[2] = 'B';
header.id[3] = 'B';
header.sizex = sizex;
header.sizey = sizey;
header.xmin = adftrans[0] + adftrans[1] / 2;
header.ymax = adftrans[3] + adftrans[5] / 2;
header.xmax = header.xmin + (static_cast<double>(sizex) - 1.0) * adftrans[1];
header.ymin = header.ymax + (static_cast<double>(sizey) - 1.0) * adftrans[5];
FILE* fout = NULL;
fopen_s(&fout, pGrdFile, "wb");
if (fout == NULL) {
printf("文件创建失败!\n");
GDALClose(pSrcDS);
return false;
}
fwrite(&header, sizeof(DSBBHeader), 1, fout);
DT_64F maxv = -1.0;
DT_64F minv = 99999.0;
DT_32F* pdata = new DT_32F[sizex];
for (int y = sizey - 1; y >= 0; --y) {
pband->RasterIO(GF_Read, 0, y, sizex, 1, pdata, sizex, 1, GDT_Float32, 0, 0, 0);
for (int x = 0; x < sizex; ++x) {
if (hasNodata && pdata[x] == nodata)
continue;
maxv = max(maxv, static_cast<DT_64F>(pdata[x]));
minv = min(minv, static_cast<DT_64F>(pdata[x]));
}
fwrite(pdata, sizeof(DT_32F), sizex, fout);
}
delete[] pdata;
header.zmin = static_cast<DT_64F>(minv);
header.zmax = static_cast<DT_64F>(maxv);
fseek(fout, 40, SEEK_SET);
fwrite(&header.zmin, sizeof(DT_64F), 1, fout);
fwrite(&header.zmax, sizeof(DT_64F), 1, fout);
PrintDSBBInfo(header);
fclose(fout);
return true;
}
测试
测试结果
输入 grd 图像:
可以看到 grd 文件的坐标原点在左下角(西南点)
原点(1,1)的地理坐标为 (113.312567, 33.917719)原图可在此处下载(里面还有一个等值线文件):
链接:https://pan.baidu.com/s/1MPuikSRIbUdi0S_pQXT0-A
提取码:fta4
保存为 tiff 图像:
输出结果:
--- 使用gdal读取grd ---
Size: 1485 x 1317
band: 1
scale: 0.000270 x -0.000270
Tie points(WN): (113.312432, 34.272779)
Tie points(ES): (113.712937, 33.917584)
Center Points(WN): (113.312567, 34.272644)
Center Points(ES): (113.712802, 33.917719)
Nodata: 99999.000
Data Type: GDT_Float32
Data range: 70.954514 ~ 3075.133301
--- Header Info ---
ID: DSBB
Size: 1485 x 1317
scale: 0.000270 x -0.000270
Center Points(WN): (113.312567, 34.272644)
Center Points(ES): (113.712802, 33.917719)
Tie points(WN): (113.312432, 34.272779)
Tie points(ES): (113.712937, 33.917584)
Data range: 70.954512 ~ 3075.133384
--- 手动读取grd ---
Size: 1485 x 1317
band: 1
scale: 0.000270 x -0.000270
Tie points(WN): (113.312567, 34.272779)
Tie points(ES): (113.713071, 33.917584)
Center Points(WN): (113.312702, 34.272644)
Center Points(ES): (113.712937, 33.917719)
Nodata: 99999.000
Data Type: GDT_Float32
Data range: 70.954514 ~ 3075.133301
--- 将gdal保存为grd ---
ID: DSBB
Size: 1485 x 1317
scale: 0.000270 x -0.000270
Center Points(WN): (113.312567, 34.272644)
Center Points(ES): (113.712802, 33.917719)
Tie points(WN): (113.312432, 34.272779)
Tie points(ES): (113.712937, 33.917584)
Data range: 70.954514 ~ 3075.133301
保存成功: C:/Users/q2799/Project/saved_res.grd
手动读取grd和gdal读取,得到的
tie point
及 分辨率都是一致的,Data range
由于float
的精度,有少许的不同(毕竟float
有效数字大概就7位)。
测试代码
main.cpp
#include <iostream>
#include "rwgrd.h"
using namespace std;
int main(int argc, char** argv)
{
cout << "Hello World" << endl;
const char* pszSrcFile = "C:/Users/q2799/Project/demData/dem overlap"
"/TH02-01AB_InSAR_20191223181856_0000003428_001_002_004_L1B.grd";
const char* pszDstFile = "C:/Users/q2799/Project/res.tiff";
DT_32F nodata = 99999.0;
GDALAllRegister();
GDALDataset* pSrcDS = grd2gdal(pszSrcFile, pszDstFile, "GTiff", &nodata);
printf("\n--- 使用gdal读取grd ---\n");
PrintDSInfo(pSrcDS);
const char* pszDstFile_manual = "C:/Users/q2799/Project/res_manual.tiff";
GDALDataset* pManDS = grd2gdal_manual(pszSrcFile, pszDstFile_manual, "GTiff", &nodata);
printf("\n--- 手动读取grd ---\n");
PrintDSInfo(pManDS);
const char * pszGrdFile = "C:/Users/q2799/Project/saved_res.grd";
printf("\n--- 将gdal保存为grd ---\n");
if (gdal2grd(pSrcDS, pszGrdFile)) {
printf("保存成功: %s\n", pszGrdFile);
}
else {
printf("保存失败!\n");
}
if (pSrcDS) {
GDALClose(pSrcDS);
}
if (pManDS) {
GDALClose(pManDS);
}
return 0;
}
rwgrd.h
#ifndef _RWGRD_H_
#define _RWGRD_H_
#include <gdal.h>
#include <gdal_priv.h>
typedef unsigned char byte;
typedef unsigned char DT_8U;
typedef unsigned short DT_16U;
typedef short DT_16S;
typedef unsigned int DT_32U;
typedef int DT_32S;
typedef float DT_32F;
typedef double DT_64F;
typedef struct {
DT_8U id[4];
DT_16U sizex;
DT_16U sizey;
DT_64F xmin;
DT_64F xmax;
DT_64F ymin;
DT_64F ymax;
DT_64F zmin;
DT_64F zmax;
} DSBBHeader;
GDALDataset* grd2gdal(const char* pszSrcFile, const char* pszDstFile, const char* pFormat, DT_32F* pNodata);
void PrintDSInfo(GDALDataset* pImgDS);
void PrintDSBBInfo(const DSBBHeader& header);
GDALDataset* grd2gdal_manual(const char* pszSrcFile, const char* pszDstFile, const char* pFormat, DT_32F* pNodata);
bool gdal2grd(GDALDataset* pSrcDS, const char* pGrdFile);
#endif
rwgrd.cpp
#include <cstdio>
#include <cmath>
#include "rwgrd.h"
using namespace std;
static const char* getDataTypeName(GDALDataType eDT)
{
switch (eDT)
{
case GDT_Byte: return "GDT_Byte";
case GDT_UInt16: return "GDT_UInt16";
case GDT_UInt32: return "GDT_UInt32";
case GDT_Int32: return "GDT_Int32";
case GDT_Float32: return "GDT_Float32";
case GDT_Float64: return "GDT_Float64";
case GDT_CInt16: return "GDT_CInt16";
case GDT_CInt32: return "GDT_CInt32";
case GDT_CFloat32: return "GDT_CFloat32";
case GDT_CFloat64: return "GDT_CFloat64";
default: return "GDT_Unknown";
}
}
// 使用 gdal 直接读取 grd 文件
GDALDataset* grd2gdal(const char *pszSrcFile, const char *pszDstFile, const char *pFormat, DT_32F *pNodata = NULL)
{
GDALAllRegister();
GDALDataset* pSrcDS = (GDALDataset*)GDALOpen(pszSrcFile, GA_ReadOnly);
if (pSrcDS == NULL) {
printf("无法打开: %s\n", pszSrcFile);
return NULL;
}
GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName(pFormat);
if (pDriver == NULL) {
printf("不支持格式: %s\n", pFormat);
GDALClose(pSrcDS);
return NULL;
}
GDALDataset* pDstDS = pDriver->CreateCopy(pszDstFile, pSrcDS, false, NULL, NULL, NULL);
if (pDstDS == NULL) {
printf("创建失败!\n");
GDALClose(pDstDS);
return NULL;
}
if (pNodata) {
pDstDS->GetRasterBand(1)->SetNoDataValue(*pNodata);
}
GDALClose(pSrcDS);
return pDstDS;
}
// 手动读取
GDALDataset* grd2gdal_manual(const char* pszSrcFile, const char* pszDstFile, const char* pFormat, DT_32F* pNodata = NULL)
{
FILE* fout = NULL;
// fout = fopen(pszSrcFile, "rb");
fopen_s(&fout, pszSrcFile, "rb");
if (fout == NULL) {
printf("无法打开文件: %s\n", pszSrcFile);
return NULL;
}
// 读取文件抬头信息
DSBBHeader header;
fread(&header, sizeof(DSBBHeader), 1, fout);
printf("\n--- Header Info ---\n");
PrintDSBBInfo(header);
// 读取数据部分
DT_32F* pdata = new DT_32F[header.sizex * header.sizey];
fread(pdata, sizeof(DT_32F), header.sizex * header.sizey, fout);
// 创建输出文件
GDALAllRegister();
GDALDriver* pDriver = GetGDALDriverManager()->GetDriverByName(pFormat);
if (pDriver == NULL) {
printf("不支持格式: %s\n", pFormat);
fclose(fout);
return NULL;
}
GDALDataset* pDstDS = pDriver->Create(pszDstFile, header.sizex, header.sizey, 1, GDT_Float32, NULL);
if (pDstDS == NULL) {
printf("创建失败!\n");
fclose(fout);
return NULL;
}
// 设置nodata
if (pNodata) {
pDstDS->GetRasterBand(1)->SetNoDataValue(*pNodata);
}
// 计算仿射参数
double adftrans[6] = {
0.0 };
adftrans[1] = (header.xmax - header.xmin) / static_cast<double>(header.sizex - 1);
adftrans[5] = (header.ymin - header.ymax) / static_cast<double>(header.sizey - 1);
adftrans[0] = header.xmin - 0.5 * adftrans[0];
adftrans[3] = header.ymax - 0.5 * adftrans[5];
pDstDS->SetGeoTransform(adftrans);
// 从最后一行开始写入数据
for (int y = header.sizey - 1; y >= 0; --y) {
pDstDS->GetRasterBand(1)->RasterIO(GF_Write, 0, y, header.sizex, 1,
pdata + header.sizex * (header.sizey - 1 - y), header.sizex, 1, GDT_Float32, 0, 0, 0);
}
fclose(fout);
return pDstDS;
}
// 将 Dataset 保存为 grd 文件
bool gdal2grd(GDALDataset *pSrcDS, const char *pGrdFile)
{
if (pSrcDS == NULL) {
printf("数据指针为空!\n");
return false;
}
DSBBHeader header;
int sizex = pSrcDS->GetRasterXSize();
int sizey = pSrcDS->GetRasterYSize();
double adftrans[6];
if (pSrcDS->GetGeoTransform(adftrans) != CE_None) {
printf("图像无仿射系数!\n");
GDALClose(pSrcDS);
return false;
}
GDALRasterBand* pband = pSrcDS->GetRasterBand(1);
if (pband->GetRasterDataType() != GDT_Float32) {
printf("图像数据类型不是 float!\n");
GDALClose(pSrcDS);
return false;
}
int hasNodata = 0;
DT_32F nodata = pband->GetNoDataValue(&hasNodata);
header.id[0] = 'D';
header.id[1] = 'S';
header.id[2] = 'B';
header.id[3] = 'B';
header.sizex = sizex;
header.sizey = sizey;
header.xmin = adftrans[0] + adftrans[1] / 2;
header.ymax = adftrans[3] + adftrans[5] / 2;
header.xmax = header.xmin + (static_cast<double>(sizex) - 1.0) * adftrans[1];
header.ymin = header.ymax + (static_cast<double>(sizey) - 1.0) * adftrans[5];
FILE* fout = NULL;
fopen_s(&fout, pGrdFile, "wb");
if (fout == NULL) {
printf("文件创建失败!\n");
GDALClose(pSrcDS);
return false;
}
fwrite(&header, sizeof(DSBBHeader), 1, fout);
DT_64F maxv = -1.0;
DT_64F minv = 99999.0;
DT_32F* pdata = new DT_32F[sizex];
for (int y = sizey - 1; y >= 0; --y) {
pband->RasterIO(GF_Read, 0, y, sizex, 1, pdata, sizex, 1, GDT_Float32, 0, 0, 0);
for (int x = 0; x < sizex; ++x) {
if (hasNodata && pdata[x] == nodata)
continue;
maxv = max(maxv, static_cast<DT_64F>(pdata[x]));
minv = min(minv, static_cast<DT_64F>(pdata[x]));
}
fwrite(pdata, sizeof(DT_32F), sizex, fout);
}
delete[] pdata;
header.zmin = static_cast<DT_64F>(minv);
header.zmax = static_cast<DT_64F>(maxv);
fseek(fout, 40, SEEK_SET);
fwrite(&header.zmin, sizeof(DT_64F), 1, fout);
fwrite(&header.zmax, sizeof(DT_64F), 1, fout);
PrintDSBBInfo(header);
fclose(fout);
return true;
}
// 打印grd header信息
void PrintDSBBInfo(const DSBBHeader & header)
{
printf("ID: %c%c%c%c\n", header.id[0], header.id[1], header.id[2], header.id[3]);
printf("Size: %d x %d\n", header.sizex, header.sizey);
double scales[2];
scales[0] = (header.xmax - header.xmin) / static_cast<double>(header.sizex - 1);
scales[1] = (header.ymin - header.ymax) / static_cast<double>(header.sizey - 1);
printf("scale: %.6f x %.6f\n", scales[0], scales[1]);
printf("Center Points(WN): (%.6f, %.6f)\n", header.xmin, header.ymax);
printf("Center Points(ES): (%.6f, %.6f)\n", header.xmax, header.ymin);
double tiePointWN[2], tiePointES[2];
tiePointWN[0] = header.xmin - 0.5 * scales[0];
tiePointWN[1] = header.ymax - 0.5 * scales[1];
tiePointES[0] = header.xmax + 0.5 * scales[0];
tiePointES[1] = header.ymin + 0.5 * scales[1];
printf("Tie points(WN): (%.6f, %.6f)\n", tiePointWN[0], tiePointWN[1]);
printf("Tie points(ES): (%.6f, %.6f)\n", tiePointES[0], tiePointES[1]);
printf("Data range: %.6f ~ %.6f\n", header.zmin, header.zmax);
}
// 打印 Dataset 信息
void PrintDSInfo(GDALDataset* pImgDS)
{
if (pImgDS == NULL) {
printf("输入图像为空!\n");
return;
}
int sizex = pImgDS->GetRasterXSize();
int sizey = pImgDS->GetRasterYSize();
int nband = pImgDS->GetRasterCount();
printf("Size: %d x %d\n", sizex, sizey);
printf("band: %d\n", nband);
double adftrans[6];
pImgDS->GetGeoTransform(adftrans);
printf("scale: %.6f x %.6f\n", adftrans[1], adftrans[5]);
printf("Tie points(WN): (%.6f, %.6f)\n", adftrans[0], adftrans[3]);
double tiePointES[2];
tiePointES[0] = adftrans[0] + adftrans[1] * sizex;
tiePointES[1] = adftrans[3] + adftrans[5] * sizey;
printf("Tie points(ES): (%.6f, %.6f)\n", tiePointES[0], tiePointES[1]);
double cenPointWN[2];
cenPointWN[0] = adftrans[0] + 0.5 * adftrans[1];
cenPointWN[1] = adftrans[3] + 0.5 * adftrans[5];
double cenPointES[2];
cenPointES[0] = tiePointES[0] - 0.5 * adftrans[1];
cenPointES[1] = tiePointES[1] - 0.5 * adftrans[5];
printf("Center Points(WN): (%.6f, %.6f)\n", cenPointWN[0], cenPointWN[1]);
printf("Center Points(ES): (%.6f, %.6f)\n", cenPointES[0], cenPointES[1]);
int isNodata = 0;
printf("Nodata: ");
pImgDS->GetRasterBand(1)->GetNoDataValue(&isNodata);
if (isNodata) {
printf("%.3f\n", pImgDS->GetRasterBand(1)->GetNoDataValue());
}
else {
printf("None\n");
}
GDALDataType eDT = pImgDS->GetRasterBand(1)->GetRasterDataType();
printf("Data Type: %s\n", getDataTypeName(eDT));
if (eDT == GDT_Float32) {
DT_32F* pdata = new DT_32F[sizex];
DT_32F maxv = -1.0f;
DT_32F minv = 99999.0f;
for (int y = 0; y < sizey; ++y) {
pImgDS->GetRasterBand(1)->RasterIO(GF_Read, 0, y, sizex, 1, pdata, sizex, 1, eDT, 0, 0, 0);
for (int x = 0; x < sizex; ++x) {
if (isNodata && pdata[x] == pImgDS->GetRasterBand(1)->GetNoDataValue()) {
continue;
}
maxv = max(pdata[x], maxv);
minv = min(pdata[x], minv);
}
}
printf("Data range: %.6f ~ %.6f\n", minv, maxv);
}
}
今天的文章gdal读取dwg文件_广联达保存的文件怎么读取不了「建议收藏」分享到此就结束了,感谢您的阅读。
版权声明:本文内容由互联网用户自发贡献,该文观点仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 举报,一经查实,本站将立刻删除。
如需转载请保留出处:https://bianchenghao.cn/88542.html