srs.py 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365
  1. """
  2. The Spatial Reference class, represents OGR Spatial Reference objects.
  3. Example:
  4. >>> from django.contrib.gis.gdal import SpatialReference
  5. >>> srs = SpatialReference('WGS84')
  6. >>> print(srs)
  7. GEOGCS["WGS 84",
  8. DATUM["WGS_1984",
  9. SPHEROID["WGS 84",6378137,298.257223563,
  10. AUTHORITY["EPSG","7030"]],
  11. TOWGS84[0,0,0,0,0,0,0],
  12. AUTHORITY["EPSG","6326"]],
  13. PRIMEM["Greenwich",0,
  14. AUTHORITY["EPSG","8901"]],
  15. UNIT["degree",0.01745329251994328,
  16. AUTHORITY["EPSG","9122"]],
  17. AUTHORITY["EPSG","4326"]]
  18. >>> print(srs.proj)
  19. +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
  20. >>> print(srs.ellipsoid)
  21. (6378137.0, 6356752.3142451793, 298.25722356300003)
  22. >>> print(srs.projected, srs.geographic)
  23. False True
  24. >>> srs.import_epsg(32140)
  25. >>> print(srs.name)
  26. NAD83 / Texas South Central
  27. """
  28. from ctypes import byref, c_char_p, c_int
  29. from enum import IntEnum
  30. from types import NoneType
  31. from django.contrib.gis.gdal.base import GDALBase
  32. from django.contrib.gis.gdal.error import SRSException
  33. from django.contrib.gis.gdal.libgdal import GDAL_VERSION
  34. from django.contrib.gis.gdal.prototypes import srs as capi
  35. from django.utils.encoding import force_bytes, force_str
  36. class AxisOrder(IntEnum):
  37. TRADITIONAL = 0
  38. AUTHORITY = 1
  39. class SpatialReference(GDALBase):
  40. """
  41. A wrapper for the OGRSpatialReference object. According to the GDAL web site,
  42. the SpatialReference object "provide[s] services to represent coordinate
  43. systems (projections and datums) and to transform between them."
  44. """
  45. destructor = capi.release_srs
  46. def __init__(self, srs_input="", srs_type="user", axis_order=None):
  47. """
  48. Create a GDAL OSR Spatial Reference object from the given input.
  49. The input may be string of OGC Well Known Text (WKT), an integer
  50. EPSG code, a PROJ string, and/or a projection "well known" shorthand
  51. string (one of 'WGS84', 'WGS72', 'NAD27', 'NAD83').
  52. """
  53. if not isinstance(axis_order, (NoneType, AxisOrder)):
  54. raise ValueError(
  55. "SpatialReference.axis_order must be an AxisOrder instance."
  56. )
  57. self.axis_order = axis_order or AxisOrder.TRADITIONAL
  58. if srs_type == "wkt":
  59. self.ptr = capi.new_srs(c_char_p(b""))
  60. self.import_wkt(srs_input)
  61. if self.axis_order == AxisOrder.TRADITIONAL and GDAL_VERSION >= (3, 0):
  62. capi.set_axis_strategy(self.ptr, self.axis_order)
  63. elif self.axis_order != AxisOrder.TRADITIONAL and GDAL_VERSION < (3, 0):
  64. raise ValueError("%s is not supported in GDAL < 3.0." % self.axis_order)
  65. return
  66. elif isinstance(srs_input, str):
  67. try:
  68. # If SRID is a string, e.g., '4326', then make acceptable
  69. # as user input.
  70. srid = int(srs_input)
  71. srs_input = "EPSG:%d" % srid
  72. except ValueError:
  73. pass
  74. elif isinstance(srs_input, int):
  75. # EPSG integer code was input.
  76. srs_type = "epsg"
  77. elif isinstance(srs_input, self.ptr_type):
  78. srs = srs_input
  79. srs_type = "ogr"
  80. else:
  81. raise TypeError('Invalid SRS type "%s"' % srs_type)
  82. if srs_type == "ogr":
  83. # Input is already an SRS pointer.
  84. srs = srs_input
  85. else:
  86. # Creating a new SRS pointer, using the string buffer.
  87. buf = c_char_p(b"")
  88. srs = capi.new_srs(buf)
  89. # If the pointer is NULL, throw an exception.
  90. if not srs:
  91. raise SRSException(
  92. "Could not create spatial reference from: %s" % srs_input
  93. )
  94. else:
  95. self.ptr = srs
  96. if self.axis_order == AxisOrder.TRADITIONAL and GDAL_VERSION >= (3, 0):
  97. capi.set_axis_strategy(self.ptr, self.axis_order)
  98. elif self.axis_order != AxisOrder.TRADITIONAL and GDAL_VERSION < (3, 0):
  99. raise ValueError("%s is not supported in GDAL < 3.0." % self.axis_order)
  100. # Importing from either the user input string or an integer SRID.
  101. if srs_type == "user":
  102. self.import_user_input(srs_input)
  103. elif srs_type == "epsg":
  104. self.import_epsg(srs_input)
  105. def __getitem__(self, target):
  106. """
  107. Return the value of the given string attribute node, None if the node
  108. doesn't exist. Can also take a tuple as a parameter, (target, child),
  109. where child is the index of the attribute in the WKT. For example:
  110. >>> wkt = 'GEOGCS["WGS 84", DATUM["WGS_1984, ... AUTHORITY["EPSG","4326"]]'
  111. >>> srs = SpatialReference(wkt) # could also use 'WGS84', or 4326
  112. >>> print(srs['GEOGCS'])
  113. WGS 84
  114. >>> print(srs['DATUM'])
  115. WGS_1984
  116. >>> print(srs['AUTHORITY'])
  117. EPSG
  118. >>> print(srs['AUTHORITY', 1]) # The authority value
  119. 4326
  120. >>> print(srs['TOWGS84', 4]) # the fourth value in this wkt
  121. 0
  122. >>> # For the units authority, have to use the pipe symbole.
  123. >>> print(srs['UNIT|AUTHORITY'])
  124. EPSG
  125. >>> print(srs['UNIT|AUTHORITY', 1]) # The authority value for the units
  126. 9122
  127. """
  128. if isinstance(target, tuple):
  129. return self.attr_value(*target)
  130. else:
  131. return self.attr_value(target)
  132. def __str__(self):
  133. "Use 'pretty' WKT."
  134. return self.pretty_wkt
  135. # #### SpatialReference Methods ####
  136. def attr_value(self, target, index=0):
  137. """
  138. The attribute value for the given target node (e.g. 'PROJCS'). The index
  139. keyword specifies an index of the child node to return.
  140. """
  141. if not isinstance(target, str) or not isinstance(index, int):
  142. raise TypeError
  143. return capi.get_attr_value(self.ptr, force_bytes(target), index)
  144. def auth_name(self, target):
  145. "Return the authority name for the given string target node."
  146. return capi.get_auth_name(
  147. self.ptr, target if target is None else force_bytes(target)
  148. )
  149. def auth_code(self, target):
  150. "Return the authority code for the given string target node."
  151. return capi.get_auth_code(
  152. self.ptr, target if target is None else force_bytes(target)
  153. )
  154. def clone(self):
  155. "Return a clone of this SpatialReference object."
  156. return SpatialReference(capi.clone_srs(self.ptr), axis_order=self.axis_order)
  157. def from_esri(self):
  158. "Morph this SpatialReference from ESRI's format to EPSG."
  159. capi.morph_from_esri(self.ptr)
  160. def identify_epsg(self):
  161. """
  162. This method inspects the WKT of this SpatialReference, and will
  163. add EPSG authority nodes where an EPSG identifier is applicable.
  164. """
  165. capi.identify_epsg(self.ptr)
  166. def to_esri(self):
  167. "Morph this SpatialReference to ESRI's format."
  168. capi.morph_to_esri(self.ptr)
  169. def validate(self):
  170. "Check to see if the given spatial reference is valid."
  171. capi.srs_validate(self.ptr)
  172. # #### Name & SRID properties ####
  173. @property
  174. def name(self):
  175. "Return the name of this Spatial Reference."
  176. if self.projected:
  177. return self.attr_value("PROJCS")
  178. elif self.geographic:
  179. return self.attr_value("GEOGCS")
  180. elif self.local:
  181. return self.attr_value("LOCAL_CS")
  182. else:
  183. return None
  184. @property
  185. def srid(self):
  186. "Return the SRID of top-level authority, or None if undefined."
  187. try:
  188. return int(self.auth_code(target=None))
  189. except (TypeError, ValueError):
  190. return None
  191. # #### Unit Properties ####
  192. @property
  193. def linear_name(self):
  194. "Return the name of the linear units."
  195. units, name = capi.linear_units(self.ptr, byref(c_char_p()))
  196. return name
  197. @property
  198. def linear_units(self):
  199. "Return the value of the linear units."
  200. units, name = capi.linear_units(self.ptr, byref(c_char_p()))
  201. return units
  202. @property
  203. def angular_name(self):
  204. "Return the name of the angular units."
  205. units, name = capi.angular_units(self.ptr, byref(c_char_p()))
  206. return name
  207. @property
  208. def angular_units(self):
  209. "Return the value of the angular units."
  210. units, name = capi.angular_units(self.ptr, byref(c_char_p()))
  211. return units
  212. @property
  213. def units(self):
  214. """
  215. Return a 2-tuple of the units value and the units name. Automatically
  216. determine whether to return the linear or angular units.
  217. """
  218. units, name = None, None
  219. if self.projected or self.local:
  220. units, name = capi.linear_units(self.ptr, byref(c_char_p()))
  221. elif self.geographic:
  222. units, name = capi.angular_units(self.ptr, byref(c_char_p()))
  223. if name is not None:
  224. name = force_str(name)
  225. return (units, name)
  226. # #### Spheroid/Ellipsoid Properties ####
  227. @property
  228. def ellipsoid(self):
  229. """
  230. Return a tuple of the ellipsoid parameters:
  231. (semimajor axis, semiminor axis, and inverse flattening)
  232. """
  233. return (self.semi_major, self.semi_minor, self.inverse_flattening)
  234. @property
  235. def semi_major(self):
  236. "Return the Semi Major Axis for this Spatial Reference."
  237. return capi.semi_major(self.ptr, byref(c_int()))
  238. @property
  239. def semi_minor(self):
  240. "Return the Semi Minor Axis for this Spatial Reference."
  241. return capi.semi_minor(self.ptr, byref(c_int()))
  242. @property
  243. def inverse_flattening(self):
  244. "Return the Inverse Flattening for this Spatial Reference."
  245. return capi.invflattening(self.ptr, byref(c_int()))
  246. # #### Boolean Properties ####
  247. @property
  248. def geographic(self):
  249. """
  250. Return True if this SpatialReference is geographic
  251. (root node is GEOGCS).
  252. """
  253. return bool(capi.isgeographic(self.ptr))
  254. @property
  255. def local(self):
  256. "Return True if this SpatialReference is local (root node is LOCAL_CS)."
  257. return bool(capi.islocal(self.ptr))
  258. @property
  259. def projected(self):
  260. """
  261. Return True if this SpatialReference is a projected coordinate system
  262. (root node is PROJCS).
  263. """
  264. return bool(capi.isprojected(self.ptr))
  265. # #### Import Routines #####
  266. def import_epsg(self, epsg):
  267. "Import the Spatial Reference from the EPSG code (an integer)."
  268. capi.from_epsg(self.ptr, epsg)
  269. def import_proj(self, proj):
  270. """Import the Spatial Reference from a PROJ string."""
  271. capi.from_proj(self.ptr, proj)
  272. def import_user_input(self, user_input):
  273. "Import the Spatial Reference from the given user input string."
  274. capi.from_user_input(self.ptr, force_bytes(user_input))
  275. def import_wkt(self, wkt):
  276. "Import the Spatial Reference from OGC WKT (string)"
  277. capi.from_wkt(self.ptr, byref(c_char_p(force_bytes(wkt))))
  278. def import_xml(self, xml):
  279. "Import the Spatial Reference from an XML string."
  280. capi.from_xml(self.ptr, xml)
  281. # #### Export Properties ####
  282. @property
  283. def wkt(self):
  284. "Return the WKT representation of this Spatial Reference."
  285. return capi.to_wkt(self.ptr, byref(c_char_p()))
  286. @property
  287. def pretty_wkt(self, simplify=0):
  288. "Return the 'pretty' representation of the WKT."
  289. return capi.to_pretty_wkt(self.ptr, byref(c_char_p()), simplify)
  290. @property
  291. def proj(self):
  292. """Return the PROJ representation for this Spatial Reference."""
  293. return capi.to_proj(self.ptr, byref(c_char_p()))
  294. @property
  295. def proj4(self):
  296. "Alias for proj()."
  297. return self.proj
  298. @property
  299. def xml(self, dialect=""):
  300. "Return the XML representation of this Spatial Reference."
  301. return capi.to_xml(self.ptr, byref(c_char_p()), force_bytes(dialect))
  302. class CoordTransform(GDALBase):
  303. "The coordinate system transformation object."
  304. destructor = capi.destroy_ct
  305. def __init__(self, source, target):
  306. "Initialize on a source and target SpatialReference objects."
  307. if not isinstance(source, SpatialReference) or not isinstance(
  308. target, SpatialReference
  309. ):
  310. raise TypeError("source and target must be of type SpatialReference")
  311. self.ptr = capi.new_ct(source._ptr, target._ptr)
  312. self._srs1_name = source.name
  313. self._srs2_name = target.name
  314. def __str__(self):
  315. return 'Transform from "%s" to "%s"' % (self._srs1_name, self._srs2_name)