[gdal-dev] Warping with GCPs at high latitudes

[Knut-Frode Dagestad]

On 24. jan. 2013 20:43, Even Rouault wrote:
> 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

I tried it now with Python API, and it works perfectly, all the way to 
the pole.

Thus the conclusion should be that GCPs can have any spatial reference 
system, but preferably one with no singularity within or near the 
destination area.

I used code like this:
# Make tranformer from GCP SRS to destination SRS
dstSRS = osr.SpatialReference()
dstSRS.ImportFromProj4(srsString)
srcSRS = osr.SpatialReference()
srcGCPProjection = srcDataset.GetGCPProjection()
srcSRS.ImportFromWkt(srcGCPProjection)
transformer = osr.CoordinateTransformation(srcSRS, dstSRS)

# Reproject all GCPs
srcGCPs = self.vrt.dataset.GetGCPs()
dstGCPs = []
for srcGCP in srcGCPs:
	(x, y, z) = transformer.TransformPoint(srcGCP.GCPX, srcGCP.GCPY, 
srcGCP.GCPZ)
	dstGCP = gdal.GCP(x, y, z, srcGCP.GCPPixel, srcGCP.GCPLine, 
srcGCP.Info, srcGCP.Id)
	dstGCPs.append(dstGCP)

# Update dataset
self.vrt.dataset.SetGCPs(dstGCPs, dstSRS.ExportToWkt())