Reproject image with LSAT projection using GDAL

I'm working near a team of data-scientist specialized in satellite imagery. One the data-scientist explain me its problem with Landsat images. Landsat images uses the UTM spatial reference and he wants to reproject all images with the Space oblique for LANDSAT projection. I have never heard about this projection before and I discovered that it is a variation of Space Oblique Mercartor Projection. This projection aims to map all image along the swath. Nasa page quotes that it is "one of the most complex projections ever devised".

Explanation

GDAL is the standard library for projection and the gdalwarp utility is usefull. It's really powerfull but you often have very complex case with unusual usage. Here it's the proj4 definition of lsat that doesn't work.

% gdalwarp -s_srs "EPSG:32630" -t_srs '+proj=lsat +lsat=5 +path=201' 24_mosaic_8bB08.sca.tif toto.tif
ERROR 1: Translating source or target SRS failed:
+proj=lsat +lsat=5 +path=201

According to the maintainer of Gdal, there is a parsing limitation in GDal! And the proj.4 string using lsat can't be converto to WKT format. He advises to uses WKT and extension and after testing it works! Thank Even Rouault to suggest the tip on the osgeo forum.

Final result

function get_srs ()
{
    name=$1
    gdalinfo $1 | grep 'AUTHORITY\["EPSG","326' | cut -d '"' -f 4
}

for f in $(ls *tif) ;
do
    epsg_code=$(get_srs "$f")
    basic_name=$(basename $f)
    gdalwarp -s_srs "EPSG:$epsg_code" -t_srs 'LOCAL_CS["lsat0",EXTENSION["PROJ4","+proj=lsat +lsat=5 +path=201 +wktext"]]' "$f" "reprojected/$f"
done
By @Romain JACQUET in
Tags :