1

i have some trouble warping one raster to different resolution. My input is a one channel tiff with. I know the geo-coordinates of x_min and y_max. I also know that the pixel resolution is 5 meters. I now want to warp the same raster to pixel resolution of 6.5 meters.

Input in the function is the numpy array and the upper left corner coordinates with the src and dst pixel resolution.

My code is:

def resample_raster(raster, raster_x0, raster_ymax, src_xres, src_yres, dst_xres, dst_yres):


    dst_height=int((dst_yres/src_yres)*raster.shape[0])
    dst_width=int((dst_xres/src_xres)*raster.shape[1])

    dstAR=np.empty((dst_height, dst_width))

    srcDS=dstDS=gdal_array.OpenArray(raster)
    srcDS.SetGeoTransform((raster_x0,src_xres,0, raster_ymax, 0, -src_yres ) )
    dstDS=gdal_array.OpenArray(dstAR)
    srcDS.SetGeoTransform((raster_x0,dst_xres,0, raster_ymax, 0, -dst_yres ) )

    gdal_warp_options=osgeo.gdal.WarpOptions(xRes=dst_xres,yRes=dst_yres, resampleAlg=osgeo.gdal.GRA_NearestNeighbour)



    osgeo.gdal.Warp(dstDS, srcDS,options=gdal_warp_options)

    new_raster=np.array(dstDS.GetRasterBand(1).ReadAsArray())

    return new_raster

The function gives me an empty array or is killed.

3 Answers 3

1

I am not certain of what you try to achieve but I have two sugestions:

  1. Follow this idea on how to resize a numpy array using scipy.
  2. Use gdal.Warp to warp a tiff to a different resolution. For example:
import gdal
infn = '/path/to/source.tif'
outfn = '/path/to/target.tif'

xres=6.5
yres=6.5
resample_alg = 'near'

ds = gdal.Warp(outfn, infn, warpoptions=dict(xRes=xres, yRes=yres, resampleAlg=resample_alg)
ds = None
0

I did this with the command line gdalwarp from my python script as follows:

# Build command 
command = f'gdalwarp -t_srs EPSG:{targetEPSG} -r cubic -tr {targetResolution} {targetResolution} -te {leftBound} {bottomBound} {rightBound} {topBound} "{inputRasterPath}" "{newPath}"'
        
# Reproject vector 
os.system(command)
0

I have downsampled a worldview-3 image (GeoTIFF format) to a WorldView-2 resolution, using gdal/osr Linux package.

First, I implement a special line to take the resolution of a reference image (image with lower resolution):

X_RES = gdalinfo -mm $PATH_TO_REFERENCE_IMAGE_WITH_LOWER_RESOLUTION | grep "Pixel Size" | grep "\([0-9]\.[0-9]\+\)\,\([-w][0-9]\.[0-9]\+\)"

Than, applied the gdalwarp routine to each image in a input folder. The complete solution:

PATH_TO_HIGHER_RESOLUTION_IMAGES=$1
PATH_TO_REFERENCE_IMAGE_WITH_LOWER_RESOLUTION=$2
PATH_TO_OUTPUT=$3

RESAMPLE_METHOD="cubic"
X_RES="0.5"

# TODO: get pixel size
#X_RES = gdalinfo -mm $PATH_TO_REFERENCE_IMAGE_WITH_LOWER_RESOLUTION | grep "Pixel Size" | grep "\([0-9]\.[0-9]\+\)\,\([-w][0-9]\.[0-9]\+\)"
#echo $X_RES

if [ -z "$X_RES" ]; then
    echo "\$X_RES resolution is empty"
fi

for image in $PATH_TO_HIGHER_RESOLUTION_IMAGES*.TIF; do
    filename=$(basename $image)

    if [[ ! -e $PATH_TO_OUTPUT ]]; then
        mkdir $PATH_TO_OUTPUT
    elif [[ ! -d $PATH_TO_OUTPUT ]]; then
        echo "$PATH_TO_OUTPUT already exists but is not a directory" 1>&2
    fi

    result_image=$PATH_TO_OUTPUT$filename

    gdalwarp -of GTiff -wo NUM_THREADS=ALL_CPUS -co TILED=YES -co TILED=YES -co COMPRESS=LZW -tr $X_RES $X_RES -r $RESAMPLE_METHOD $image $result_image
done

and usage:

./resample.sh ../folder-with-higher-resolution/ ../folder-with-lower-resolution/REFERENCE_IMAGE.TIF ../output-resampled-images/

** Observe the -tr and -ts flags. First one is for pixel size, where another is for desired image dimension.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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