[gdal-dev] Warping with GCPs at high latitudes

[Even Rouault]

Le jeudi 24 janvier 2013 13:18:50, Knut-Frode Dagestad a écrit :
> Hi,
> 
> I am having trouble to perform accurate warping of images with GCPs
> (lon,lat) in the high Arctic (to "+proj=stere +lat_0=90 +lon_0=0");
> problems starting already around 80 degrees of latitude.
> Thin Plate Spline (-tps) performs better than polynomial fits, but
> output is often severely dragged at one or two edges, and sometimes it
> looks like the destination image is mirrored across a line.
> Trimming images at the edges and adjusting warp options SAMPLE_STEPS and
> SAMPLE_GRID helps a little, but not sufficiently.
> 
> For the polynomial fits I assume the problem is related to precision due
> to the fact that lon and lat become singular at the pole. However, it
> should then probably be better if the algorithm could:
> - convert first GCPs to the output projection (in my case stereographic,
> which has no singularity for the relevant images)
> - make the polynomial fit GCP_stereo <-> pixel, line
> - use this polynomial for the warping
> 
> I guess most GDAL users have a long wishlist of switches, but how about
> a switch to control whether GCP polynomial should rather be made in
> destination projection? Or a switch to gdal_translate to change the
> GCPProjection and GCP coordinates. A simpler and more pragmatic solution
> could be to make a separate Python script gdal_translate_gcps.py?
> 
> Are there other suggestion to deal with these problems?

In your use case, the source image is more or less very close to its final 
shape in Stereographics projection (the warping doesn't warp much the image), 
right ?

Did you try to do your suggestion at hand to actually check that your idea 
leads to the expected result ? This is basically a matter of running :

1) gdaltransform -s_srs EPSG:4326 -t_srs "+proj=stere +lat_0=90 +lon_0=0" on 
the coordinates of the GCPs of the source image
2) gdal_translate -a_srs  "+proj=stere +lat_0=90 +lon_0=0"  [-gcp x y X Y]* 
src.tif src_reprojected_gcp.tif where x y come from the source image and X Y 
are the reprojected coordinates
3) gdalwarp src_reprojected_gcp.tif out.tif

> 
> 
> Best regards from Knut-Frode
> 
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev