Different CRS and their impact on how you work with maps

[JosefWN]

Background

I made this ticket to try to explain, to the best of my abilities, how I use flutter_map and some of the challenges it has posed. I probably should have made it sooner but I knew it would be a lengthy explanation and I've been a bit busy the past year.

To me, choosing which CRS to use is similar to choosing which programming language to use; what's downright impossible in one CRS could be a breeze in another. If you switch from C to Dart, your first thought likely wouldn't be to look into disabling the garbage collector, but rather how to leverage garbage collection and other language features to the greatest extent. Similarly, in leveraging the strengths and avoiding the pitfalls of the CRS you have chosen, you will work with your coordinates differently than you would if you had chosen another CRS.

Since there aren't that many families of projections to choose from however, you don't necessarily have the luxury of choice. So when you hear people complaining about seemingly exotic CRS, it's not likely that they wanted to work with them, but simply that they had to, which I will explain further below.

Spurred by: https://github.com/fleaflet/flutter_map/pull/1750 Relates to: https://github.com/fleaflet/flutter_map/issues/1386, https://github.com/fleaflet/flutter_map/pull/1295, https://github.com/fleaflet/flutter_map/issues/191 / https://github.com/fleaflet/flutter_map/issues/211

Technical background

To be able to make sense of latitudes and longitudes they need to belong to some geometry resembling that of the Earth's. This geometry, called a geoid, needs a spherical coordinate system defining the origin of latitude/longitude and how the axes are oriented. The geoid, along with the definition of the coordinate system, is collectively called a datum.

As the Earth formed it was hotter and softer than it is today. Rotating around its own axis, the planets mass was pulled toward the equator, flattening the regions around the poles, like a spinning pizza dough. Part of WGS84 is a datum that defines the coordinate system of the Earth on an "oblate spheroid", a flattened sphere more generally referred to as an ellipsoid. The WGS84 datum is what GNSS/GPS systems and most projections use to understand how they should relate to the latitudes and longitudes.

This alignment on the WGS84 datum is convenient, because it ensures that when you give me a latitude and a longitude, it means exactly the same thing to you as it does to me, regardless of CRS. You could use a spherical datum/approximation to speed up calculations because of the symmetry (at the cost of precision), but the projection of the source data is pretty much always derived from the WGS84 datum.

The second, trickier part, is how we unfold this ellipsoid so it forms a 2-dimensional map, especially if we want a rectangular map, as we're used to thinking of 2D maps. WGS84 not only defines a datum, but also a projection from that datum onto a plane.

The projection is called Equirectangular projection, although in the odd case that the projection is omitted by the CRS, the projection library can choose a suitable projection (typically Equirectangular projection or Mercator).

The combination of the WGS84 datum and the projection forms a CRS which is confusingly also referred to as WGS84, or EPSG:4326. Common CRS and their building blocks, such as the datum, are identified by EPSG codes. This modular system makes it possible to use or even roll your own CRS without having to define all parameters from scratch.

You can see from the example in this repo that EPSG:4326 looks increasingly distorted the farther you go from the equator (more examples below).

If the Earth was a cylinder you could unfold it like you would a rolled up map, without distortion, but how do you peel an orange in a way that forms a rectangle?

Mercator is a popular cylindrical projection, where latitudes and longitudes are projected onto a cylinder centered on/containing the WGS84 geoid, the cylinder is then unfolded to form a rectangular map. The projection effectively stretches the pointy ends of the peel, closing any gaps in the map. There will be no distortion at the equator, where the geoid and the cylinder intersect, but the farther you go from the equator, the greater the distortions (the projection error). Take Greenland for example, which is not larger than Africa.

Whether it's a problem or not depends not only on which region your application intends to cover, but also how you are planning on using your map. The distortion caused by Mercator is obviously a problem for both Europe and North America if you want to compare the size/area of US states or European countries to countries on the equator, for example.

Just like different mapping solutions have different requirements, different projections exist to fulfill them; some preserve angles or area, some are global but with heavily distorted regions (or limited to exclude these regions), and some are local, showing regions which are not properly supported by global projections and/or minimizing the projection error in specific regions of interest. While there are projections like the Dymaxion which attempt to find a middle ground between accurate local projections and global projections, they are not very intuitive (or what do you think?). Simply put, no one size fits all.

