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

Ahmet Temiz ahmettemiz88 at gmail.com
Tue May 30 10:59:15 PDT 2023


Thank you..
I solved the problem.  Cutting polygon has to be converted to pixel
coordinates.  I have missed that point. now it works well ..

Regards

30 May 2023 Sal 16:20 tarihinde Rahkonen Jukka <
jukka.rahkonen at maanmittauslaitos.fi> şunu yazdı:

> Hi,
>
>
>
> I can’t really read the code, but based on what happens, maybe your
> cutline is outside the image. So you
>
>    - create an output band
>    - fill it with zeroes
>    - and then you write pixels within the cutline into the output
>
> Because the output stays all zeroes, and there are no errors, then I guess
> that the pixels within your cutline are all nodata.
>
> Maybe printing the WKT out from “geom” could prove or disprove my guess.
>
>
>
> -Jukka Rahkonen-
>
>
>
> *Lähettäjä:* gdal-dev <gdal-dev-bounces at lists.osgeo.org> *Puolesta *Ahmet
> Temiz
> *Lähetetty:* tiistai 30. toukokuuta 2023 15.57
> *Vastaanottaja:* Norman Vine <nhv at meganet.net>
> *Kopio:* gdal-dev <gdal-dev at lists.osgeo.org>
> *Aihe:* Re: [gdal-dev] wrong with clipping a tif file with a polygon
>
>
>
> Thank you ,
>
>
>
> I already worked on that code. But clipping was not succesful. output.tif
> is blank on qgis.
>
>
>
> On Tue, May 30, 2023 at 3:24 PM Norman Vine <nhv at meganet.net> wrote:
>
> 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
>
>
>
>
> --
>
> 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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20230530/366785ab/attachment-0001.htm>


More information about the gdal-dev mailing list