2

(This question was asked on another forum 3 months ago but received no response. I'd really welcome some help getting over the problem)

My input map is a JPG scan of a not-to-scale map. I would like to warp it to EPSG:3857 so that it can forced to scale. The result is not warping as I'd expect. Se images below.

The approach taken is to add 9 widely dispersed GCPs (ground control points) but no projection to a new TIF file. The GCPs have been added as X, Y, Longitide, Latitude, Z=0. gdal_edit does this step (from Python) without error and gdalinfo in.tif now shows:

 Size is 3421, 1925 Coordinate System is `' GCP Projection = GCP[ 0]: Id=1, Info= (462,111) -> (-1.541,50.759,0) GCP[ 1]: Id=2, Info= (732,657) -> (-1.496,50.704,0) GCP[ 2]: Id=3, Info= (1998,1764) -> (-1.285,50.588,0) GCP[ 3]: Id=4, Info= (1911,687) -> (-1.291,50.701,0) GCP[ 4]: Id=5, Info= (1896,108) -> (-1.3,50.762,0) GCP[ 5]: Id=6, Info= (2691,360) -> (-1.158,50.732,0) GCP[ 6]: Id=7, Info= (2739,1098) -> (-1.149,50.659,0) GCP[ 7]: Id=8, Info= (3108,756) -> (-1.084,50.689,0) GCP[ 8]: Id=9, Info= (2163,1320) -> (-1.254,50.635,0) 

The subsequent warp command is supposed to introduce a spatial reference system:

gdalwarp -of GTiff -t_srs EPSG:3857 in.tif out-warped-3857.tif 

It runs, but the resulting map is stretched wrongly, from 3421x1925 to 3852x1381. Changing the output projection to EPSG:4326, the result is identical which I wouldn't expect.

There are three images to illustrate:

Reference map (OpenTopoMap) enter image description here

Input map enter image description here

Output map enter image description here

Questions: Is this the best way to add projection and force a map to scale? Why is the vertical dimension shrinking so much?

2
  • If you want to change the projection you must tell what is the original coordinate system in which the ground control points are. gdalwarp -of GTiff -s_srs epsg:4326 -t_srs EPSG:3857 in.tif out-warped-3857.tif should make difference. Commented Jan 4, 2021 at 20:35
  • @user30184. It seems to work - thank you. The GCPs came from Google Earth (epsg:4326) but, given that the source map was not projected as epsg:4326, I didn't realise I needed to supply that as a parameter. Commented Jan 5, 2021 at 9:36

1 Answer 1

2

When ground control points are added into a GeoTIFF with GDAL there is no information about the coordinate system yet. It can be noticed from the gdalinfo report that shows the Coordinate System as empty.

Coordinate System is `' GCP Projection = GCP[ 0]: Id=1, Info= (462,111) -> (-1.541,50.759,0) 

What is known is just that pixel (462,111) is located in coordinates (-1.541,50.759,0) in some coordinate system.

Successful re-projection with gdalwarp requires that both the source CRS (-s_srs) and target CRS (-t_srs) are known. If CRS is stored into the source data gdalwarp reads it automatically and therefore the -s_srs parameter is not mandatory. If source data has no CRS and user does not give -s_srs then it is assumed that -s_srs=-t_srs. That's why the result with only -t_srs parameter looked identical for both EPSG:4326 and EPSG:3857.

The gdalwarp command for converting a GeoTIFF with ground control points measured in EPSG:4326 into normal georeferenced image in EPSG:3857 is:

gdalwarp -of GTiff -s_srs epsg:4326 -t_srs EPSG:3857 in.tif out-warped-3857.tif 

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.