When I have referred to the weaknesses of global WGS84-type projections in other posts, I'm somewhat sloppily referring to the broad family of spherical projections like WGS84 (and by extension cylindrical projections like Mercator).

One way of avoiding this issue is of course to not project to 2D to begin with, but to present the WGS84 geoid (or a spherical approximation) in 3D as-is, which Mapbox and Google Maps have done. In theory this approach should be able to solve most of these problems, but both Mapbox and Google still struggle with singularities about four times the size of Spains total land area around the two poles collectively (very sad for me 🥲).

The case for marrying a CRS

I often use the polar stereographic projections, a family of cartesian planar projections that preserve angles. For a map of the central Arctic, the plane could be a tangent to the Earth's surface at the North Pole.

You can imagine a source of light shining from the South Pole toward the North Pole, projecting the surface of the Earth onto the projection plane (top left shows a rotated view of the projection onto the plane to illustrate the concept):

If we are planning on being further south, say at 80°N, we can move the plane so it intersects the Earth at this latitude. The intersection between the Earth and the plane is called the latitude of true scale, and as the name implies, there is no distortion at this intersection. The farther we move away from the latitude of true scale, either north or south, the greater the distortions in the map. Moving north is slightly more forgiving than moving south, because it's flatter.

Scales are typically well preserved throughout the central Arctic, but the map gets increasingly distorted as it approaches the equator. Because a given (x, y) coordinate in the plane can only correspond to one point on the geoid, the projection is limited to either the northern or the southern hemisphere. In practice these maps exist to cover the gaps in popular projections such as those mentioned earlier, and as such, they rarely extend below 60°N or above 60°S, where they will not perform well anyway.

One example of a polar stereographic CRS is EPSG:3413, for which there is an example in this repo.

This CRS is also using the WGS84 datum, so it's painless to transition between local coordinates and lat/lon. This is something users of cartesian coordinate systems are reluctant to do however.

For one, the CRS is already cartesian/planar, so it doesn't need to be projected, all you need to do is scale your points to match the zoom level in flutter_map before you plot them.

As an aside: if your input projection is the same as the output projection that you render on the map, you technically don't even need to declare a CRS unless you want to translate backwards from local points to LatLng. Being able to bypass the CRS is especially useful in situations like mapping floor plans which might not correspond to an existing CRS: https://github.com/fleaflet/flutter_map/pull/1750#issuecomment-1884042084

Let's say we want to measure the distance between two LatLng. Rather than using the very complex/accurate Vincency method, we use the simple Haversine method (from latlong2):

double distance(final LatLng p1, final LatLng p2) {
  final sinDLat = math.sin((p2.latitudeInRad - p1.latitudeInRad) / 2);
  final sinDLng = math.sin((p2.longitudeInRad - p1.longitudeInRad) / 2);

  // Sides
  final a = sinDLat * sinDLat +
      sinDLng *
          sinDLng *
          math.cos(p1.latitudeInRad) *
          math.cos(p2.latitudeInRad);
  final c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a));

  return equatorRadius * c;
}

In a polar CRS we would just use the Pythagorean theorem:

var dist = const Point(1, 2).distanceTo(const Point(3, 4));

The projection has little distortion in the central Arctic, so this is fairly accurate even over longer distances.

Calculating the bearing between two points using latlong2 involves calls to: sin, cos, tan, sin, cos, atan2. In a cartesian space? atan2. By constraining the CRS you get the best of both worlds, low computational cost/high speed as well as high accuracy, but only in certain regions which are relevant to you.

The fact that virtually anything you do with the points come at a fraction of the computational cost and programmatic complexity makes it more feasible to do heavy lifting on the fly, for example, observing and modeling/prognosticating behavior of vehicles, vessels or sensors as their locations are continuously fed to your app.

In other words, your app will be riddled with Points. You will end up relying mostly on dart:math to work with the coordinates, simply because it's intuitive, blazingly fast and native to your CRS. Like @ignatz put it, I seem married to a particular CRS, and I am, but not without reason.

