0

I have the need to export georeferenced images from Leaflet.js on the client side. Exporting an image from Leaflet is not a problem as there are plenty of existing plugins for this, but I'd like to include a world file with the export so the resulting image can be read into GIS software. I have a working script fort his, but I can't seem to nail down the correct parameters for my world file such that the resulting georeferenced image is located exactly correctly.

Here's my current script

// map is a Leaflet map object
let bounds = map.getBounds();   // Leaflet LatLngBounds
let topLeft = bounds.getNorthWest();
let bottomRight = bounds.getSouthEast();
let width_deg = bottomRight.lng - topLeft.lng;
let height_deg = topLeft.lat - bottomRight.lat;
let width_px = $(map._container).width()    // Width of the map in px
let height_px = $(map._container).height()  // Height of the map in px
let scaleX = width_deg / width_px;
let scaleY = height_deg / height_px;   

let jgwText = `${scaleX}
0
0
-${scaleY}
${topLeft.lng}
${topLeft.lat}`

This seems to work well at large scales (ie zoomed in to city-level or so), but at smaller scales there is some distortion along the y-axis. One thing I noticed is that all examples of world files I can find (and those produced from QGIS or ArcMap) all have the x-scale and y-scale parameters being exactly equal (oppositely signed). In my calculations, these terms are different unless you are sitting right on the equator.

Example world file produced from QGIS

0.08984380916303301   // x-scale (size of px in x direction)
0                     // rotation parameter 1
0                     // rotation parameter 2
-0.08984380916303301  // y-scale (size of px in y direction)
-130.8723208723141056 // x-coord of top left px
51.73651369984968085  // y-coord of top left px

Example world file produced from my calcs

0.021972656250000017
0
0
-0.015362443783773333
-130.91308593750003
51.781435604431195

Example of produced image using my calcs with correct state boundaries overlaid:

Does anyone have any idea what I'm doing wrong here?

IvanSanchez
  • 18,272
  • 3
  • 30
  • 45
  • Are we talking EPSG:3857 or EPSG:4326? – IvanSanchez Sep 17 '20 at 17:44
  • @IvanSanchez, I've been working in 4326 just because it's simpler to get geographic lat lngs from Leaflet. I can certainly use a projected coordinate system though if it will help. The example world files shown are in 4326 – Nic Williams Sep 17 '20 at 21:59
  • What I mean is, are your *images* in EPSG:4326? – IvanSanchez Sep 17 '20 at 23:29
  • I'm not sure I understand completely. I'm producing normal, non-georeferenced images by exporting from Leaflet using Leaflet.Export (https://github.com/Flexberry/Leaflet.Export). I'm then attempting to manually create a world file for this image by grabbing parameters from the Leaflet window on export (pixel size in X and Y direction, coordinate location of the top left px). The coordinates used in this world file will dictate the CRS of the final image right? I've been trying in WGS84 to avoid coordinate transformations so far. Sorry if there's a fundamental munderstanding on my end – Nic Williams Sep 18 '20 at 00:01
  • 1
    @IvanSanchez, problem solved! I realized in my attempts to use EPSG:3857 I'd calculated the width and height of the map using Leaflet's L.map.distance() function, which returned correct distances using EPSG:4326. I switched to calculating the width and height by projecting the corner points to 3857 and simply subtracting X and Y values. This gave me equal x and y scales and a perfect fit. Thanks for the help! – Nic Williams Sep 18 '20 at 00:29

1 Answers1

1

Problem was solved by using EPSG:3857 for the worldfile, and ensuring the width and height of the map bounds was also measured in this coordinate system. I had tried using EPSG:3857 for the worldfile, but measured the width and height of the map bounds using Leaflet's L.map.distance() function. To solve the problem, I instead projected corner points of the map bounds to EPSG:3857 using L.CRS.EPSG3857.project(), the simply subtracted the X,Y values.

Corrected code is shown below, where map is a Leaflet map object (L.map)

// Get map bounds and corner points in 4326
let bounds = map.getBounds();
let topLeft = bounds.getNorthWest();
let bottomRight = bounds.getSouthEast();
let topRight = bounds.getNorthEast();

// get width and height in px of the map container
let width_px = $(map._container).width()
let height_px = $(map._container).height()

// project corner points to 3857
let topLeft_3857 = L.CRS.EPSG3857.project(topLeft)
let topRight_3857 = L.CRS.EPSG3857.project(topRight)
let bottomRight_3857 = L.CRS.EPSG3857.project(bottomRight)

// calculate width and height in meters using epsg:3857
let width_m = topRight_3857.x - topLeft_3857.x
let height_m = topRight_3857.y - bottomRight_3857.y

// calculate the scale in x and y directions in meters (this is the width and height of a single pixel in the output image)
let scaleX_m = width_m / width_px
let scaleY_m = height_m / height_px

// worldfiles need the CENTRE of the top left px, what we currently have is the TOPLEFT point of the px.
// Adjust by subtracting half a pixel width and height from the x,y
let topLeftCenterPxX = topLeft_3857.x - (scaleX / 2)
let topLeftCenterPxY = topLeft_3857.y - (scaleY / 2)

// format the text of the worldfile
let jgwText = `
${scaleX_m}
0
0
-${scaleY_m}
${topLeftCenterPxX}
${topLeftCenterPxY}
`

For anyone else with this problem, you'll know things are correct when your scale-x and scale-y values are exactly equal (but oppositely signed)!

Thanks @IvanSanchez for pointing me in the right direction :)