crs_fin <- 3067
crs_ll <- 4326R, leaflet, and MML tiles: Changing the output CRS
maps
Task
Make an interactive map (HTML widget) where the Coordinate Reference System (CRS) of the map is the official CRS used by Finnish institutes, TM35FIN / EPSG:3067.
Define the target CRS and Mercator EPSG codes to be used later.
We will use these libraries.
library(sf)
library(leaflet)
library(tmap)
library(geofi)Preliminaries
Example data: Polygons
To have some data to depict and check the output, we fetch the municipalities of Finland using the geofi package.
mun <- geofi::get_municipalities()Requesting response from: http://geo.stat.fi/geoserver/wfs?service=WFS&version=1.0.0&request=getFeature&typename=tilastointialueet%3Akunta4500k_2024Warning: Coercing CRS to epsg:3067 (ETRS89 / TM35FIN)Data is licensed under: Attribution 4.0 International (CC BY 4.0)Example basemap that supports the target CRS
Most likely the data is to plotted on top of some kind of a basemap, also known as a theme map or layer. A typical basemap service is of the type Web Map Service (WMS) or Web Map Tile Service (WMTS), and provides pre-rendered pixel-images known as tiles. These tiles are fast to download and are to be stiched together as the background map by the rendering engine. The images usually come in the unprojected lon-lat, Mercator CRS.
An example of a WMTS service provider that supports the Finnish goverment CRS is the National Land Survey of Finland, or Maanmittauslaitos (MML). Details of their open data services can be found at https://www.maanmittauslaitos.fi/karttakuvapalvelu. The English side is not complete, so try to parse the details from the Finnish pages.
To use this WMTS provider, a free API key is needed. I’ve stored my key in the environmental variable MML_API_KEY, defined in the local ~/.Rprofile . See Sys.setenv help page for doing this.
To use the custom WMTS, we define a url-template string that describes the service in a way that the leaflet-engine can request the tiles.
mml_key <- Sys.getenv("MML_API_KEY") # use your own here.
# url template; From an example
umml_fin <- "https://avoin-karttakuva.maanmittauslaitos.fi/avoin/wmts/1.0.0/taustakartta/default/ETRS-TM35FIN/{z}/{y}/{x}.png"
# Augment the template with the api key
umml <- sprintf("%s?api-key=%s",umml_fin, mml_key)Note the “ETRS-TM35FIN” string. This is where we request the tiles to be in the target projection. The other interesting string is “taustakartta”, which defines the layer we want. “selkokartta” and “maastokartta” would be alternatives, see 1.
Note: I think this is not clean WMTS, more like WMS, but I’ve not found a way to make leaflet-package support the full WMTS 🤷
Making the CRS understandable by leaflet
leaflet assumes that all input data are in EPSG:3867 Web Mercator CRS 2. The list of basemap provides all default to this CRS, and will not work with the target projection.
Other projections can be used, but the CRS need to be specified manually. Let us define the target CRS in a format that leaflet (and derivatives) understands. The specification is entagled with the properties of the tile service, particularly the resolution and the origin of the tile mosaics, so adjust accordingly. We need to
- Use the Proj4-string syntax. This string we get from https://epsg.io/3067 (scroll down therein).
- The resolution (meters-to-pixels) and the origin (top-left corner) of the tile service
The CRS information, to be used with the MML tiles, is then compiled like so:
lcrs <- leafletCRS(crsClass = "L.Proj.CRS",
code = "EPSG:3067" , # just a name tag
proj4def = "+proj=utm +zone=35 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs",
origin = c(-548576, 8388608), # This is MML basemap related
resolutions = 2^(13:-1) # This is MML basemap related
)The origin and resolution details were found in a technical document provide by MML 3.
Interactive map using leaflet
The projection information is now ready to be applied to a map. First though: leaflet expects all input data to be in lon-lat, so we need to first transform the data to be shown. The back-transformation to the target CRS will take place inside leaflet code when rendered.
mun_ll <- sf::st_transform(mun, crs_ll)We can now produce the map. The target CRS for output is provided in the options parameter.
leaflet(mun_ll,
options = leafletOptions(crs = lcrs) ## Here is the main thing, will project data
) |>
# add the basemap, now already in the target CRS
addTiles( urlTemplate = umml ) |>
# then draw the example data
addPolygons(fillColor = ~colorFactor("viridis", domain = 1:16)(ely_code),
popup = ~name,
weight = 2,
highlightOptions = highlightOptions(color = "red", bringToFront = TRUE)
) |>
setView(lng = 25, lat = 68, zoom = 3) # lon-lat again
The projection is now correct.
tmap example, uses leaflet
If the above example for leaflet worked, we should be able to also use tmap, which provides a higher level functionality for producing maps.
# Set interactive mode
tmap_mode("view")
tm_shape(mun_ll) +
tm_view(projection = lcrs, # here is the main thing, will project data
set.view = c(lon = 25, lat = 68, zoom = 3) # lon-lat again
) +
tm_basemap(server = umml) +
tm_borders() +
tm_fill(col = "ely_name_fi", alpha = .8, legend.show = FALSE,
id = "name",
popup.vars = c("name", "ely_name_fi")
) 
Note that while writing this example the tm_shape data needed to be in lat-lon too if using a custom CRS basemap provider. I think this is because when manipulating the CRS using tm_view, the data is automatically re-projected from assumed lon-lat regardless of the data CRS.