1# (C) British Crown Copyright 2014 - 2018, Met Office 2# 3# This file is part of pyepsg. 4# 5# pyepsg is free software: you can redistribute it and/or modify it under 6# the terms of the GNU Lesser General Public License as published by the 7# Free Software Foundation, either version 3 of the License, or 8# (at your option) any later version. 9# 10# pyepsg is distributed in the hope that it will be useful, 11# but WITHOUT ANY WARRANTY; without even the implied warranty of 12# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13# GNU Lesser General Public License for more details. 14# 15# You should have received a copy of the GNU Lesser General Public License 16# along with pyepsg. If not, see <http://www.gnu.org/licenses/>. 17""" 18Provides simple access to http://epsg.io/. 19 20The entry point for this package is the :func:`get()` function. 21 22""" 23from __future__ import print_function 24 25import sys 26import weakref 27import xml.etree.ElementTree as ET 28 29import requests 30 31 32EPSG_IO_URL = 'http://epsg.io/' 33 34GML_NS = '{http://www.opengis.net/gml/3.2}' 35XLINK_NS = '{http://www.w3.org/1999/xlink}' 36 37 38_cache = weakref.WeakValueDictionary() 39 40 41class EPSG(object): 42 """Parent class of all objects returned by pyepsg.""" 43 def __init__(self, element): 44 self.element = element 45 46 @property 47 def identifier(self): 48 """The official URN for this object.""" 49 return self.element.find(GML_NS + 'identifier').text 50 51 52class UOM(EPSG): 53 """A unit of measure.""" 54 @property 55 def name(self): 56 """The human-readable name.""" 57 return self.element.find(GML_NS + 'name').text 58 59 60class Axis(EPSG): 61 """A single coordinate axis.""" 62 63 @property 64 def direction(self): 65 """A description of the orientation of this axis.""" 66 return self.element.find(GML_NS + 'axisDirection').text 67 68 @property 69 def uom(self): 70 """The name of the unit of measure used on this axis.""" 71 uom = self.element.attrib['uom'] 72 return get(uom).name 73 74 def __repr__(self): 75 return '<Axis: {self.direction} / {self.uom}>'.format(self=self) 76 77 78class CartesianCS(EPSG): 79 """A 1-, 2-, or 3-dimensional cartesian coordinate system.""" 80 @property 81 def axes(self): 82 """An ordered list of :class:`Axis` objects describing X and Y.""" 83 axes = self.element.findall(GML_NS + 'axis') 84 return [Axis(axis.find(GML_NS + 'CoordinateSystemAxis')) for 85 axis in axes] 86 87 @property 88 def name(self): 89 """The human-readable name.""" 90 return self.element.find(GML_NS + 'name').text 91 92 @property 93 def remarks(self): 94 """Human-readable comments.""" 95 return self.element.find(GML_NS + 'remarks').text 96 97 def __repr__(self): 98 name = self.name 99 if len(name) > 38: 100 name = name[:38] + '..' 101 return '<CartesianCS: {name}>'.format(name=name) 102 103 104class CRS(EPSG): 105 """ 106 Abstract parent class for :class:`GeodeticCRS`, :class:`ProjectedCRS` 107 and :class:`CompoundCRS`. 108 109 """ 110 111 @property 112 def id(self): 113 """The EPSG code for this CRS.""" 114 id = self.element.attrib[GML_NS + 'id'] 115 code = id.split('-')[-1] 116 return code 117 118 @property 119 def name(self): 120 """The human-readable name.""" 121 return self.element.find(GML_NS + 'name').text 122 123 @property 124 def scope(self): 125 """A human-readable description of the intended usage for this CRS.""" 126 return self.element.find(GML_NS + 'scope').text 127 128 def as_esri_wkt(self): 129 """ 130 Return the ESRI WKT string which corresponds to the CRS. 131 132 For example:: 133 134 >>> print(get(27700).as_esri_wkt()) # doctest: +ELLIPSIS 135 PROJCS["OSGB_1936_British_National_Grid",GEOGCS["GCS_OSGB 19... 136 137 """ 138 url = '{prefix}{code}.esriwkt?download'.format(prefix=EPSG_IO_URL, 139 code=self.id) 140 return requests.get(url).text 141 142 def as_html(self): 143 """ 144 Return the OGC WKT which corresponds to the CRS as HTML. 145 146 For example:: 147 148 >>> print(get(27700).as_html()) # doctest: +ELLIPSIS 149 <div class="syntax"><pre><span class="gh">PROJCS</span><span... 150 151 """ 152 url = '{prefix}{code}.html?download'.format(prefix=EPSG_IO_URL, 153 code=self.id) 154 return requests.get(url).text 155 156 def as_proj4(self): 157 """ 158 Return the PROJ.4 string which corresponds to the CRS. 159 160 For example:: 161 162 >>> print(get(21781).as_proj4()) 163 +proj=somerc +lat_0=46.95240555555556 +lon_0=7.439583333333333 \ 164+k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel \ 165+towgs84=674.4,15.1,405.3,0,0,0,0 +units=m +no_defs 166 167 """ 168 url = '{prefix}{code}.proj4?download'.format(prefix=EPSG_IO_URL, 169 code=self.id) 170 return requests.get(url).text.strip() 171 172 def as_wkt(self): 173 """ 174 Return the OGC WKT string which corresponds to the CRS. 175 176 For example:: 177 178 >>> print(get(27700).as_wkt()) # doctest: +ELLIPSIS 179 PROJCS["OSGB 1936 / British National Grid",GEOGCS["OSGB 1936... 180 181 """ 182 url = '{prefix}{code}.wkt?download'.format(prefix=EPSG_IO_URL, 183 code=self.id) 184 return requests.get(url).text 185 186 def domain_of_validity(self): 187 """ 188 Return the domain of validity for this CRS as: 189 (west, east, south, north). 190 191 For example:: 192 193 >>> print(get(21781).domain_of_validity()) 194 [5.96, 10.49, 45.82, 47.81] 195 196 197 """ 198 # TODO: Generalise interface to return a polygon? (Can we find 199 # something that uses a polygon instead?) 200 domain = self.element.find(GML_NS + 'domainOfValidity') 201 domain_href = domain.attrib[XLINK_NS + 'href'] 202 url = '{prefix}{code}.gml?download'.format(prefix=EPSG_IO_URL, 203 code=domain_href) 204 xml = requests.get(url).content 205 gml = ET.fromstring(xml) 206 207 def extract_bound(tag): 208 ns = '{http://www.isotc211.org/2005/gmd}' 209 xpath = './/{ns}EX_GeographicBoundingBox/{ns}{tag}/'.format( 210 ns=ns, 211 tag=tag) 212 bound = gml.find(xpath) 213 return float(bound.text) 214 215 tags = ('westBoundLongitude', 'eastBoundLongitude', 216 'southBoundLatitude', 'northBoundLatitude') 217 bounds = [extract_bound(tag) for tag in tags] 218 return bounds 219 220 221class GeodeticCRS(CRS): 222 """ 223 Represents a single geodetic CRS. 224 225 """ 226 227 def __repr__(self): 228 return '<GeodeticCRS: {self.id}, {self.name}>'.format(self=self) 229 230 231class ProjectedCRS(CRS): 232 """ 233 Represents a single projected CRS. 234 235 """ 236 @property 237 def base_geodetic_crs(self): 238 """The :class:`GeodeticCRS` on which this projection is based.""" 239 base = self.element.find(GML_NS + 'baseGeodeticCRS') 240 href = base.attrib[XLINK_NS + 'href'] 241 return get(href) 242 243 @property 244 def cartesian_cs(self): 245 """The :class:`CartesianCS` which describes the coordinate axes.""" 246 cs = self.element.find(GML_NS + 'cartesianCS') 247 href = cs.attrib[XLINK_NS + 'href'] 248 return get(href) 249 250 def __repr__(self): 251 return '<ProjectedCRS: {self.id}, {self.name}>'.format(self=self) 252 253 254class CompoundCRS(CRS): 255 """ 256 Represents a single compound CRS. 257 258 """ 259 def __repr__(self): 260 return '<CompoundCRS: {self.id}, {self.name}>'.format(self=self) 261 262 263def get(code): 264 """ 265 Return an object that corresponds to the given EPSG code. 266 267 Currently supported object types are: 268 - :class:`GeodeticCRS` 269 - :class:`ProjectedCRS` 270 - :class:`CartesianCS` 271 - :class:`UOM` 272 273 For example:: 274 275 >>> print(get(27700)) 276 <ProjectedCRS: 27700, OSGB 1936 / British National Grid> 277 >>> print(get('4400-cs')) 278 <CartesianCS: Cartesian 2D CS. Axes: easting, northi..> 279 >>> print(get(5973)) 280 <CompoundCRS: 5973, ETRS89 / UTM zone 33 + NN2000 height> 281 282 """ 283 instance = _cache.get(code) 284 if instance is None: 285 url = '{prefix}{code}.gml?download'.format(prefix=EPSG_IO_URL, 286 code=code) 287 xml = requests.get(url).content 288 root = ET.fromstring(xml) 289 class_for_tag = { 290 GML_NS + 'CartesianCS': CartesianCS, 291 GML_NS + 'GeodeticCRS': GeodeticCRS, 292 GML_NS + 'ProjectedCRS': ProjectedCRS, 293 GML_NS + 'CompoundCRS': CompoundCRS, 294 GML_NS + 'BaseUnit': UOM, 295 } 296 if root.tag in class_for_tag: 297 instance = class_for_tag[root.tag](root) 298 else: 299 raise ValueError('Unsupported code type: {}'.format(root.tag)) 300 _cache[code] = instance 301 return instance 302 303 304if __name__ == '__main__': 305 import doctest 306 failure_count, test_count = doctest.testmod() 307 sys.exit(failure_count) 308