The difference projection makes

Below you can see the same image, of the Arctic ocean, in three different projections.

Polar stereographic projection produces a circular map, but since its four circular arcs have been cropped to make the image a rectangle in EPSG:3996, it produces waves at the bottom edges of the EPSG:4326 and EPSG:3857 images which are expected, and not a projection error:

EPSG:4326 (WGS84):

EPSG:3857 (Web Mercator*):

The Web Mercator's distortions are not as eye-catching as the Equirectangular projection (WGS84), but that doesn't mean that the distortions are not there, in fact, they are extreme. Note that Greenland is about 4-5 times the width of Norway & Sweden for Web Mercator, but in IBCAO Polar Stereographic, it is just about twice the width. This is in a region where the polar projection is accurate, and the Web Mercator projector struggles the most. The reason Web Mercator doesn't show as much of the map as the other two projections is because it doesn't support going beyond 85°N.

EPSG:3996 (IBCAO Polar Stereographic):

* Web Mercator is the default for Google Maps, Mapbox etc. when rendering in 2D.

Where flutter_map comes in

One could argue that flutter_map already supports cartesian coordinate systems, but it doesn't support a cartesian-centric approach to development (which is where you will end up if your cartesian app is anything beyond a trivial example).

There is CrsSimple (likely inspired by Leaflet):

// Custom CRS for non geographical maps
@immutable
class CrsSimple extends _CrsWithStaticTransformation {
  const CrsSimple()
      : super(
          code: 'CRS.SIMPLE',
          transformation: const _Transformation(1, 0, -1, 0),
          projection: const _LonLat(), // <---
          infinite: false,
          wrapLat: null,
          wrapLng: null,
        );
}

https://leafletjs.com/examples/crs-simple/crs-simple.html

This is assuming that you will treat the LatLon somewhat awkwardly as if they were cartesian coordinates.

To me it feels like shoe-horning in support for simple cartesian CRS:

• We have bounds on lat/lon which means you might not be working in your native, but scaled, coordinates.
• Semantically it's completely incorrect, which complicates the library users code and can lead to misunderstandings
  • Point is intentionally ambiguous, but LatLng is very specifically in our context referring to spherical coordinates
  • By treating very different types as if they were the same you effectively don't have any type safety (throw two SimpleCrs LatLngs at Vincenty above and it would gladly give you a completely incorrect distance measurement)

Another way of working around the limitations is to take your cartesian coordinates, unproject them to LatLng, only to let the layers project them back into cartesian coordinates again.

The cost of passing a list of, say, 10k points to a layer (just getting the coordinates into the layer), is suddenly the cost of (un)projecting 10k points and then reprojecting them back to the very same points again. There is also a slight loss in precision in jumping back and forth between CRS. It's a very expensive per-frame cost if you compare it to what passing a list to a method costs.

Possible solutions

When we project LatLng inside our components, we effectively turn these coordinates into cartesian coordinates. In other words, most of flutter_map is cartesian in nature. This is what makes working in cartesian coordinate systems a bit frustrating, because the support for these coordinates is there already, it's just that it's hidden behind the API:s.

In a (normal) map app there will always be two geographical coordinate systems, spherical LatLng and projected/planar Point, it would make sense if the user could choose, at least to some degree, in which system they want to work.

For me the bulk of the work in forking portions of flutter_map has been in the layer code. There were parts of the code where my projection unavoidably necessitated a re-write, but there were also pretty large parts where the code could easily have supported my use-case if only the API was a bit more clever. I also believe this is where users of local CRS tend to get stuck, wondering how to get their non-LatLng into their layers, especially when there is no clear mapping between Point and LatLng.

It's a tricky problem and I don't have a good solution off-hand, but I haven't mulled on it for too much either yet. Just having a way of getting Points into the layers would get us a long way. This could also benefit people that are plotting large quantities of points, because there is no need to project the points from LatLng to Point for every frame. If you have a static list of 1 million points and you project them when your app starts, that projected list can last the lifetime of your app.

[ignatz]

