私はGDALライブラリを使用しています。現在、左上の点と右上の点を取り込んで、元の画像から切り出すことができます。私が今やりたいことは、2 つの WKT ポイントを取り、X、Y 座標に変換して同じことを行うことです。GeoTransform とそれが使用している座標系 (WGS84) を知っていれば、これが可能かどうか疑問に思っていました。
3 に答える
私も以前にこれに出くわしましたが、これは座標変換を行うのに十分な方法です。
GDAL ドキュメントからの注意:
GDALDataset::GetProjectionRef() によって返される座標系は、GDALDataset::GetGeoTransform() によって返されるアフィン地理参照変換によって暗示された地理参照座標を記述します。
これをOGRCoordinateTransformationと一緒に使用して、変換を行うことができます。
基本的に、コードは次のようになります。
// Load up some dataset.
dataset = (GDALDataset *) GDALOpen( mapfile, GA_ReadOnly );
// Define Geographic coordinate system - set it to WGS84.
OGRSpatialReference *poSRS_Geog = new OGRSpatialReference();
poSRS_Geog->importFromEPSG( 4326 ); // WGS84
// Define Projected coordinate system - set to the GeoTransform.
const char *sProj = dataset->GetProjectionRef();
OGRSpatialReference *poSRS_Proj = new OGRSpatialReference( sProj );
// Set up the coordinate transform (geographic-to-projected).
OGRCoordinateTransformation *poCT_Geog2Proj;
poCT_Geog2Proj = OGRCreateCoordinateTransformation( poSRS_Geog, poSRS_Proj );
// Now everything is set up and we set transforming coordinates!
// Pass Lon/Lat coordinates to the Transform function:
double x = lon;
double y = lat;
poCT_Geog2Proj->Transform( 1, &x, &y );
// Now x and y variables will contain the X/Y pixel coordinates.
これが、経度/緯度とピクセル座標の間で変換する方法です。Transform()
で配列を使用し、複数の座標をまとめて変換できることに注意してください。最初の引数は変換する座標ペアの数で、2 番目と 3 番目の引数は x と y へのポインターです。ここで 1 組だけ変換します。
逆変換を設定するのも同様に簡単であることに注意してください。
// Set up the coordinate transform (projected-to-geographic).
OGRCoordinateTransformation *poCT_Proj2Geog;
poCT_Proj2Geog = OGRCreateCoordinateTransformation( poSRS_Proj, poSRS_Geog );
画像のアフィン変換を使用して、サンプルの緯度/経度を計算しました。私が抱えていた唯一の問題は、画像が北向きの場合、緯度/経度を計算するときに geoTransform[2] と geTransform[4] をゼロにする必要があることでした。
x = (int)Math.Abs(Math.Round((Latitude - geotransform[0]) / geotransform[1]));
y = (int)Math.Abs(Math.Round((Longitude - geotransform[3]) / geotransform[5]));
ブルートフォースしたい場合は、次のようにすることができます(私はこれを実行し、機能しましたが、これは単なる擬似コードです):
//Get the Pixel for the length and width, this portion is for the full image
pixelXSize = AbsoluteValue((latitudeAt(Zero)-(latitudeAt(Length)-1))/imageLength);
pixelYSize = AbsoluteValue((longitudeAt(Zero)-(LongitudeAt(Width)-1))/imageWidth);
//Calculate the x,y conversion for the points you want to calculate
x = AbsoluteValue((latitudeToConvert-latitudeAt(Zero))/pixelXSize);
y = AbsoluteValue((longitudeToConvert-longitudteAt(Zero))/pixelYSize);
この答えは、1 ~ 2 ピクセルずれている可能性があります。geotransform を使用している場合、North Facing 変数が返される回答を台無しにする可能性があります。これまでのところ、北向きの画像に対してのみこれをテストしました。