Create a GeoTIFF from a GeoTIFF using only one of it's values

Question

Does anyone have a method for extracting pixels from a GeoTIFF that are of a particular value. So for example if a band had values from 0 to 10 and I wanted to make an image just using the zero values and to leave the other pixel values as unknown.

My GeoTIFF is 19 GB large so I'd prefer to avoid using QGIS. I've already tried GRASS and it would freeze on me. I would prefer a solution that can use the GDAL tools but Python would be good too.

My eventual goal is to extract the polygons that surround regions of interest using the monochrome image generated from the above.

Answer 1

I've not tried this but I think the example in the gdal_calc manual page

gdal_calc.py -A input.tif --outfile=result.tif --calc="A*(A>0)" --NoDataValue=0

could be changed to

gdal_calc.py -A input.tif --outfile=result.tif --calc="A*(A==1)" --NoDataValue=0

should work to extract value 1.

Answer 2

this is one way to do what you want using GDAL and NUMPY.

from osgeo import gdal
import numpy as np

rasterFile = 'C:/Users/pepito/yourinputraster.tif'
in_ds = gdal.Open(rasterFile)
in_band = in_ds.GetRasterBand(1)
in_data = in_band.ReadAsArray()

#this is the line that makes the reclass
reclas = np.where(in_data == 0, 1, 0)

# output as TIFF
gtiff_driver = gdal.GetDriverByName('GTiff')
out_ds = gtiff_driver.Create('C:/Users/pepito/youroutputraster.tif',in_band.XSize, in_band.YSize, 1, in_band.DataType)
out_ds.SetProjection(in_ds.GetProjection())
out_ds.SetGeoTransform(in_ds.GetGeoTransform())
out_band = out_ds.GetRasterBand(1)
out_band.SetNoDataValue(0)
out_band.WriteArray(reclas)
out_ds.FlushCache()

However if your raster is 19GB you will probably have to process it in chunks.