I have a map file in a tiff format that is not geo-referenced, but I know the following projection information: +proj=lcc +lon_0=10.856 +lat_0=51.322 +lat_1=42 +lat_2=57 +datum=WGS84. I would like to transform it into EPSG:3857.
In my first attempt I tried to fix four arbitrary points and perform warping
gdal_translate -a_srs "+proj=lcc +lon_0=10.856 +lat_0=51.322 +lat_1=42 +lat_2=57 +datum=WGS84" -of GTiff -gcp 19725 11865 1 56 -gcp 103755 12990 23 56 -gcp 112755 80400 23 46 -gcp 8925 78990 1 46 original.tif translated.tif
But the coordinates of the corners of the warped image where not as expected and the map was distorted.
Corner Coordinates: Upper Left ( 1208480.629, 6678543.890) ( 10d51'21.48"E, 51d19'21.08"N) Lower Left ( 1208480.629, 6678541.123) ( 10d51'21.48"E, 51d19'21.03"N) Upper Right ( 1208486.180, 6678543.890) ( 10d51'21.66"E, 51d19'21.08"N) Lower Right ( 1208486.180, 6678541.123) ( 10d51'21.66"E, 51d19'21.03"N) Center ( 1208483.404, 6678542.507) ( 10d51'21.57"E, 51d19'21.06"N)
After reading many posts, here and elsewhere, I came to this:
- Converting the UL and LR corners to x and y
cs2cs +proj=latlong +lon_0=10.856 +lat_0=51.322 +lat_1=42 +lat_2=57 +datum=WGS84 +to +init=epsg:3857 -3.63 57.28 > -404089.75 7817564.89 23.76 44.39 > 2644951.10 5525995.28
- Geo-referencing the tiff
gdal_translate -a_srs "+proj=lcc +lon_0=10.856 +lat_0=51.322 +lat_1=42 +lat_2=57 +datum=WGS84" -a_ullr -404089.75 7817564.89 2644951.10 5525995.28 original.tif translated.tif
The result of gdalinfo translated.tif is as follows:
Size is 14735, 11333 Coordinate System is: PROJCS["unnamed", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]], PROJECTION["Lambert_Conformal_Conic_2SP"], PARAMETER["standard_parallel_1",42], PARAMETER["standard_parallel_2",57], PARAMETER["latitude_of_origin",51.322], PARAMETER["central_meridian",10.856], PARAMETER["false_easting",0], PARAMETER["false_northing",0], UNIT["metre",1, AUTHORITY["EPSG","9001"]]] Origin = (-404089.750000000000000,7817564.889999999664724) Pixel Size = (206.925066168985410,-202.203265684284787) Metadata: AREA_OR_POINT=Area TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch) TIFFTAG_XRESOLUTION=144 TIFFTAG_YRESOLUTION=144 Image Structure Metadata: INTERLEAVE=PIXEL Corner Coordinates: Upper Left ( -404089.750, 7817564.890) (146d18'53.31"E, 73d27'53.56"N) Lower Left ( -404089.750, 5525995.280) (158d45'33.92"W, 88d 1'24.64"N) Upper Right ( 2644951.100, 7817564.890) (172d26'19.88"W, 64d26'59.68"N) Lower Right ( 2644951.100, 5525995.280) (138d13'57.20"E, 73d22'12.76"N) Center ( 1120430.675, 6671780.085) (161d52'17.63"W, 79d37'35.75"N)
However, the coordinates of the corners above are completely off. When I attempted to warp it, the image could become so huge that I had to stop the process after a couple of minutes.
gdalwarp -overwrite -t_srs EPSG:3857 -r near -co COMPRESS=LZW translated.tif warped.tif
The result of gdalinfo warped.tif is as follows:
Corner Coordinates: Upper Left (-20037508.339,29704755.766) (180d 0' 0.00"W, 88d54'44.28"N) Lower Left (-20037508.339, 9464937.842) (180d 0' 0.00"W, 64d26'59.58"N) Upper Right (20037503.600,29704755.766) (179d59'59.85"E, 88d54'44.28"N) Lower Right (20037503.600, 9464937.842) (179d59'59.85"E, 64d26'59.58"N) Center ( -2.370,19584846.804) ( 0d 0' 0.08"W, 84d41'15.52"N)
What did I miss?