tags:

views:

333

answers:

3

Hi all, hope you are well

I am trying to convert lat/long coordinates to OSGB36 x and y using the proj.4 library.

Has anyone else successfully done this? I need to fill the srcPrj4String and destPrj4String variables, e.g.

string srcPrj4String = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs";
string destPrj4String = "+proj=utm +zone=11 +ellps=GRS80 +datum=NAD83 +units=m";

but I can't figure out the what the destPrj4String should be with OSGB36 - i know the datum should be +datum=OSGB36, but everything I try, doesn't work

Any ideas?

Many thanks in advance

leddy

+2  A: 

got it:

string srcPrj4String = "+proj=longlat +ellps=WGS84 +towgs84=0,0,0 +no_ defs"; (remove the space between no_ and defs - for some reason it's making the text italic)

string destPrj4String = "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894 +units=m +no_defs";

cheers!

leddy
+1  A: 

Googling turns up this from Dr John Stevenson, a Manchester University Earth Science academic - who should get it right if anyone does. Here's a quote.


The problem was that going to OSGB36 requires both a projection and a datum conversion. Prior to October 2007, proj was only carrying out the projection, thus resulting in the large offset. You can check if you have the new version by running 'proj -v' or by looking at your epsg file:

cat /usr/share/proj/epsg | grep -A 1 "British National Grid" 

# OSGB 1936 / British National Grid 
<27700> +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 
+y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs  <>

The new versions have +datum=OSGB36.

If you have an old version, you can correct it by replacing the line with:

+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.999601 +x_0=400000 +y_0=-100000 
+ellps=airy 
+towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894 +units=m 
+no_defs <>

A complication is that OSGB36 is slightly distorted with respect to GPS projections (such as WGS84 and ETRS89). This offset is small, and is only important for higher precision surveying. Many searches about OSGB36 offsets bring up pages relating to this. If you want to compensate for this too, you can download a nadgrid file and use it. For my data, this moved the points by about 1 m.

MarkJ
A: 

EPSG:27700 on spatailreference.org gives various strings for defining this, including one for proj4.

Here is example code in ruby using the proj4 bindings :

#!/usr/bin/ruby
require 'rubygems'
require 'proj4'

#Some example WGS84 lat lon coordinates to convert:
lon = -0.10322
lat = 51.52237

srcPoint = Proj4::Point.new(Math::PI * lon.to_f / 180, 
                            Math::PI * lat.to_f / 180)

srcPrj  = Proj4::Projection.new("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs") 
destPrj = Proj4::Projection.new("+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.999601 +x_0=400000 +y_0=-100000 +ellps=airy +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894 +units=m +no_defs <>")

point = srcPrj.transform(destPrj, srcPoint)

puts "http://www.openstreetmap.org/?mlat=" + lat.to_s + "&mlon=" + lon.to_s + "&zoom=16"
puts "Converts to:";
puts "http://streetmap.co.uk/grid/" + point.x.round.to_s + "_" + point.y.round.to_s + "_106"

The output:

http://www.openstreetmap.org/?mlat=51.52237&amp;mlon=-0.10322&amp;zoom=16
Converts to:
http://streetmap.co.uk/grid/531691_182089_106

So this is working accurately now. Originally I was trying just the 'destPrj' string, and calling the 'forward' method, but this refused to do the datum conversion, resulting everything being 100m out. It seemed to be necessary to use the 'srcPrj' string and the 'transform' method, to get the datum conversion happening.

See also my blog post: Ruby code for converting to UK Ordnance Survey coordinate systems from WGS84? which includes a pure ruby version (not proj4) for doing the same

Harry Wood