GDALワープC ++ APIを使用したラスターファイルのトリミング

Question

DEMラスターファイルをシェープファイルでトリミングする必要があります。このシェープファイルには多くのポリゴンがあります。トリミング操作後にいくつかの変数を計算する必要があります。

アルゴリズムはC ++で書かれているので、GDAL C ++ APIを呼び出したいと思います。ワープしたファイルをディスクに書き出す必要はありません。テストに合格したら、VRTまたはMEMのいずれかを使用する予定です。

すでに各ポリゴンの座標を取得しており、cutline引数を使用してポリゴンとして渡しましたが、今のところ操作は効果がありません。

            OGRLinearRing poExteriorRing;
            for (i = 0; i < nPt; i++) 
              {
                dX = (*iIterator).vVertex.at(i).dX;
                dY = (*iIterator).vVertex.at(i).dY;
                double dDummy1 = (dX - dX_origin) / dResolution_elevation;
                long lColumn_index = long(round(dDummy1));
                double dDummy2 = (dY_origin - dY) / dResolution_elevation;
                long lRow_index = long(round(dDummy2));     
                poExteriorRing.addPoint(lColumn_index, lRow_index);
              }
            poExteriorRing.closeRings();
            OGRPolygon polygon;
            polygon.addRing(&poExteriorRing);

            GDALDataset  *poVRTDS;
            poDS_elevation->GetGeoTransform(adfGeoTransform);
            poVRTDS = poDriver->CreateCopy( "./test.tiff", poDS_elevation, FALSE, NULL, NULL, NULL );

            GDALWarpOptions *psWarpOptions = GDALCreateWarpOptions();
            ...
            psWarpOptions->hCutline = &polygon;                
            GDALWarpOperation oOperation;
            oOperation.Initialize( psWarpOptions );
            oOperation.ChunkAndWarpImage( 0, 0,
                            GDALGetRasterXSize( poDS_elevation ),
                            GDALGetRasterYSize( poDS_elevation ) );
            GDALClose( poVRTDS );

次のように、ポリゴンの頂点から画像のインデックスを作成する必要があると思います。CUTLINE：これにはカットラインのWKTジオメトリが含まれている場合があります。GDALWarpOperation :: Initialize（）によってジオメトリに変換され、GDALWarpOptionshCutlineフィールドに割り当てられます。座標は、ソースピクセル/ライン座標で表す必要があります。注：これは、gdalwarpユーティリティの-cutlineオプションに対して行われた仮定とは異なります。https://gdal.org/api/gdalwarp_cpp.html

Answer 1

何度か試みた後、いくつかのブログ投稿に従って問題を修正することができました。 http://www.jeepxie.net/article/293437.html

コア部分は、座標を変換する必要があるということです。

eError = TransformCutlineToSource( hWrkSrcDS, hCutline,
                                  &(psWO->papszWarpOptions),
                                  psOptions->papszTO );

また、hCutlineとcutlineの両方を同時に使用しないでください。