0

I'm trying to write a method (in VB.NET) to transform datum from EPSG:3785 to WGS84. I'm currently making use of ProjNet in my solution so I used it to perform this transformation. I thought everything was going swimmingly at first; I plugged in my data and got my transformed data points out. When I checked the values the results were not exactly what I hoped.

To check the results I used QGIS to compare the data sets (before and after the ProjNet transformation) and noticed that there is a significant error term appearing in the y-axis.

Valid XHTML.

In terms of what my code does: firstly there isn't a lot of documentation for Proj.Net but there was some helpful posts I used here:

  1. How to perform transformations using ProjNet?
  2. https://projnet.codeplex.com/wikipage?title=CreateProjection&referringTitle=FAQ

Similar to the other post, my call to the ProjNet transformation looks like so:

    Dim csf = new ProjNet.CoordinateSystems.CoordinateSystemFactory()
    Dim ctf = new ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory()

    Dim sys3785 as ICoordinateSystem = csf.CreateFromWkt(_EPSG3785)
    Dim sys84 as ICoordinateSystem = csf.CreateFromWkt(_WGS84)

    Dim transformToWgs84 As ICoordinateTransformation = ctf.CreateFromCoordinateSystems(sys3785, sys84)

    Dim epsg3785Points = New Double() {_x, _y}
    Dim wgs84Points As Double() = transformToWgs84.MathTransform.Transform(epsg3785Points)

Attributes _x and _y are data points in EPSG3785. _EPSG3785 and _WGS84 are constants defined from the spatialreference.org resources as:

 Private Const _EPSG3785 = "PROJCS[""Popular Visualisation CRS / Mercator"", GEOGCS[""Popular Visualisation CRS"", DATUM[""D_Popular_Visualisation_Datum"", SPHEROID[""Popular_Visualisation_Sphere"",6378137,0]], PRIMEM[""Greenwich"",0], UNIT[""Degree"",0.017453292519943295]],PROJECTION[""Mercator""],PARAMETER[""central_me idian"",0], PARAMETER[""scale_factor"",1],PARAMETER[""false_easting"",0], PARAMETER[""false_northing"",0], PARAMETER[""latitude_of_origin"",0], UNIT[""Meter"",1]]"
 Private Const _WGS84 = "GEOGCS[""GCS_WGS_1984"",DATUM[""D_WGS_1984"",SPHEROID[""WGS_1984"",6378137,298.257223563]],PRIMEM[""Greenwich"",0],UNIT[""Degree"",0.0174532925199433]]"

Now ProjNet was throwing an error at first with the spacial reference WKT because EPSG:3785 didn't have defined "latitude_of_origin". I'm not too sure why so I appended it with a value of 0 as shown above. That's pretty much as far as I got with it before I got this transformation error.

Out of interest I decided to write a little test to see how the y-values differ between QGIS and the ProjNet transformations. I found that the difference across the y-axis was fairly linear the further the point moves from the origin (x-axis was always near exact). Valid XHTML.

Now I could write a little correction factor to the y-axis to re-align the data but that doesn't seem to solve the issue here. To me the problem exists either as an error in the projection strings or within ProjNet's transformation.

Is there anything apparent as to why this transformation error is appearing?

  • Are you sure that the other projection is EPSG:3785 and not 3857? – user30184 Sep 12 '16 at 18:54
  • Can you give a sample point? or try x=0, y=5000000. For EPSG:3857, lat=40.9162744707. For EPSG:3395, lat=41.1067450395. 3857 and 3785 should give identical results. Most issues with y offsets and "web Mercator" are due to ellipsoidal/spherical equations being used. 3395 is on an ellipsoid. – mkennedy Sep 12 '16 at 22:45
  • You might also take a look at alastaira.wordpress.com/2011/01/23/…. Adding the EPSG code in the last line of the WKT should solve the issue. – AndreJ Sep 13 '16 at 5:51
  • Reason I believe the projection is in 3785 is due to the data alignment shown in the QGIS image. The "datums" layer when loaded as 3785 lies on the network where expected (the network is in WGS84). I tried using 3857 and I get fractionally better results in the ProjNet transformation but there still exists an error on the y-axis. QGIS is projecting the independent layers correctly but when I replicate it in ProjNet I get this error. – GetReQ Sep 13 '16 at 8:52
  • @mkennedy Sure: the item shown in the OP is read in as {849289.734113, 1007603.928995}. ProjNet's reprojection into WGS84 gives {7.629299488, 9.073781567} whereas QGIS's export into WGS84 gives {7.629299488, 9.014043796}. – GetReQ Sep 13 '16 at 9:00
1

Thanks to @AndreJ for his reference, it was an interesting read and it contained a reference to solving this projection issue:

http://projnet.codeplex.com/discussions/83818

As suggested by tanoshimi's post in the link above, I added the parameter to the EPSG:3785 WKT and got a far more accurate WGS84 projection through ProjNet:

 PARAMETER["semi_minor",6378137]

The issue does sit within ProjNets transformation as discussed within it's CodePlex board.

Did the same y-comparison between ProjNet and QGIs and got a far more consistent transformation through ProjNet.

Valid XHTML.

This is now far more visually apparent (green dot compared with yellow).

Valid XHTML.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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