How to recreate coordinate transformation from vesselfinder (EPSG:4326 to EPSG:3857)?

Question

I'm currently trying to figure out how vesselfinder.com calculates its Box Boundaries (bbox) which they use to query data from their backend.

Given an input like: lat, lon = 59.8230, 22.9586

They fetch data by using this bbox: 13761899,35886447,13779795,35898097

If I try to get a similar bbox by using bboxfinder.com, I get the following values, which aren't even close to what I was expecting: 2553560.4710,8358928.9331,2556565.4293,8360514.8411

The website above is using EPSG:4326 (WGS 84) to EPSG:3857 (WHS 84 / Pseudo-Mercator) by default. I tried to verify in the JS code of vesselfinder that they're using this conversion as well.

    var c = new s.geom.MultiLineString(t);
    return c.transform('EPSG:4326', 'EPSG:3857'),

There are also the following ones mentioned, but I'm pretty sure, that it has to be the upper shown transformation.

it = [
  new $('EPSG:3857'),
  new $('EPSG:102100'),
  new $('EPSG:102113'),
  new $('EPSG:900913'),

The questions now are: What am I doing wrong? / Where do I think wrong?

I also tried using Python for the conversion and even tried out the other mentioned EPSG:XXXXXX types, but haven't got the desired result. I also changed the order of both EPSG types when creating the Transformer, but again, not the desired results.

from pyproj import Transformer

TRAN_4326_TO_3857 = Transformer.from_crs("EPSG:4326", "EPSG:3857", always_xy=True)

lat = 59.823002
lon = 22.958583
expansion = 2000

res = TRAN_4326_TO_3857.transform(lng, lat)
bbox = (round(res[0]-expansion), round(res[1]-expansion), round(res[0]+expansion), round(res[1]+expansion))

print(bbox)
# (2455738, 8260436, 2655738, 8460436)

This one is close to the one I got from bboxfinder, but is again not even close to the bbox vesselfinder is using.

For what it's worth, I get the same results as you for EPSG:3857, so the answer to what you're doing wrong is probably "Assuming that the CRS is EPSG:3857". This raises the obvious follow-up question "So what is it, then?" to which, unfortunately, I have no answer. — Ture Pålsson, Aug 13 '22 at 18:26
Exactly, it's just an assumption made from the findings in their JS code. But glad to hear, that's most likely not an implementation issue. Hope someone else has an idea. Otherwise, I'll try to bruteforce all CRS and figure out which one comes close to the required bbox. — , Aug 13 '22 at 19:56
I don't know what vesselfinder is doing, but your Python code mixed up lon and lat. `always_xy` in Transformer means longitude first, latitude second, so it should be `transform(lon, lat)`. This gives the result close to bboxfinder.com, (2553738, 8358436, 2557738, 8362436). — Michael Entin, Aug 13 '22 at 21:29
This is what I also got when playing around with the order of lat / lon. I'm now d'accord with bboxfinder, but again not with vesselfinder. Anyway, thanks for the hint. — , Aug 14 '22 at 08:10
Try the: ``EPSG:4326 EPSG:6933`` also you can refer https://stackoverflow.com/questions/72393359/how-to-add-point-to-the-linestring — R. Baraiya, Aug 24 '22 at 07:51

score 1 · Answer 1 · answered Aug 21 '22 at 19:00

https://gis.stackexchange.com/a/370496 seems to have the math.

convertCoordinates(lon, lat) {
    var x = (lon * 20037508.34) / 180;
    var y = Math.log(Math.tan(((90 + lat) * Math.PI) / 360)) /
                  (Math.PI / 180);
    y = (y * 20037508.34) / 180;
    return [x, y];
}

Or, in C# ( https://gis.stackexchange.com/a/325551 )

public static double[] latLonToMeters(double lat, double lon)
{
    //Debug.Log(Mathd.Tan((90d + lat) * Mathd.PI / 720));
    //Debug.Log(Mathd.Tan((90d + lat) * Mathd.PI / 360d));
    double[] meters = new double[2];
    meters[0] = lat * originShift / 180d;
    meters[1] = Mathd.Log(Mathd.Tan((90d+lon) * Mathd.PI / 360d)) /
                   (Mathd.PI / 180d);
    //meters[1] = Mathd.Log(5) / (Mathd.PI / 180d);
    meters[1] = meters[1] * originShift / 180d;
    return meters;
}

In any case, note the website those came from; that may be a better place to get the algorithm. (Then if you need help converting to your preferred language, come back here.)

How to recreate coordinate transformation from vesselfinder (EPSG:4326 to EPSG:3857)?

1 Answers1