[gdal-dev] wrong with clipping a tif file with a polygon

[Norman Vine]

perhaps this will help
https://stackoverflow.com/questions/70664043/clip-raster-with-polygon-with-gdal-c

> On May 30, 2023, at 8:12 AM, Ahmet Temiz <ahmettemiz88 at gmail.com> wrote:
> 
> Hi,
> I am trying to clip a tif file with a polygon. But resulting output.tif file shows nothing.
> 
> Where am I doing wrong ?
> 
> "  #include <gdal/gdal.h>
> #include <gdal/gdal_priv.h>
> #include <gdal/gdalwarper.h>
> #include <gdal/ogrsf_frmts.h>
> #include <gdal/gdalwarper.h>
> #include <string> // for string class
> using namespace std;
> 
> int main() {
>   const char *inputPath = "input.tif";
>   const char *outputPath = "output.tif";
> 
>   // clipper Polygon
>   // THIS FILE MUST BE IN PIXEL/LINE COORDINATES or otherwise one should
>   // copy the function gdalwarp_lib.cpp:TransformCutlineToSource()
>   // from GDAL's sources
>   // It is expected that it contains a single polygon feature
>   const char *read_filenamePoly = "input.geojson";
> 
>   GDALDataset *hSrcDS;
>   GDALDataset *hDstDS;
> 
>   GDALAllRegister();
>   auto poDriver = GetGDALDriverManager()->GetDriverByName("GTiff");
>   hSrcDS = (GDALDataset *)GDALOpen(inputPath, GA_ReadOnly);
> 
>   hDstDS = (GDALDataset *)poDriver->CreateCopy(
>   outputPath, hSrcDS, 0, nullptr, nullptr, nullptr);
>   // Without this step the cutline is useless - because the background
>   // will be carried over from the original image
>   hDstDS->GetRasterBand(1)->Fill(0);
> 
>   hDstDS->SetProjection(hSrcDS->GetProjectionRef());
> 
>   const char *src_srs = hSrcDS->GetProjectionRef();
>   const char *dst_srs = hDstDS->GetProjectionRef();
>   printf("%s ",src_srs);
>   // clipper Layer
>   GDALDataset *poDSClipper;
>   poDSClipper = (GDALDataset *)GDALOpenEx(
>     read_filenamePoly, GDAL_OF_UPDATE, NULL, NULL, NULL);
> 
>   auto poLayerClipper = poDSClipper->GetLayer(0);
>   auto geom = poLayerClipper->GetNextFeature()->GetGeometryRef();
>   //char* wkt = nullptr;
>   if (geom != nullptr && geom->IsValid())
>   { 
>       printf("\n geom is valid \n");
>   }
>   // setup warp options
>   GDALWarpOptions *psWarpOptions = GDALCreateWarpOptions();
>   psWarpOptions->hSrcDS = hSrcDS;
>   psWarpOptions->hDstDS = hDstDS;
>   psWarpOptions->nBandCount = 1;
>   psWarpOptions->panSrcBands =
>     (int *)CPLMalloc(sizeof(int) * psWarpOptions->nBandCount);
>   psWarpOptions->panSrcBands[0] = 1;
>   psWarpOptions->panDstBands =
>     (int *)CPLMalloc(sizeof(int) * psWarpOptions->nBandCount);
>   psWarpOptions->panDstBands[0] = 1;
>   psWarpOptions->pfnProgress = GDALTermProgress;
>   psWarpOptions->hCutline = geom;
> 
>   // Establish reprojection transformer.
>   psWarpOptions->pTransformerArg = GDALCreateGenImgProjTransformer(
>     hSrcDS, src_srs, hDstDS, dst_srs, TRUE, 1000, 1);
>   psWarpOptions->pfnTransformer = GDALGenImgProjTransform;
>   GDALWarpOperation oOperation;
>   oOperation.Initialize(psWarpOptions);
>   oOperation.ChunkAndWarpImage(
>     0, 0, GDALGetRasterXSize(hDstDS), GDALGetRasterYSize(hDstDS));
>   GDALDestroyGenImgProjTransformer(psWarpOptions->pTransformerArg);
>   GDALDestroyWarpOptions(psWarpOptions);
>   GDALClose(hDstDS);
>   GDALClose(hSrcDS);
> }
> "
> 
> -- 
> Ahmet Temiz
> Jeoloji Müh.
> Afet ve Acil Durum Yönetimi Başkanlığı
> Deprem  Dairesi Başkanlığı
> 
> 
> ________________________
> 
> Ahmet Temiz
> Geological Eng.
> 
> Disaster and Emergency Management
> of Presidency
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev