-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathlat_lon_to_bng.py
executable file
·98 lines (67 loc) · 2.59 KB
/
lat_lon_to_bng.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
#!/usr/bin/env python
import numpy as np
import sys
import string
import argparse
import netCDF4 as nc
import matplotlib.pyplot as plt
import datetime as dt
import mpl_toolkits.basemap.pyproj as pyproj
################################################################################
################################################################################
#
# Define the class
#
################################################################################
################################################################################
class llbng:
def __init__(self):
# projection from:
# http://spatialreference.org/ref/epsg/27700/proj4/
self.osgb36=pyproj.Proj('+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' ) # UK Ordnance Survey, 1936 datum
# http://spatialreference.org/ref/epsg/4258/
self.etrs89=pyproj.Proj('+proj=longlat +ellps=GRS80 +no_defs')
def ll_to_bng(self,lon,lat):
x, y = pyproj.transform(self.etrs89, self.osgb36, lon, lat )
return x, y
def read_lat_lon(self,fname):
f=open(fname,'r')
sites=[]
for lin in f.readlines():
if lin.strip()[0] != '#' :
els=lin.strip().split(',')
if len(els)>0:
sites.append({})
sites[-1]['sitename']=els[0]
sites[-1]['lat']=float(els[1])
sites[-1]['lon']=float(els[2])
f.close()
return sites
################################################################################
# Parse the input
################################################################################
def parse_input(self):
parser=argparse.ArgumentParser(description='Get BNG square(s) containing pair(s) of lat/lon coords')
# positional
parser.add_argument('infname',help='File containing sitename, lat, lon elev')
parser.add_argument('outfname',help='Output file')
args=parser.parse_args()
return args.infname, \
args.outfname
################################################################################
################################################################################
#
# Start the main routine
#
################################################################################
################################################################################
if __name__=='__main__':
l=llbng()
infname, outfname = l.parse_input()
sites=l.read_lat_lon(infname)
for site in sites:
site['x'],site['y']=l.ll_to_bng(site['lon'],site['lat'])
of=open(outfname,'w')
for site in sites:
of.write('%s %f %f %f xx\n'%(site['sitename'].replace(' ','_').replace("'",''),site['x'],site['y'],0.0))
of.close()