Coordinate System – Convert ETRS89 Portugal TM06 Coordinate to Google Maps Coordinate

coordinate systemepsggoogle maps

I am trying convert an coordinate (point = {-26841,7714704367;42953,3965652831}) from shapefile with follow prj file: 

PROJCS["ETRS89_Portugal_TM06",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137,298.257222101]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",39.66825833333333],PARAMETER["central_meridian",-8.133108333333334],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]

To a Google Maps coordenate. I am reading about it and I found that I need to transform this coordenate into EPSG 3857. I am using the follow code to do that:

const string latlongwkt = "PROJCS[\"Popular Visualisation CRS / Mercator\", GEOGCS[\"Popular Visualisation CRS\", DATUM[\"WGS84\", SPHEROID[\"WGS84\", 6378137.0, 298.257223563, AUTHORITY[\"EPSG\",\"7059\"]], AUTHORITY[\"EPSG\",\"6055\"]], PRIMEM[\"Greenwich\", 0, AUTHORITY[\"EPSG\", \"8901\"]], UNIT[\"degree\", 0.0174532925199433, AUTHORITY[\"EPSG\", \"9102\"]], AXIS[\"E\", EAST], AXIS[\"N\", NORTH], AUTHORITY[\"EPSG\",\"4055\"]], PROJECTION[\"Mercator\"], PARAMETER[\"semi_minor\",6378137], PARAMETER[\"False_Easting\", 0], PARAMETER[\"False_Northing\", 0], PARAMETER[\"Central_Meridian\", 0], PARAMETER[\"Latitude_of_origin\", 0], UNIT[\"metre\", 1, AUTHORITY[\"EPSG\", \"9001\"]], AXIS[\"East\", EAST], AXIS[\"North\", NORTH], AUTHORITY[\"EPSG\",\"3785\"]]";
IProjectedCoordinateSystem googleProjectedCoordinates = CoordinateSystemWktReader.Parse(latlongwkt) as IProjectedCoordinateSystem;

var wktstring = "PROJCS[\"ETRS89_Portugal_TM06\",GEOGCS[\"GCS_ETRS_1989\",DATUM[\"D_ETRS_1989\",SPHEROID[\"GRS_1980\",6378137,298.257222101]],PRIMEM[\"Greenwich\",0],UNIT[\"Degree\",0.017453292519943295]],PROJECTION[\"Transverse_Mercator\"],PARAMETER[\"latitude_of_origin\",39.66825833333333],PARAMETER[\"central_meridian\",-8.133108333333334],PARAMETER[\"scale_factor\",1],PARAMETER[\"false_easting\",0],PARAMETER[\"false_northing\",0],UNIT[\"Meter\",1]]";
IProjectedCoordinateSystem sourceCS = CoordinateSystemWktReader.Parse(wktstring) as IProjectedCoordinateSystem;

CoordinateTransformationFactory ctFactory = new CoordinateTransformationFactory();
ICoordinateTransformation transformer1 = ctFactory.CreateFromCoordinateSystems(sourceCS, googleProjectedCoordinates);

double[] toPoint = transformer1.MathTransform.Transform(fromPoint);

But it still give me a wrong coordinate.

I also tried to replace IProjectedCoordinateSystem googleProjectedCoordinates to:ICoordinateSystem geographicCoordinates = ProjNet.CoordinateSystems.GeographicCoordinateSystem.WGS84

And I got and lat/log coordinate but with wrong location in Google Map.

I already spent a lot of time and I'm stuck on it. Could anyone give me some help? I really can't figure out what is wrong and why it didn't work.

I am using the Proj.Net (http://projnet.codeplex.com/) to apply the projection.

I am using the wkt (that I think is compatible with Google Maps coordinates, correct me if I am wrong) provided in this article (https://alastaira.wordpress.com/2011/01/23/the-google-maps-bing-maps-spherical-mercator-projection/)

I made some changes in code, I changed proj.net to DotSpatial that it is based in Proj.4 and it was recently updated (Thu Apr 28, 2016 at 8:00 AM).

So the code is the follow:

Shapefile nycBoroughs = Shapefile.OpenFile(filename);
nycBoroughs.Reproject(KnownCoordinateSystems.Geographic.Europe.ETRS1989);

and now I got the coordinates as you said:
coordenate = {(-8,44768822459892, 40,0546873641873)}

I also apply the formulas and I got the coordinates also as you said:
x = -940392.35141178139, y = 4873892.4828391289

But now I can't use this coordinates directly in Google Maps, I need to convert again to ETRS1989 to get the coordinates in degree and use it with google maps?

Sorry, probably it is a simple question, but I am a programmer and i don't have enough know-how in system coordinates to understand it easily.

Best Answer

The Proj.Net homepage notes that the code is rather old, from 2009. The development of Proj.4 has improved since then, including proper use of the Google Mercator projection.

If you can't use the Proj.4/GDAL libraries, I suggest to reproject from TM06 to ETRS89/WGS84 coordinates (should be around (-8.4477 40.0547)),then convert to Google mercator using the formulas in the linked article:

double x = lon * 20037508.34 / 180;
double y = Math.Log(Math.Tan((90 + lat) * Math.PI / 360)) / (Math.PI / 180) * 20037508.34 / 180;

Note that Math.Log is the logarithmus naturalis, not decimalis.

The result should be (-940392 4873892).

Related Question