[GIS] Clip raster using mask other raster using Python GDAL

gdalpython

I want to clip raster using for mask other raster and I have find this question and I try to follow this code:

import gdal
from gdalconst import GA_ReadOnly

data = gdal.Open('img_mask.tif', GA_ReadOnly)
geoTransform = data.GetGeoTransform()
minx = geoTransform[0]
maxy = geoTransform[3]
maxx = minx + geoTransform[1] * data.RasterXSize
miny = maxy + geoTransform[5] * data.RasterYSize
call('gdal_translate -projwin ' + ' '.join([str(x) for x in [minx, maxy, maxx, miny]]) + ' -of GTiff img_orig.tif img_out.tif', shell=True)

But that code doesn't working the output TIFF have some extent with original TIFF against the clip doesn't work. Any idea?

Best Answer

You should use gdal.Translate method. Something like the following should work.

from osgeo import gdal
from gdalconst import GA_ReadOnly

maskDs = gdal.Open('myMask.tif', GA_ReadOnly)# your mask raster
projection=maskDs.GetProjectionRef()
geoTransform = maskDs.GetGeoTransform()
minx = geoTransform[0]
maxy = geoTransform[3]
maxx = minx + geoTransform[1] * maskDs.RasterXSize
miny = maxy + geoTransform[5] * maskDs.RasterYSize


data=gdal.Open('myRaster.tif', GA_ReadOnly) #Your data the one you want to clip
output='output.tif' #output file
gdal.Translate(output,data,format='GTiff',projWin=[minx,maxy,maxx,miny],outputSRS=projection)

Related Solutions

[GIS] Python GDAL: Save array as raster with projection from other file

Eddy, try a different approach to write the array as a raster. The following code works for me, and for your case you would only need to change the corner coordinate, pixel size and number of pixels for the reprojected values, as you did in your code. So the first part would be replaced by some of your procedures.

from osgeo import gdal

def array_to_raster(array):
    """Array > Raster
    Save a raster from a C order array.

    :param array: ndarray
    """
    dst_filename = '/a_file/name.tiff'


    # You need to get those values like you did.
    x_pixels = 16  # number of pixels in x
    y_pixels = 16  # number of pixels in y
    PIXEL_SIZE = 3  # size of the pixel...        
    x_min = 553648  
    y_max = 7784555  # x_min & y_max are like the "top left" corner.
    wkt_projection = 'a projection in wkt that you got from other file'

    driver = gdal.GetDriverByName('GTiff')

    dataset = driver.Create(
        dst_filename,
        x_pixels,
        y_pixels,
        1,
        gdal.GDT_Float32, )

    dataset.SetGeoTransform((
        x_min,    # 0
        PIXEL_SIZE,  # 1
        0,                      # 2
        y_max,    # 3
        0,                      # 4
        -PIXEL_SIZE))  

    dataset.SetProjection(wkt_projection)
    dataset.GetRasterBand(1).WriteArray(array)
    dataset.FlushCache()  # Write to disk.
    return dataset, dataset.GetRasterBand(1)  #If you need to return, remenber to return  also the dataset because the band don`t live without dataset.

[GIS] Clip raster by mask without change values

The reason your output image is 8-bit is because you explicitly ensure the image is 8-bit the step before you write it.

# Clip the image using the mask
clip = gdalnumeric.numpy.choose(mask, (clip,0)).astype(gdalnumeric.numpy.uint8)

# Save clipping as tiff
gdalnumeric.SaveArray(clip, "%s.tif" % output, format="GTiff", prototype=raster)

Best Answer

Related Solutions

[GIS] Python GDAL: Save array as raster with projection from other file

[GIS] Clip raster by mask without change values

Related Question