Thanks for the write-up. I certainly hadn't expected the MapBox issues you encountered but just briefly looking at it, I guess that's the kind of continuity issue you can expected when you jump between projections.

Anyway back to your FlutterMap pains. I feel like you're mixing a few points, which is no problem, but let me try to pull them apart.

The "performance" arc

I'll ignore the distance measurement argument since errors and precision requirements are highly domain specific. I don't think any generic library should fall back to Pythagorean theorem. I take this more as point towards the second arc, you're operating in planar space already thus you can save some cycles for some points in the right places.

As for re-projecting on every frame, I agree with you that it adds to the cost. That said, I don't think that means that users should pre-project all points. Users who use LatLng would also benefit from not re-projecting. Thus breaking up the crs into a projetion and transformation stage would make a lot of sense. E.g. have the polygon layer project all points on construction and then just transform them when needed :thumbs up:

The "I don't wanna deal with projections" arc

I can certainly see the two somewhat different use-cases:

• You're doing mapping that doesn't cover the globe, e.g. in an office building
• Or your use case: while LatLng would work, I'd rather work in planar space either because that's how points are supplied to me or I just wanna do simpler math and am ok with the errors.

FWIW:

One could argue that flutter_map already supports cartesian coordinate systems, but it doesn't support a cartesian-centric approach to development (which is where you will end up if your cartesian app is anything beyond a trivial example).

Personally, I wouldn't argue that. I'd consider the SimpleCrs a hack at best.

Another way of working around the limitations is to take your cartesian coordinates, unproject them to LatLng, only to let the layers project them back into cartesian coordinates again.

Not for the office-building-case but for yours I could see this to be a reasonable or even good approach especially if the re-projection concerns would be addressed. I assume it would at least reduce the size of your fork.

As an aside: if your input projection is the same as the output projection that you render on the map, you technically don't even need to declare a CRS unless you want to translate backwards from local points to LatLng.

In the LatLng-centric state FM is in, this would be a non-feature but I see that this would make sense in world where flutter map supports 1st-class planar CRSs.

... My take-aways:

• We should definitely get rid of unnecessary re-projections, because why not?.
• It remains an open discussion in FlutterMap should/could support cartesian coordinates.
• Many optimizations/trade-offs are domain-specific. Pre-projected points could be allowed by some or all layers if we break up the Crs and add a secondary constructor

[JosefWN]

Thanks for the write-up. I certainly hadn't expected the MapBox issues you encountered but just briefly looking at it, I guess that's the kind of continuity issue you can expected when you jump between projections.

The only jumping I am doing is from LatLng to Globe in Mapbox (one projection), it's a simpler frontend used for a slightly different purpose but in the same region. I don't know the specifics about that particular projection but they have struggled a lot with Globe (and with code quality in general) the past 1-2 years. A whole lot of bugs in their releases.

I'll ignore the distance measurement argument since errors and precision requirements are highly domain specific. I don't think any generic library should fall back to Pythagorean theorem. I take this more as point towards the second arc, you're operating in planar space already thus you can save some cycles for some points in the right places.

Of course, that was my point. It would be a terrible idea to use it broadly, but for me it makes a lot of sense. The point of that section is intended to motivate why I am having Point everywhere in my code (which from an outsiders perspective might seem a bit crazy).

I made this ticket to try to explain, to the best of my abilities, how I use flutter_map and some of the challenges it has posed.

As I stated in the very beginning of the post ^

In the LatLng-centric state FM is in, this would be a non-feature but I see that this would make sense in world where flutter map supports 1st-class planar CRSs.

Yes, some level of stringency is required. From my point of view I always use a CRS, but I can certainly feel for the people that don't. Part of the purpose of the post is to highlight that the GIS-world is pretty diverse and that there is more than one way to skin a cat, but perhaps it was superfluous.

[JaffaKetchup]

Thanks for this write up @JosefWN (and discussing @ignatz), very helpful for CRS-beginners like myself and I'm sure plenty of other people :D. We are paying attention to this thread, we've just got other higher priority things to focus on at the moment.

[JaffaKetchup]

I'm converting this to a discussion, since I think it is better placed there and can be long-lasting.