I am trying to do a raster calculations on a 28 GB GeoTiff using gdal_calc.py (3.0.4) but the process fails

gdal_calc.py -A DEM_0.tif --outfile=DEM_1.tif --calc="A+0.33" --format=GTiff --type=Float32 --co="COMPRESS=LZW"

when the size of the output file becomes ~650 MB the command keeps running but the file stops increasing in size. It happens if I run the process on the SSD internally but also on an external SSD. Furthermore, I noticed that the rate of increase of the filesize is only 70MB / min (which is rather slow I think, although gdal_calc only uses one CPU).

Thanks to a comment below I also created another geotiff with tilesize 64x64. I ran gdal_calc.py overnight and pretty much the same thing happened - it was still running after more than 7 hours so I killed the process.

GDAL or other software solutions welcome.

shortened gdalinfo:

Size is 146140, 86076


Band 1 Block=256x256 Type=Float32, ColorInterp=Gray NoData Value=-32767 Overviews: 73070x43038, 36535x21519, 18268x10760, 9134x5380, 4567x2690, 2284x1345, 1142x673, 571x337, 286x169, 143x85

Unit Type: metre

  • These big raster jobs can be tough to resolve once you start suspecting memory management. Out of curiosity how much RAM does your system have? Windows or..? – elrobis Aug 4 '20 at 14:56
  • Linux mint with 24 GB RAM – EOF Aug 4 '20 at 15:01
  • Interesting. The image's demand just barely exceeds your available RAM. I wonder if you could, temporarily, increase your swap space to the full extend of your external drive to see if it'll get through the task? – elrobis Aug 4 '20 at 18:12
  • Thanks. I have 15 GB swap now but don't have space to increase to the full extend of the external drive (250GB). It wouldn't be a good permanent solution either. – EOF Aug 4 '20 at 18:39
  • 2
    Tiled output should be better because otherwise GDAL must write out 146140 pixels wide rows. Original image has 246x256 sized tiles and thus one row is covered by 571 tiles. For getting data for one row 571x256x256 pixels must be read, which means 40098475 pixels. If input and output are both tiled then at least in theory GDAL should be able to work tile by tile. – user30184 Aug 4 '20 at 18:52

Although it isn't a true gdal_calc answer it solves the issue of working with a 28 GB raster using GRASS through R.


dem  <- raster('DEM.tif) 
bb   <- st_bbox(dem)

initGRASS(gisBase = "/usr/lib/grass78", 
          home = '/home/user/tmp', # to have sufficient disk space
          gisDbase = '/home/user/tmp', 
          location = "wadden", 
          mapset = "PERMANENT",
          override = TRUE)

execGRASS("g.proj", flags = c("c", "quiet"), 
          proj4 = st_crs(dem)$proj4string)

execGRASS("g.region", flags = c("quiet"), 
          n = as.character(bb["ymax"]), s = as.character(bb["ymin"]), 
          e = as.character(bb["xmax"]), w = as.character(bb["xmin"]), 
          res = paste(res(dem)[1])) 

execGRASS("r.in.gdal", flags=c("o", "overwrite"), parameters=list(input='DEM.tif'), output="dem"))
execGRASS("r.mapcalc", expression="dem_fixed = dem + 0.33")
execGRASS("r.out.gdal", parameters=list(
                          createopt=c("COMPRESS=LZW",  "BIGTIFF=YES")),
                        flags=c("f", "overwrite") )

The issue might have been caused by an older version of gdal_calc.py Versions starting from 3.3 have useful new options, so it is good to try with the newest gdal again. Options like TILED=YES or COMPRESS=DEFLATE can be added with the syntax: --co="TILED=YES" one can define many options like this so I use always --co='COMPRESS=DEFLATE' --co='PREDICTOR=YES' --co='BIGTIFF=IF_SAFER' --co='TILED=YES'

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.