1

Full disc Earth

How do I gdal_translate an orthographic PNG image of the Earth centered at Lat=+50N Lon=+5E to proper GeoTIFFf format ?

I cannot work out how to specify the limits of a full disc. (diameter of Earth in metres? 180 degrees sphere? 90 degrees from centre?) Which argument do I add to the Spatial Reference system string? How do I specify GCPs to a full disc if required?

Image Ortho-5_50-clean.png is png 2048px x 2048px full disc Earth with transparency enabled.

gdal_translate -of Gtiff  -a_srs "+proj=ortho +ellps=WGS84 +datum=WGS84 +lat_0=50 +lon_0=5"   Ortho-5_50-clean.png Ortho-5_50-clean.tif

This command to reproject as longlat fails with an incomprehensible error:

gdalwarp -t_srs "+proj=longlat +datum=WGS84 +ellps=WGS84 +units=m" test_ortho.tif test_EQ.tif

ERROR 1: The transformation is already "north up" or a transformation between pixel/line and georeferenced coordinates cannot be computed for Ortho-5_50-clean.tif. There is no affine transformation and no GCPs. Specify transformation option SRC_METHOD=NO_GEOTRANSFORM to bypass this check.


Version -- GDAL 3.4.1, released 2021/12/27
Version -- ImageMagick (convert) 6.9.11-60
Version -- Linux 6.5.0-14-generic #14~22.04.1-Ubuntu


gdalinfo Ortho-5_50-clean.tif :

Driver: GTiff/GeoTIFF
Files: Ortho-5_50-clean.tif
Size is 2048, 2048
Coordinate System is:
PROJCRS["unknown",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["Orthographic",
        METHOD["Orthographic",
            ID["EPSG",9840]],
        PARAMETER["Latitude of natural origin",85,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",-30,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["False easting",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]]]
Data axis to CRS axis mapping: 1,2
  Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 2048.0)
Upper Right ( 2048.0,    0.0)
Lower Right ( 2048.0, 2048.0)
Center      ( 1024.0, 1024.0)
Band 1 Block=2048x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA 
Band 2 Block=2048x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA 
Band 3 Block=2048x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA 
Band 4 Block=2048x1 Type=Byte, ColorInterp=Alpha

2 Answers 2

3

Let's assume a spherical datum. In the orthographic projection the disk measures the diameter of the sphere.

Let's assume a sphere of 6378137 meters in radius, so the georeferenced coordinates of the image vertices are:

  • Upper Left X: -6378137
  • Upper Left Y: 6378137
  • Lower Right X: 6378137
  • Lower Right Y: -6378137

That's all we need, let's create the translate:

$ gdal_translate -of Gtiff -a_srs "+proj=ortho +R=6378137 +lat_0=50 +lon_0=5" -a_ullr -6378137 6378137 6378137 -6378137 Ortho-5_50-clean.png Ortho-5_50-clean.tif

Input file size is 2048, 2048
0...10...20...30...40...50...60...70...80...90...100 - done.

$ gdalinfo Ortho-5_50-clean.tif

Driver: GTiff/GeoTIFF
Files: Ortho-5_50-clean.tif
Size is 2048, 2048
Coordinate System is:
PROJCRS["unknown",
    BASEGEOGCRS["unknown",
        DATUM["unknown",
            ELLIPSOID["unknown",6378137,0,
                LENGTHUNIT["metre",1,
                    ID["EPSG",9001]]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CONVERSION["Orthographic",
        METHOD["Orthographic",
            ID["EPSG",9840]],
        PARAMETER["Latitude of natural origin",50,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",5,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["False easting",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]]]
Data axis to CRS axis mapping: 1,2
Origin = (-6378137.000000000000000,6378137.000000000000000)
Pixel Size = (6228.649414062500000,-6228.649414062500000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
ERROR 1: Point outside of projection domain
Upper Left  (-6378137.000, 6378137.000) 
ERROR 1: Point outside of projection domain
Lower Left  (-6378137.000,-6378137.000) 
ERROR 1: Point outside of projection domain
Upper Right ( 6378137.000, 6378137.000) 
ERROR 1: Point outside of projection domain
Lower Right ( 6378137.000,-6378137.000) 
Center      (   0.0000000,   0.0000000) (  5d 0' 0.00"E, 50d 0' 0.00"N)
Band 1 Block=2048x1 Type=Byte, ColorInterp=Red
Band 2 Block=2048x1 Type=Byte, ColorInterp=Green
Band 3 Block=2048x1 Type=Byte, ColorInterp=Blue
0

To translate a full disc PNG image of the Earth in orthographic or perspective to a proper geotiff you must include the geographic limits in metres... This is the radius of the Earth (~6367444.65 metres). I have no idea how to calculate the precise measurements that will vary depending on the position of the camera. The distances will vary as the equatorial bulge shifts with the latitude of the view. minimum at Lat=0, max at Lat=90.

My example command line assumes a spherical Earth:

gdal_translate -of Gtiff -a_ullr -6367444.65 6367444.65 6367444.65 -6367444.65 -a_srs "+proj=ortho +ellps=WGS84 +datum=WGS84 +lat_0=50 +lon_0=5" source_disc.png target_disc.tif

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.