GDAL數據集寫入空間坐標參考

目錄

• 1. 概述
• 2. 柵格數據
• 3. 矢量數據

1. 概述

可以通過GDAL給地理數據寫入空間參考信息，不過要注意的是GDAL給矢量數據和柵格數據寫入空間坐標參考的接口不太一樣。

2. 柵格數據

實現代碼如下：

#include <iostream>
#include <gdal_priv.h>
#include <string>

using namespace std;

int main()
{
	GDALAllRegister();
	CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");  //支持中文路徑
	CPLSetConfigOption("SHAPE_ENCODING", "");  //解決中文亂碼問題	
	CPLSetConfigOption("GDAL_DATA", "D:/Work/GDALBuild/gdal-2.4.2/install/data");

	GDALDriver *pDriver = GetGDALDriverManager()->GetDriverByName("GTIFF"); //圖像驅動
	char** ppszOptions = NULL;
	ppszOptions = CSLSetNameValue(ppszOptions, "BIGTIFF", "IF_NEEDED"); //配置圖像信息
	const char* dstPath = "dst.tif";
	GDALDataset* dst = pDriver->Create(dstPath, 256, 256, 3, GDT_Byte, ppszOptions);
	if (dst == nullptr)
	{
		printf("Can't Write Image!");
		return false;
	}

	//空間參考
	OGRSpatialReference spatialReference;
	spatialReference.importFromEPSG(4326);				//wgs84地理坐標系
	char *pszWKT = nullptr;
	spatialReference.exportToWkt(&pszWKT);

	dst->SetProjection(pszWKT);

	CPLFree(pszWKT);
	pszWKT = nullptr;

	//坐標信息
	double padfTransform[6] = {
		114.0,		//左上角點坐標X
		0.000001,		//X方向的分辨率
		0,		//旋轉系數，如果為0，就是標准的正北向圖像
		34.0,			//左上角點坐標Y
		0,			//旋轉系數，如果為0，就是標准的正北向圖像
		0.000001,			//Y方向的分辨率
	};
	dst->SetGeoTransform(padfTransform);

	GDALClose(dst);
}

這里創建了一個wgs84地理坐標系空間參考的柵格數據，通過OGRSpatialReference類導出了描述空間參考的wkt字符串，寫入到GDAL數據集中。

3. 矢量數據

實現代碼如下：

#include <iostream>
#include <gdal_priv.h>
#include <ogrsf_frmts.h>

using namespace std;

int main()
{
	GDALAllRegister();
	CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");  //支持中文路徑
	CPLSetConfigOption("SHAPE_ENCODING", "");  //解決中文亂碼問題	
	CPLSetConfigOption("GDAL_DATA", "D:/Work/GDALBuild/gdal-2.4.2/install/data");

	//空間參考
	OGRSpatialReference spatialReference;
	spatialReference.importFromEPSG(4326);				//wgs84地理坐標系
	
	//創建
	GDALDriver* driver = GetGDALDriverManager()->GetDriverByName("ESRI Shapefile");
	if (!driver)
	{
		printf("Get Driver ESRI Shapefile Error！\n");
		return false;
	}

	GDALDataset* dataset = driver->Create("dst.shp", 0, 0, 0, GDT_Unknown, NULL);
	OGRLayer* poLayer = dataset->CreateLayer("houseType", &spatialReference, wkbPolygon, NULL);

	//創建屬性字段
	{
		// 字符串
		OGRFieldDefn oField1("名稱", OFTString);
		oField1.SetWidth(8);
		if (poLayer->CreateField(&oField1) != OGRERR_NONE) {
			printf("Creating Name field failed.\n"); return FALSE;
		}

		// 浮點數
		OGRFieldDefn oField2("面積", OFTReal);
		oField2.SetPrecision(3);
		if (poLayer->CreateField(&oField2) != OGRERR_NONE) {
			printf("Creating Name field failed.\n"); return FALSE;
		}

		// 整型
		OGRFieldDefn oField3("結點數", OFTInteger);
		if (poLayer->CreateField(&oField3) != OGRERR_NONE) {
			printf("Creating Name field failed.\n"); return FALSE;
		}
	}

	//創建特征
	OGRFeature *poFeature = new OGRFeature(poLayer->GetLayerDefn());

	OGRLinearRing ogrring;
	int pNum = 4;
	ogrring.setNumPoints(pNum);
	ogrring.setPoint(0, 114.0, 34.0, 0.0);
	ogrring.setPoint(1, 115.0, 34.0, 0.0);
	ogrring.setPoint(2, 115.0, 35.0, 0.0);
	ogrring.setPoint(3, 114.0, 35.0, 0.0);
	   
	OGRPolygon polygon;
	polygon.addRing(&ogrring);
	poFeature->SetGeometry(&polygon);

	poFeature->SetField("名稱", "多邊形");
	poFeature->SetField("面積", polygon.get_Area());
	poFeature->SetField("結點數", pNum);

	if (poLayer->CreateFeature(poFeature) != OGRERR_NONE)
	{
		printf("Failed to create feature in shapefile.\n");
		return false;
	}

	//釋放
	GDALClose(dataset);
	dataset = nullptr;
}

與寫入到柵格數據不同，空間參考信息寫入到矢量數據是寫入到GDAL數據集的圖層類中的，並且直接傳入OGRSpatialReference類即可。