8

I have read pixel values of an raster image using GDAL libraries in visual studio 2010(vc++).

Now I have to crop the image (subset) according to the grid given in shape file.

I have got a code for cropping image on internet as given below, but its having some error in it. Please review it and help .

GDALDataset* crop(GDALDataset* poDS, int xoff, int yoff, int xsize, int ysize) { VRTDataset *poVDS = (VRTDataset*)VRTCreate(xsize, ysize); // Copy dataset info const char* pszProjection = poDS->GetProjectionRef(); if (pszProjection != NULL && strlen(pszProjection) > 0) poVDS->SetProjection(pszProjection); double adfGeoTransform[6]; if (poDS->GetGeoTransform(adfGeoTransform) == CE_None) { // Adapt the geotransform matrix to the subarea adfGeoTransform[0] += xoff * adfGeoTransform[1] + yoff * adfGeoTransform[2]; adfGeoTransform[3] += xoff * adfGeoTransform[4] + yoff * adfGeoTransform[5]; poVDS->SetGeoTransform(adfGeoTransform); } poVDS->SetMetadata(poDS->GetMetadata()); // Here I also copy metadata from my own domain char **papszMD; papszMD = poDS->GetMetadata(MD_DOMAIN_MSAT); if (papszMD != NULL) poVDS->SetMetadata(papszMD, MD_DOMAIN_MSAT); for (int i = 0; i < poDS->GetRasterCount(); ++i) { GDALRasterBand* poSrcBand = poDS->GetRasterBand(i + 1); GDALDataType eBandType = poSrcBand->GetRasterDataType(); poVDS->AddBand(eBandType, NULL); VRTSourcedRasterBand* poVRTBand = (VRTSourcedRasterBand*)poVDS- >GetRasterBand(i + 1); poVRTBand->AddSimpleSource(poSrcBand, xoff, yoff, xsize, ysize, 0, 0, xsize, ysize); poVRTBand->CopyCommonInfoFrom(poSrcBand); // Again, I copy my own metadata papszMD = poSrcBand->GetMetadata(MD_DOMAIN_MSAT); if (papszMD != NULL) poVRTBand->SetMetadata(papszMD, MD_DOMAIN_MSAT); } return poVDS; }
Improve this question
1

3 Answers 3

Reset to default
12

Unprojected:

gdal_translate -srcwin xstart ystart xstop ystop input.raster output.raster

Projected:

gdal_translate -projwin ulx uly lrx lry input.raster output.raster

If you are using gdal on your system, gdal_translate is installed.

Improve this answer
5
  • Hello jay laura: thank you for the reply but the lines you have suggested works for python i guess. i only uses c++ code. So i dont think gdal_translate work in c++ coding . Commented Jan 8, 2013 at 17:01
  • 1
    @hardik, gdal_translate is written in C++ and operates independently of the gdal python bindings, here's the code: trac.osgeo.org/gdal/browser/trunk/gdal/apps/gdal_translate.cpp. It's part of the core utility set and a compiled binary is generally installed by default (unless you're not installing and only building from source). Commented Mar 22, 2013 at 16:22
  • 1
    The unprojected command seems inconsistent with the documentation gdal.org/gdal_translate.html, which says -srcwin xoff yoff xsize ysize. Commented Feb 17, 2016 at 18:44
  • Can you explain what you mean by "projected" and "unprojected". My raster has a projection, but I don't want to reproject it — which form is correct? Commented Dec 2, 2019 at 18:59
  • @colllin "projected" means you can specify the coordinates using the coordinate format of the projection (e.g. for example if using EPSG:4326 you can provide 90/180 lat/lon coordinates) and "unprojected" means you specify coordinates as x/y pixel values of the raster. For projected, you can specify the SRS using -projwin_srs. It does not actually reproject anything in the output. Commented May 21, 2021 at 8:08
4

See the GDALRasterBand Class. Carefully using the nXOff and nXSize (etc.) parameters of GDALRasterBand::RasterIO will read a "cropped" part of the raster.

Improve this answer
3

Thank you Mike. I have solved the problem of cropping the image using VC++ with gdal libraries. I have created VRTDataset of my desired size of raster to be cropped and then save it using CreateCopy().

Improve this answer

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.