Edit after input below and e.g. here: leaflet has to make a compromise between "understandable for non-gis-experts" and "correct gis-handling". What I understand now, leaflet focussed on easy handling, thus, projections cannot really be implemented.
Originial question: My goal is to display "real" densities correctly, that is in events/m^2. At the moment I am struggeling with leaflet's projection's. What do I have to do, to display all crsClasses correctly? My data generates uniform random points on the shere.
I guess I miss to transform the data?
library(leaflet)
n <- 10000
z <- 2*runif(n) - 1
phi <- 2*pi*runif(n) - pi
x <- sin(phi)*sqrt(1 - z^2)
y <- cos(phi)*sqrt(1 - z^2)
theta <- acos(z)
# lat / lng
lat <- theta*180/pi - 90
lng <- phi*180/pi
# Working fine:
# plot3D::polygon3D(x, y, z)
rgl::plot3d(x, y, z)
# rgl::plot3d(sin(phi)*sin(theta), cos(phi)*sin(theta), cos(theta)) # the same...
# EPSG:3857, also known as "Google Mercator" or "Web Mercator", the first in the following list
crsClasses <- c("L.CRS.EPSG3857", "L.CRS.EPSG4326", "L.CRS.EPSG3395", "L.CRS.Simple", "L.Proj.CRS")
epsg3857 <- leafletCRS(crsClass = crsClasses[1])
leaflet(options = leafletOptions(
crs = epsg3857, worldCopyJump = FALSE)) %>%
addTiles() %>%
addProviderTiles(providers$OpenStreetMap, group = "OSM") %>%
addProviderTiles(providers$Stamen.TonerLite, group = "Toner Lite") %>%
addLayersControl(baseGroups = c("OSM", "Toner Lite")) %>%
addCircleMarkers(lat = lat,
lng = lng,
radius = 1,
stroke = FALSE,
fillOpacity = 0.7)
Using crsClasses[1]
... looks ok, but I don't understand why this corresponds to worldCopyJump = FALSE
.
Using crsClasses[2]
... looks strange
crsClasses[3...6]
also look strange. It is really hard to find good tutorials for a GIS-newbie
References:
Comment after input from IvanSanchez: Is the following strategy correct?
- Using crsClass L.CRS.EPSG4326, the map changed, but the data remained unchanged. Thus, I also have to transform the data.
This means I have to use
new <- sp::spTransform(x, CRSobj, ...)
. It looks like- I have to figure out the EPSG code, here you already told me (and you could guess) it is EPSG:4326. Is there a systemtic way to get that result?
- I need to transform (lat, lng) to an a proper object. From
??spTransform
: "x ... object to be transformed" and "CRSobj ... object of class CRS". The help leads tosp::CRS(projargs, doCheckCRSArgs=TRUE)
. But then I get stuck... - Call the new lat/lng in addCircleMarkers().
- As you pointed out, I should look for equal-area projection for densities
- From here I would look for a Lambert azimuthal equal-area (LAEA), as they recommend. From here, I would guess, I take the EPSG:3575(??) for Europe, but perhaps the "WGS 84 / North Pole LAEA Europe" is wrong?
- The I am lost: ...find (or build) base maps for that projection. I think this means for the map and the data?
Can you provide a guidline through that jungle or at least to one of the cases above? Or do you have a good reference to learn how to do that properly?
From the Leaflet docs:
crs ... Coordinate Reference System to use. Don't change this if you're not sure what it means.
But unfortunatelly, they don't tell you, where to look if you need to change.