1

I would like change a raster dataset from EPSG: 3035 (ETRS89/LAEA Europa) to EPSG: 6204 (Macedonia State Coordinate System)

My WKT for the first coordinate system is:

PROJCS["ETRS89 / LAEA Europe",
    GEOGCS["ETRS89",
        DATUM["European_Terrestrial_Reference_System_1989",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            TOWGS84[0,0,0,0,0,0,0],
            AUTHORITY["EPSG","6258"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4258"]],
    PROJECTION["Lambert_Azimuthal_Equal_Area"],
    PARAMETER["latitude_of_center",52],
    PARAMETER["longitude_of_center",10],
    PARAMETER["false_easting",4321000],
    PARAMETER["false_northing",3210000],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","3035"]]

and for the Second SR:

PROJCS["unnamed",
    GEOGCS["Bessel 1841",
        DATUM["unknown",
            SPHEROID["bessel",6377397.155,299.1528128],
            TOWGS84[682,-203,480,0,0,0,0]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",21],
    PARAMETER["scale_factor",0.9999],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["Meter",1],
    AUTHORITY["epsg","6204"]]

Since I want to reproject my data, and the projections have different datums, how this process would look like using gdalwarp?

According to this post (How to reproject a raster file in QGIS with datum transformation?) gdalwarp is the way to go, in the following form:

gdalwarp -s_srs "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs" -t_srs "+proj=tmerc +lat_0=0 +lon_0=21 +k=0.9999 +x_0=500000 +y_0=0 +ellps=bessel +towgs84=682,-203,480,0,0,0,0 +units=m +no_defs" input.tif output.tif

Is that a valid approach? I can't see the datum information in the proj4 strings I am providing for both map projections. Is there anything I have to take care of when projecting a raster datasets from a map projection which uses a different datum to my target Spatial Reference?

3
  • When I run the command line tool from gdal, it throws the error: ERROR1: Too many points (10201 out of 10201) failed to transform, unable to compute output bounds. Dec 18, 2019 at 6:52
  • The +towgs84= parameter is the datum transformation. Both definitions are using WGS84 geodetic datum as a pivot. Dec 18, 2019 at 13:55
  • About the error provided in you previous comment, the geometries of all the features lies outside the valid geodetic coordinates (latitude modules greater than 90 or longitude modules greater than 360) in the transformation. I don't know if there are a bounding box validation in the target projection, but seems like a problem in the input coordinates. Can you provide some point coordinate to test? Dec 18, 2019 at 14:02

0

Your Answer

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

Browse other questions tagged or ask your own question.