4
votes

I need to convert Eastings and Northings OSGB36 coordinates into latitude and longitude coordinates using Go. Therefore I'm wondering if there is any package that do just this. Writing one from scratch doesn't seem trivial. I'm running the code in a sandboxed VM so the code must be pure Go.

Input:

  • Northing - Distance in metres north of National Grid origin.
  • Easting - Distance in metres east of National Grid origin.

Example:

 348356,862582

Output (Decimal degrees -DDD):

  • Latitude
  • Longitude

Example:

41.40338, 2.17403
3
What kind of co-ordinates do you have and what do you want - you need to be a bit more specific as there are dozens of different types of co-ordinates in use. An example would be good.Nick Craig-Wood
I've updated the question. Hopefully it makes sense now.themihai
@JeremyBanks I'm trying to use "go". As the question seamed unclear I provided the equivalent in PHP. However it seems that leads to confusion so now I just removed the PHP example.themihai
You might try having a look at this project.Dmitri Goldring
@dommage it doesn't seem to work. I have only Northing, Easting coordinatesthemihai

3 Answers

3
votes

You could try the go-proj-4 library which is a wrapper to the comprehensive PROJ.4 - Cartographic Projections Library.

I haven't tried either but it looks like it should do the job.

Alternatively you could port the code on this page

2
votes

The code below is very naively ported from the javascript (© 2005-2014 Chris Veness) at the link Nick Craig-Wood provided. I've only ported the OsGridToLatLong() function, but the others shouldn't be too difficult. Also, this treats all values as float64. You may want to treat northing and easting as int.

package main

import (
    "fmt"
    "math"
)

const (
    radToDeg = 180 / math.Pi
    degToRad = math.Pi / 180

    a    = 6377563.396
    b    = 6356256.909  // Airy 1830 major & minor semi-axes
    f0   = 0.9996012717 // NatGrid scale factor on central meridian
    lat0 = 49 * degToRad
    lon0 = -2 * degToRad // NatGrid true origin
    n0   = -100000.0
    e0   = 400000.0        // northing & easting of true origin, metres
    e2   = 1 - (b*b)/(a*a) // eccentricity squared
    n    = (a - b) / (a + b)
    n2   = n * n
    n3   = n * n * n
)

func OsGridToLatLong(northing, easting float64) (float64, float64) {
    lat := lat0
    m := 0.0
    for northing-n0-m >= 1e-5 { // until < 0.01mm
        lat = (northing-n0-m)/(a*f0) + lat
        ma := (1 + n + (5/4)*n2 + (5/4)*n3) * (lat - lat0)
        mb := (3*n + 3*n*n + (21/8)*n3) * math.Sin(lat-lat0) * math.Cos(lat+lat0)
        mc := ((15/8)*n2 + (15/8)*n3) * math.Sin(2*(lat-lat0)) * math.Cos(2*(lat+lat0))
        md := (35 / 24) * n3 * math.Sin(3*(lat-lat0)) * math.Cos(3*(lat+lat0))
        m = b * f0 * (ma - mb + mc - md) // meridional arc
    }
    cosLat := math.Cos(lat)
    sinLat := math.Sin(lat)
    nu := a * f0 / math.Sqrt(1-e2*sinLat*sinLat)                 // transverse radius of curvature
    rho := a * f0 * (1 - e2) / math.Pow(1-e2*sinLat*sinLat, 1.5) // meridional radius of curvature
    eta2 := nu/rho - 1
    tanLat := math.Tan(lat)
    tan2lat := tanLat * tanLat
    tan4lat := tan2lat * tan2lat
    tan6lat := tan4lat * tan2lat
    secLat := 1 / cosLat
    nu3 := nu * nu * nu
    nu5 := nu3 * nu * nu
    nu7 := nu5 * nu * nu
    vii := tanLat / (2 * rho * nu)
    viii := tanLat / (24 * rho * nu3) * (5 + 3*tan2lat + eta2 - 9*tan2lat*eta2)
    ix := tanLat / (720 * rho * nu5) * (61 + 90*tan2lat + 45*tan4lat)
    x := secLat / nu
    xi := secLat / (6 * nu3) * (nu/rho + 2*tan2lat)
    xii := secLat / (120 * nu5) * (5 + 28*tan2lat + 24*tan4lat)
    xiia := secLat / (5040 * nu7) * (61 + 662*tan2lat + 1320*tan4lat + 720*tan6lat)
    de := easting - e0
    de2 := de * de
    de3 := de2 * de
    de4 := de2 * de2
    de5 := de3 * de2
    de6 := de4 * de2
    de7 := de5 * de2
    lat = lat - vii*de2 + viii*de4 - ix*de6
    lon := lon0 + x*de - xi*de3 + xii*de5 - xiia*de7
    return lat * radToDeg, lon * radToDeg
}

func main() {
    lat, lon := OsGridToLatLong(348356.0, 862582.0)
    fmt.Printf("Latitude: %fN\nLongitude: %fE\n", lat, lon)
}

Produces:

Latitude: 52.833026N
Longitude: 4.871525E

These are OSGB-36 coordinates.

From Chris Veness's original post:

"The Ordnance Survey uses ‘OSGB-36’, based on an elliptical model of the earth’s surface which is a good fit to the UK. GPS systems generally use the world-wide ‘WGS-84’, based on an elliptical model which is a best approximation to the entire earth. At Greenwich, these differ by about 126m (they coincide somewhere in the Atlantic ocean; there’s more on Wikipedia)"

I have no idea if this is where the OSGB-36 coords for Northing: 348356, Easting: 862582 are meant to be, but the code above, and the javascript code (JSFiddle), put it somewhere in the northern Netherlands. (Although the coords shown on the map are WGS-84, and not OSGB-36. See converting between OSGB-36 & WGS-84 for more details.)

Probably completely wrong

This code hasn't been properly tested (It produces the same output as the original javascript code, but that's no guarantee of correctness), and may have hideous bugs, but it should point you in the right direction (no pun intended).

Playground

1
votes

Hi I had the same problem. I was porting our codebase from python. We used this library - (UTM python version) by Tobias Bieniek. It is public on github.

I didn't find something good and lightweight. I just ported this library on golang - UTM golang library. It's pure go which contain only two methods. You can find badge-link to godoc in repository. I use it in a production.

For speed use:

go get github.com/im7mortal/UTM