Updated MapProjections

This commit is contained in:
Clemens 2022-01-21 00:17:01 +01:00
parent c6c597a142
commit 544ac32ee3
5 changed files with 129 additions and 112 deletions

View file

@ -30,8 +30,24 @@ namespace MapControl.Projections
private double scaleFactor;
private string bboxFormat;
public IMathTransform LocationToMapTransform { get; private set; }
public IMathTransform MapToLocationTransform { get; private set; }
public GeoApiProjection(string wkt = null)
{
if (wkt != null)
{
WKT = wkt;
}
}
/// <summary>
/// Gets or sets an OGC Well-known text representation of a coordinate system,
/// i.e. a PROJCS[...] or GEOGCS[...] string as used by https://epsg.io or http://spatialreference.org.
/// Setting this property updates the CoordinateSystem property with an ICoordinateSystem created from the WKT string.
/// </summary>
public string WKT
{
get { return CoordinateSystem?.WKT; }
set { CoordinateSystem = new CoordinateSystemFactory().CreateFromWkt(value); }
}
/// <summary>
/// Gets or sets the ICoordinateSystem of the MapProjection.
@ -85,16 +101,9 @@ namespace MapControl.Projections
}
}
/// <summary>
/// Gets or sets an OGC Well-known text representation of a coordinate system,
/// i.e. a PROJCS[...] or GEOGCS[...] string as used by https://epsg.io or http://spatialreference.org.
/// Setting this property updates the CoordinateSystem property with an ICoordinateSystem created from the WKT string.
/// </summary>
public string WKT
{
get { return CoordinateSystem?.WKT; }
set { CoordinateSystem = new CoordinateSystemFactory().CreateFromWkt(value); }
}
public IMathTransform LocationToMapTransform { get; private set; }
public IMathTransform MapToLocationTransform { get; private set; }
public override bool IsNormalCylindrical
{
@ -106,15 +115,6 @@ namespace MapControl.Projections
get { return isWebMercator; }
}
public GeoApiProjection(string wkt = null)
{
if (wkt != null)
{
WKT = wkt;
}
}
public override Point LocationToMap(Location location)
{
if (LocationToMapTransform == null)

View file

@ -11,40 +11,90 @@ namespace MapControl.Projections
{
public class GeoApiProjectionFactory : MapProjectionFactory
{
public const int WorldMercator = 3395;
public const int WebMercator = 3857;
public const int Etrs89UtmNorthFirst = 25828;
public const int Etrs89UtmNorthLast = 25838;
public const int Wgs84UtmNorthFirst = 32601;
public const int Wgs84UtmNorthLast = 32660;
public const int Wgs84UpsNorth = 32661;
public const int Wgs84UtmSouthFirst = 32701;
public const int Wgs84UtmSouthLast = 32760;
public const int Wgs84UpsSouth = 32761;
private readonly Dictionary<int, string> wkts = new Dictionary<int, string>();
private readonly HttpClient httpClient = new HttpClient();
public override MapProjection CreateProjection(string projectionDefinition)
public override MapProjection CreateProjection(string crsId)
{
var projection = base.CreateProjection(projectionDefinition);
MapProjection projection = null;
if (projection == null &&
projectionDefinition.StartsWith("EPSG:") &&
int.TryParse(projectionDefinition.Substring(5), out int epsgCode))
if (crsId.StartsWith("EPSG:") && int.TryParse(crsId.Substring(5), out int epsgCode))
{
if (epsgCode >= 32601 && epsgCode <= 32660)
switch (epsgCode)
{
projection = new UtmProjection(epsgCode - 32600, true);
}
else if (epsgCode == 32661)
{
projection = new UpsNorthProjection();
}
else if (epsgCode >= 32701 && epsgCode <= 32760)
{
projection = new UtmProjection(epsgCode - 32700, false);
}
else if (epsgCode == 32761)
{
projection = new UpsSouthProjection();
}
else
{
projection = new GeoApiProjection(GetWkt(epsgCode));
case WorldMercator:
projection = new WorldMercatorProjection();
break;
case WebMercator:
projection = new WebMercatorProjection();
break;
case int c when c >= Etrs89UtmNorthFirst && c <= Etrs89UtmNorthLast:
projection = new GeoApiProjection(GetEtrs89UtmWkt(epsgCode));
break;
case int c when c >= Wgs84UtmNorthFirst && c <= Wgs84UtmNorthLast:
projection = new UtmProjection(epsgCode - Wgs84UtmNorthFirst + 1, true);
break;
case int c when c >= Wgs84UtmSouthFirst && c <= Wgs84UtmSouthLast:
projection = new UtmProjection(epsgCode - Wgs84UtmSouthFirst + 1, false);
break;
case Wgs84UpsNorth:
projection = new UpsNorthProjection();
break;
case Wgs84UpsSouth:
projection = new UpsSouthProjection();
break;
default:
projection = new GeoApiProjection(GetWkt(epsgCode));
break;
}
}
return projection;
return projection ?? base.CreateProjection(crsId);
}
private static string GetEtrs89UtmWkt(int epsgCode)
{
const string etrs89UtmWktFormat
= "PROJCS[\"ETRS89 / UTM zone {1}N\","
+ "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[\"Transverse_Mercator\"],"
+ "PARAMETER[\"latitude_of_origin\",0],"
+ "PARAMETER[\"central_meridian\",{2}],"
+ "PARAMETER[\"scale_factor\",0.9996],"
+ "PARAMETER[\"false_easting\",500000],"
+ "PARAMETER[\"false_northing\",0],"
+ "UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],"
+ "AXIS[\"Easting\",EAST],"
+ "AXIS[\"Northing\",NORTH],"
+ "AUTHORITY[\"EPSG\",\"{0}\"]]";
int centralMeridian = 6 * (epsgCode - Etrs89UtmNorthFirst) - 15;
return string.Format(etrs89UtmWktFormat, epsgCode, epsgCode - 25800, centralMeridian);
}
private string GetWkt(int epsgCode)

View file

@ -16,15 +16,15 @@ namespace MapControl.Projections
/// </summary>
public class PolarStereographicProjection : MapProjection
{
public static double ConvergenceTolerance = 1e-6;
public static int MaxIterations = 10;
public static double ConvergenceTolerance { get; set; } = 1e-6;
public static int MaxIterations { get; set; } = 10;
private readonly bool north;
private readonly double scaleFactor;
private readonly double falseEasting;
private readonly double falseNorthing;
public PolarStereographicProjection(string crsId, bool north, double scaleFactor = 1d, double falseEasting = 0d, double falseNorthing = 0d)
public PolarStereographicProjection(string crsId, bool north, double scaleFactor, double falseEasting, double falseNorthing)
{
CrsId = crsId;
this.north = north;

View file

@ -2,47 +2,20 @@
// © 2022 Clemens Fischer
// Licensed under the Microsoft Public License (Ms-PL)
using System;
using ProjNet.CoordinateSystems;
using System;
namespace MapControl.Projections
{
public class UtmProjection : GeoApiProjection
{
private string zone;
public string Zone
{
get { return zone; }
set
{
if (zone != value)
{
if (string.IsNullOrEmpty(value))
{
throw new ArgumentException("Invalid UTM zone.");
}
var hemisphere = value[value.Length - 1];
if ((hemisphere != 'N' && hemisphere != 'S') ||
!int.TryParse(value.Substring(0, value.Length - 1), out int zoneNumber))
{
throw new ArgumentException("Invalid UTM zone.");
}
SetZone(zoneNumber, hemisphere == 'N');
}
}
}
public UtmProjection()
{
}
public UtmProjection(int zoneNumber, bool north)
public UtmProjection(int zone, bool north)
{
SetZone(zoneNumber, north);
SetZone(zone, north);
}
public UtmProjection(Location location)
@ -50,27 +23,21 @@ namespace MapControl.Projections
SetZone(location);
}
public void SetZone(int zoneNumber, bool north)
{
if (zoneNumber < 1 || zoneNumber > 60)
{
throw new ArgumentException("Invalid UTM zone number.", nameof(zoneNumber));
}
var zoneName = zoneNumber.ToString() + (north ? "N" : "S");
if (zone != zoneName)
{
zone = zoneName;
CoordinateSystem = ProjectedCoordinateSystem.WGS84_UTM(zoneNumber, north);
}
}
public void SetZone(Location location)
{
var zoneNumber = Math.Min((int)(Location.NormalizeLongitude(location.Longitude) + 180d) / 6 + 1, 60);
SetZone(zoneNumber, location.Latitude >= 0d);
}
public void SetZone(int zone, bool north)
{
if (zone < 1 || zone > 60)
{
throw new ArgumentException("Invalid UTM zone number.", nameof(zone));
}
CoordinateSystem = ProjectedCoordinateSystem.WGS84_UTM(zone, north);
}
}
}

View file

@ -17,24 +17,24 @@ namespace MapControl.Projections
{
public WorldMercatorProjection()
{
WKT = "PROJCS[\"WGS 84 / World Mercator\"," +
"GEOGCS[\"WGS 84\"," +
"DATUM[\"WGS_1984\"," +
"SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]]," +
"AUTHORITY[\"EPSG\",\"6326\"]]," +
"PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]]," +
"UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]]," +
"AUTHORITY[\"EPSG\",\"4326\"]]," +
"PROJECTION[\"Mercator_1SP\"]," +
"PARAMETER[\"latitude_of_origin\",0]," +
"PARAMETER[\"central_meridian\",0]," +
"PARAMETER[\"scale_factor\",1]," +
"PARAMETER[\"false_easting\",0]," +
"PARAMETER[\"false_northing\",0]," +
"UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]]," +
"AXIS[\"Easting\",EAST]," +
"AXIS[\"Northing\",NORTH]," +
"AUTHORITY[\"EPSG\",\"3395\"]]";
WKT = "PROJCS[\"WGS 84 / World Mercator\","
+ "GEOGCS[\"WGS 84\","
+ "DATUM[\"WGS_1984\","
+ "SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],"
+ "AUTHORITY[\"EPSG\",\"6326\"]],"
+ "PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],"
+ "UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],"
+ "AUTHORITY[\"EPSG\",\"4326\"]],"
+ "PROJECTION[\"Mercator_1SP\"],"
+ "PARAMETER[\"latitude_of_origin\",0],"
+ "PARAMETER[\"central_meridian\",0],"
+ "PARAMETER[\"scale_factor\",1],"
+ "PARAMETER[\"false_easting\",0],"
+ "PARAMETER[\"false_northing\",0],"
+ "UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],"
+ "AXIS[\"Easting\",EAST],"
+ "AXIS[\"Northing\",NORTH],"
+ "AUTHORITY[\"EPSG\",\"3395\"]]";
}
public override Vector GetRelativeScale(Location location)