CoordinateSharp icon indicating copy to clipboard operation
CoordinateSharp copied to clipboard

Add Support For EPSG:3857 WGS 84 / Pseudo-Mercator

Open Reag opened this issue 2 years ago • 1 comments

I work with geo-data that is in UTM 36N, but I need to sync this data to a map that uses EPSG:3857 (Web Mercator projection). EPSG:3857 is the projection that backs Google and Bing Maps.

Would it be possible to add this format for conversion in the library? I already use CoordinateSharp for latlon conversion to ECEF in another part of the project, and having a single lib for all of this work would be very helpful!

For reference, the model contains data in the format ~(72800.000, 3650000.000) and I need to convert it to ~(3945000.00, 3892000.00)

Reag avatar Jul 14 '22 08:07 Reag

I think this is doable. We will investigate feasibility and performance impacts.

References: https://www.iogp.org/wp-content/uploads/2019/09/373-07-02.pdf https://lyzidiamond.github.io/posts/4326-vs-3857

Tronald avatar Jul 16 '22 21:07 Tronald

EXAMPLE PROGRAM WRITTEN FOR IMPLEMENTATION REFERENCE

Credit Shan Siddiqui for assistance with the math.

Using CoordinateSharp;
Using CoordinateSharp.Formatters;

static void Main(string[] args)
{
   To_Web_Mercator();
}

/// <summary>
/// Convert to web mercator. Fails near poles, so logic must be built in to truncate +- 85.06 per standard once implemented.
/// </summary>
static void To_Web_Mercator()
{
            //E = Easting
            //FE = False Easting
            //a= Semi-Major Axis
            //lng = Longitude
            //lngNO = Longitude of Natural Origin
            //N = Northing
            //FN = False Northing
            //ln = NATURAL LOGARITHM

            //Easting = FE + a * (lng - lngNO)
            //Northing = FN + a * ln[tan(pi/4*latitude/2]

            Coordinate c = Coordinate.Parse("24°22'54.433\" 100°20'0\"");

            double FE = 0;
            double FN = 0;
            double a = 6378137.0;

            double E = FE + a * (c.Longitude.ToRadians() - 0);
            double N = FN + a * Math.Log(Math.Tan(Math.PI / 4 + c.Latitude.ToRadians() / 2)); 

            Console.WriteLine($"{E} {N}");

            From_Web_Mercator(E, N);

}

static void From_Web_Mercator(double easting, double northing)
{
            //E = Easting
            //FE = False Easting
            //a= Semi-Major Axis
            //lng = Longitude
            //lngNO = Longitude of Natural Origin
            //N = Northing
            //FN = False Northing

            //D = -(FN - N) / a = (FN-n)/a
            //lat = pi/2 - 2 * atan(e^D) where e=base of natural log 2.7182818
            //lng = [(E-FE)/a] +lngNO

            double N = northing;
            double E = easting;
            double FE = 0;
            double FN = 0;
            
            double a = 6378137.0;
            double lngNO = 0;
            double D = -(N - FN) / a;

            double lat = Math.PI / 2 - 2 * Math.Atan(Math.Pow(Math.E,D));
            double lng = ((E - FE) / a) + lngNO;
            Coordinate c = new Coordinate(Extensions.ToDegrees(lat), Extensions.ToDegrees(lng));
            Console.WriteLine(c);

}

NOTES

Class nomenclature should elude to "Web Mercator" projection with default ESPG CRS code 3857. Update documentation for clarity when implementing (to include existing geodetic system).

Determine if Web Mercator should exist in base Coordinate class, or exist outside with Eager Loading options. Determination should be based on performance.

Translation Comparison Tool: https://epsg.io/transform#s_srs=3857&t_srs=4326&x=11169055.5762583&y=2800000.0031361

Reference Pages 44-45 for conversion math documentation and example: https://www.iogp.org/wp-content/uploads/2019/09/373-07-02.pdf

Tronald avatar Aug 14 '22 23:08 Tronald

STATUS AS OF 8/17/2022

WebMercator Implemented as new class with separate eager load options to handle EPSG:3857 conversions for mapping applications in sandbox. False Easting/Northing values will exist in math formulas, but will not be exposed to users unless requested. Feature if requested will be released with subsequent updates.

TO DO

-Update Parsers -Create Unit Tests -Update Documentation -Update Web Guides -Update Eager Loading Guides -Update Benchmarks -Push to Git -Release BETA to Nuget

Hoping to have this complete in the next few weeks.

Tronald avatar Aug 18 '22 00:08 Tronald

Thank you for working on this! Will make my life a lot easier if I can do all my conversions from a single library

Reag avatar Aug 18 '22 04:08 Reag

Pre-Release version is out if you wish to get a jump and test. There is no online documentation yet, but the operations is identical to other coordinate systems.

Just note the conversions become inaccurate above/below +/- 85.06 Latitude. This is a documented limitation of the system itself. Coordinate.WebMercator will throw a FormatException if called with any datum other than WGS84 (as WebMercator requires WGS84). If you are not changing datums, you shouldn't have issue.

https://www.nuget.org/packages/CoordinateSharp/2.14.1.1-beta

Coordinate c = new Coordinate(45, 64);
Console.WriteLine(c.WebMercator.ToString()); //7124447.411mE 5621521.486mN
		
c = Coordinate.Parse("3624447.411 2721521.486");
Console.WriteLine(c.ToString()); //N 23º 44' 17.003" E 32º 33' 32.274"

Tronald avatar Aug 28 '22 01:08 Tronald

Fantastic! Thank you for your work on this feature!

Reag avatar Aug 28 '22 04:08 Reag

Feature added and complete with v2.14.2.1 published today to nuget.org.

No user feedback from BETA testing period. Minor changes were made to tests, documentation and licensing since pre-release so version has been updated.

Thank you for the recommendation.

Tronald avatar Sep 24 '22 17:09 Tronald