spatial_extent.py 56 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557155815591560156115621563156415651566156715681569157015711572157315741575157615771578157915801581158215831584158515861587158815891590159115921593159415951596159715981599160016011602160316041605160616071608160916101611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680168116821683168416851686168716881689169016911692169316941695169616971698169917001701170217031704170517061707170817091710171117121713171417151716171717181719172017211722172317241725172617271728172917301731173217331734173517361737173817391740174117421743174417451746174717481749175017511752175317541755175617571758175917601761176217631764176517661767176817691770177117721773177417751776177717781779178017811782178317841785178617871788178917901791179217931794179517961797179817991800180118021803180418051806180718081809181018111812181318141815181618171818181918201821182218231824182518261827182818291830183118321833183418351836183718381839184018411842184318441845184618471848184918501851185218531854
  1. """
  2. Spatial extents classes for map layer and space time datasets
  3. Usage:
  4. .. code-block:: python
  5. >>> import grass.temporal as tgis
  6. >>> tgis.init()
  7. >>> extent = tgis.RasterSpatialExtent(
  8. ... ident="raster@PERMANENT", north=90, south=90, east=180, west=180,
  9. ... top=100, bottom=-20)
  10. >>> extent = tgis.Raster3DSpatialExtent(
  11. ... ident="raster3d@PERMANENT", north=90, south=90, east=180, west=180,
  12. ... top=100, bottom=-20)
  13. >>> extent = tgis.VectorSpatialExtent(
  14. ... ident="vector@PERMANENT", north=90, south=90, east=180, west=180,
  15. ... top=100, bottom=-20)
  16. >>> extent = tgis.STRDSSpatialExtent(
  17. ... ident="strds@PERMANENT", north=90, south=90, east=180, west=180,
  18. ... top=100, bottom=-20)
  19. >>> extent = tgis.STR3DSSpatialExtent(
  20. ... ident="str3ds@PERMANENT", north=90, south=90, east=180, west=180,
  21. ... top=100, bottom=-20)
  22. >>> extent = tgis.STVDSSpatialExtent(
  23. ... ident="stvds@PERMANENT", north=90, south=90, east=180, west=180,
  24. ... top=100, bottom=-20)
  25. (C) 2012-2013 by the GRASS Development Team
  26. This program is free software under the GNU General Public
  27. License (>=v2). Read the file COPYING that comes with GRASS
  28. for details.
  29. :authors: Soeren Gebbert
  30. """
  31. from __future__ import print_function
  32. from .base import SQLDatabaseInterface
  33. from .core import init
  34. from datetime import datetime
  35. class SpatialExtent(SQLDatabaseInterface):
  36. """This is the spatial extent base class for all maps and space time datasets
  37. This class implements a three dimensional axis aligned bounding box
  38. and functions to compute topological relationships
  39. Usage:
  40. .. code-block:: python
  41. >>> init()
  42. >>> extent = SpatialExtent(table="raster_spatial_extent",
  43. ... ident="soil@PERMANENT", north=90, south=90, east=180, west=180,
  44. ... top=100, bottom=-20)
  45. >>> extent.id
  46. 'soil@PERMANENT'
  47. >>> extent.north
  48. 90.0
  49. >>> extent.south
  50. 90.0
  51. >>> extent.east
  52. 180.0
  53. >>> extent.west
  54. 180.0
  55. >>> extent.top
  56. 100.0
  57. >>> extent.bottom
  58. -20.0
  59. >>> extent.print_info()
  60. +-------------------- Spatial extent ----------------------------------------+
  61. | North:...................... 90.0
  62. | South:...................... 90.0
  63. | East:.. .................... 180.0
  64. | West:....................... 180.0
  65. | Top:........................ 100.0
  66. | Bottom:..................... -20.0
  67. >>> extent.print_shell_info()
  68. north=90.0
  69. south=90.0
  70. east=180.0
  71. west=180.0
  72. top=100.0
  73. bottom=-20.0
  74. """
  75. def __init__(self, table=None, ident=None, north=None, south=None,
  76. east=None, west=None, top=None, bottom=None, proj="XY"):
  77. SQLDatabaseInterface.__init__(self, table, ident)
  78. self.set_id(ident)
  79. self.set_spatial_extent_from_values(north, south, east, west, top,
  80. bottom)
  81. self.set_projection(proj)
  82. def overlapping_2d(self, extent):
  83. """Return True if this (A) and the provided spatial extent (B) overlaps
  84. in two dimensional space.
  85. Code is lend from wind_overlap.c in lib/gis
  86. Overlapping includes the spatial relations:
  87. - contain
  88. - in
  89. - cover
  90. - covered
  91. - equivalent
  92. .. code-block:: python
  93. >>> A = SpatialExtent(north=80, south=20, east=60, west=10)
  94. >>> B = SpatialExtent(north=80, south=20, east=60, west=10)
  95. >>> A.overlapping_2d(B)
  96. True
  97. :param extent: The spatial extent to check overlapping with
  98. :return: True or False
  99. """
  100. if self.get_projection() != extent.get_projection():
  101. self.msgr.error(_("Projections are different. Unable to compute "
  102. "overlapping_2d for spatial extents"))
  103. return False
  104. N = extent.get_north()
  105. S = extent.get_south()
  106. E = extent.get_east()
  107. W = extent.get_west()
  108. # Adjust the east and west in case of LL projection
  109. if self.get_projection() == "LL":
  110. while E < self.get_west():
  111. E += 360.0
  112. W += 360.0
  113. while W > self.get_east():
  114. E -= 360.0
  115. W -= 360.0
  116. if(self.get_north() <= S):
  117. return False
  118. if(self.get_south() >= N):
  119. return False
  120. if self.get_east() <= W:
  121. return False
  122. if self.get_west() >= E:
  123. return False
  124. return True
  125. def overlapping(self, extent):
  126. """Return True if this (A) and the provided spatial
  127. extent (B) overlaps in three dimensional space.
  128. Overlapping includes the spatial relations:
  129. - contain
  130. - in
  131. - cover
  132. - covered
  133. - equivalent
  134. Usage:
  135. .. code-block:: python
  136. >>> A = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  137. >>> B = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  138. >>> A.overlapping(B)
  139. True
  140. :param extent: The spatial extent to check overlapping with
  141. :return: True or False
  142. """
  143. if not self.overlapping_2d(extent):
  144. return False
  145. T = extent.get_top()
  146. B = extent.get_bottom()
  147. if self.get_top() <= B:
  148. return False
  149. if self.get_bottom() >= T:
  150. return False
  151. return True
  152. def intersect_2d(self, extent):
  153. """Return the two dimensional intersection as spatial_extent
  154. object or None in case no intersection was found.
  155. :param extent: The spatial extent to intersect with
  156. :return: The intersection spatial extent
  157. """
  158. if not self.overlapping_2d(extent):
  159. return None
  160. eN = extent.get_north()
  161. eS = extent.get_south()
  162. eE = extent.get_east()
  163. eW = extent.get_west()
  164. N = self.get_north()
  165. S = self.get_south()
  166. E = self.get_east()
  167. W = self.get_west()
  168. # Adjust the east and west in case of LL projection
  169. if self.get_projection() == "LL":
  170. while eE < W:
  171. eE += 360.0
  172. eW += 360.0
  173. while eW > E:
  174. eE -= 360.0
  175. eW -= 360.0
  176. # Compute the extent
  177. nN = N
  178. nS = S
  179. nE = E
  180. nW = W
  181. if W < eW:
  182. nW = eW
  183. if E > eE:
  184. nE = eE
  185. if N > eN:
  186. nN = eN
  187. if S < eS:
  188. nS = eS
  189. new = SpatialExtent(north=nN, south=nS, east=nE, west=nW,
  190. top=0, bottom=0, proj=self.get_projection())
  191. return new
  192. def intersect(self, extent):
  193. """Return the three dimensional intersection as spatial_extent
  194. object or None in case no intersection was found.
  195. Usage:
  196. .. code-block:: python
  197. >>> A = SpatialExtent(north=80, south=20, east=60, west=10,
  198. ... bottom=-50, top=50)
  199. >>> B = SpatialExtent(north=80, south=20, east=60, west=10,
  200. ... bottom=-50, top=50)
  201. >>> C = A.intersect(B)
  202. >>> C.print_info()
  203. +-------------------- Spatial extent ----------------------------------------+
  204. | North:...................... 80.0
  205. | South:...................... 20.0
  206. | East:.. .................... 60.0
  207. | West:....................... 10.0
  208. | Top:........................ 50.0
  209. | Bottom:..................... -50.0
  210. >>> B = SpatialExtent(north=40, south=30, east=60, west=10,
  211. ... bottom=-50, top=50)
  212. >>> C = A.intersect(B)
  213. >>> C.print_info()
  214. +-------------------- Spatial extent ----------------------------------------+
  215. | North:...................... 40.0
  216. | South:...................... 30.0
  217. | East:.. .................... 60.0
  218. | West:....................... 10.0
  219. | Top:........................ 50.0
  220. | Bottom:..................... -50.0
  221. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  222. ... bottom=-50, top=50)
  223. >>> C = A.intersect(B)
  224. >>> C.print_info()
  225. +-------------------- Spatial extent ----------------------------------------+
  226. | North:...................... 40.0
  227. | South:...................... 30.0
  228. | East:.. .................... 60.0
  229. | West:....................... 30.0
  230. | Top:........................ 50.0
  231. | Bottom:..................... -50.0
  232. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  233. ... bottom=-30, top=50)
  234. >>> C = A.intersect(B)
  235. >>> C.print_info()
  236. +-------------------- Spatial extent ----------------------------------------+
  237. | North:...................... 40.0
  238. | South:...................... 30.0
  239. | East:.. .................... 60.0
  240. | West:....................... 30.0
  241. | Top:........................ 50.0
  242. | Bottom:..................... -30.0
  243. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  244. ... bottom=-30, top=30)
  245. >>> C = A.intersect(B)
  246. >>> C.print_info()
  247. +-------------------- Spatial extent ----------------------------------------+
  248. | North:...................... 40.0
  249. | South:...................... 30.0
  250. | East:.. .................... 60.0
  251. | West:....................... 30.0
  252. | Top:........................ 30.0
  253. | Bottom:..................... -30.0
  254. :param extent: The spatial extent to intersect with
  255. :return: The intersection spatial extent
  256. """
  257. if not self.overlapping(extent):
  258. return None
  259. new = self.intersect_2d(extent)
  260. eT = extent.get_top()
  261. eB = extent.get_bottom()
  262. T = self.get_top()
  263. B = self.get_bottom()
  264. nT = T
  265. nB = B
  266. if B < eB:
  267. nB = eB
  268. if T > eT:
  269. nT = eT
  270. new.set_top(nT)
  271. new.set_bottom(nB)
  272. return new
  273. def union_2d(self, extent):
  274. """Return the two dimensional union as spatial_extent
  275. object or None in case the extents does not overlap or meet.
  276. :param extent: The spatial extent to create a union with
  277. :return: The union spatial extent
  278. """
  279. if not self.overlapping_2d(extent) and not self.meet_2d(extent):
  280. return None
  281. return self.disjoint_union_2d(extent)
  282. def disjoint_union_2d(self, extent):
  283. """Return the two dimensional union as spatial_extent.
  284. :param extent: The spatial extent to create a union with
  285. :return: The union spatial extent
  286. """
  287. eN = extent.get_north()
  288. eS = extent.get_south()
  289. eE = extent.get_east()
  290. eW = extent.get_west()
  291. N = self.get_north()
  292. S = self.get_south()
  293. E = self.get_east()
  294. W = self.get_west()
  295. # Adjust the east and west in case of LL projection
  296. if self.get_projection() == "LL":
  297. while eE < W:
  298. eE += 360.0
  299. eW += 360.0
  300. while eW > E:
  301. eE -= 360.0
  302. eW -= 360.0
  303. # Compute the extent
  304. nN = N
  305. nS = S
  306. nE = E
  307. nW = W
  308. if W > eW:
  309. nW = eW
  310. if E < eE:
  311. nE = eE
  312. if N < eN:
  313. nN = eN
  314. if S > eS:
  315. nS = eS
  316. new = SpatialExtent(north=nN, south=nS, east=nE, west=nW,
  317. top=0, bottom=0, proj=self.get_projection())
  318. return new
  319. def union(self, extent):
  320. """Return the three dimensional union as spatial_extent
  321. object or None in case the extents does not overlap or meet.
  322. :param extent: The spatial extent to create a union with
  323. :return: The union spatial extent
  324. """
  325. if not self.overlapping(extent) and not self.meet(extent):
  326. return None
  327. return self.disjoint_union(extent)
  328. def disjoint_union(self, extent):
  329. """Return the three dimensional union as spatial_extent .
  330. Usage:
  331. .. code-block:: python
  332. >>> A = SpatialExtent(north=80, south=20, east=60, west=10,
  333. ... bottom=-50, top=50)
  334. >>> B = SpatialExtent(north=80, south=20, east=60, west=10,
  335. ... bottom=-50, top=50)
  336. >>> C = A.disjoint_union(B)
  337. >>> C.print_info()
  338. +-------------------- Spatial extent ----------------------------------------+
  339. | North:...................... 80.0
  340. | South:...................... 20.0
  341. | East:.. .................... 60.0
  342. | West:....................... 10.0
  343. | Top:........................ 50.0
  344. | Bottom:..................... -50.0
  345. >>> B = SpatialExtent(north=40, south=30, east=60, west=10,
  346. ... bottom=-50, top=50)
  347. >>> C = A.disjoint_union(B)
  348. >>> C.print_info()
  349. +-------------------- Spatial extent ----------------------------------------+
  350. | North:...................... 80.0
  351. | South:...................... 20.0
  352. | East:.. .................... 60.0
  353. | West:....................... 10.0
  354. | Top:........................ 50.0
  355. | Bottom:..................... -50.0
  356. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  357. ... bottom=-50, top=50)
  358. >>> C = A.disjoint_union(B)
  359. >>> C.print_info()
  360. +-------------------- Spatial extent ----------------------------------------+
  361. | North:...................... 80.0
  362. | South:...................... 20.0
  363. | East:.. .................... 60.0
  364. | West:....................... 10.0
  365. | Top:........................ 50.0
  366. | Bottom:..................... -50.0
  367. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  368. ... bottom=-30, top=50)
  369. >>> C = A.disjoint_union(B)
  370. >>> C.print_info()
  371. +-------------------- Spatial extent ----------------------------------------+
  372. | North:...................... 80.0
  373. | South:...................... 20.0
  374. | East:.. .................... 60.0
  375. | West:....................... 10.0
  376. | Top:........................ 50.0
  377. | Bottom:..................... -50.0
  378. >>> B = SpatialExtent(north=40, south=30, east=60, west=30,
  379. ... bottom=-30, top=30)
  380. >>> C = A.disjoint_union(B)
  381. >>> C.print_info()
  382. +-------------------- Spatial extent ----------------------------------------+
  383. | North:...................... 80.0
  384. | South:...................... 20.0
  385. | East:.. .................... 60.0
  386. | West:....................... 10.0
  387. | Top:........................ 50.0
  388. | Bottom:..................... -50.0
  389. >>> A = SpatialExtent(north=80, south=20, east=60, west=10,
  390. ... bottom=-50, top=50)
  391. >>> B = SpatialExtent(north=90, south=80, east=70, west=20,
  392. ... bottom=-30, top=60)
  393. >>> C = A.disjoint_union(B)
  394. >>> C.print_info()
  395. +-------------------- Spatial extent ----------------------------------------+
  396. | North:...................... 90.0
  397. | South:...................... 20.0
  398. | East:.. .................... 70.0
  399. | West:....................... 10.0
  400. | Top:........................ 60.0
  401. | Bottom:..................... -50.0
  402. :param extent: The spatial extent to create a disjoint union with
  403. :return: The union spatial extent
  404. """
  405. new = self.disjoint_union_2d(extent)
  406. eT = extent.get_top()
  407. eB = extent.get_bottom()
  408. T = self.get_top()
  409. B = self.get_bottom()
  410. nT = T
  411. nB = B
  412. if B > eB:
  413. nB = eB
  414. if T < eT:
  415. nT = eT
  416. new.set_top(nT)
  417. new.set_bottom(nB)
  418. return new
  419. def is_in_2d(self, extent):
  420. """Return True if this extent (A) is located in the provided spatial
  421. extent (B) in two dimensions.
  422. ::
  423. _____
  424. |A _ |
  425. | |_| |
  426. |_____|B
  427. :param extent: The spatial extent
  428. :return: True or False
  429. """
  430. if self.get_projection() != extent.get_projection():
  431. self.msgr.error(_("Projections are different. Unable to compute "
  432. "is_in_2d for spatial extents"))
  433. return False
  434. eN = extent.get_north()
  435. eS = extent.get_south()
  436. eE = extent.get_east()
  437. eW = extent.get_west()
  438. N = self.get_north()
  439. S = self.get_south()
  440. E = self.get_east()
  441. W = self.get_west()
  442. # Adjust the east and west in case of LL projection
  443. if self.get_projection() == "LL":
  444. while eE < W:
  445. eE += 360.0
  446. eW += 360.0
  447. while eW > E:
  448. eE -= 360.0
  449. eW -= 360.0
  450. if W <= eW:
  451. return False
  452. if E >= eE:
  453. return False
  454. if N >= eN:
  455. return False
  456. if S <= eS:
  457. return False
  458. return True
  459. def is_in(self, extent):
  460. """Return True if this extent (A) is located in the provided spatial
  461. extent (B) in three dimensions.
  462. Usage:
  463. .. code-block:: python
  464. >>> A = SpatialExtent(north=79, south=21, east=59, west=11,
  465. ... bottom=-49, top=49)
  466. >>> B = SpatialExtent(north=80, south=20, east=60, west=10,
  467. ... bottom=-50, top=50)
  468. >>> A.is_in(B)
  469. True
  470. >>> B.is_in(A)
  471. False
  472. :param extent: The spatial extent
  473. :return: True or False
  474. """
  475. if not self.is_in_2d(extent):
  476. return False
  477. eT = extent.get_top()
  478. eB = extent.get_bottom()
  479. T = self.get_top()
  480. B = self.get_bottom()
  481. if B <= eB:
  482. return False
  483. if T >= eT:
  484. return False
  485. return True
  486. def contain_2d(self, extent):
  487. """Return True if this extent (A) contains the provided spatial
  488. extent (B) in two dimensions.
  489. Usage:
  490. .. code-block:: python
  491. >>> A = SpatialExtent(north=80, south=20, east=60, west=10)
  492. >>> B = SpatialExtent(north=79, south=21, east=59, west=11)
  493. >>> A.contain_2d(B)
  494. True
  495. >>> B.contain_2d(A)
  496. False
  497. :param extent: The spatial extent
  498. :return: True or False
  499. """
  500. return extent.is_in_2d(self)
  501. def contain(self, extent):
  502. """Return True if this extent (A) contains the provided spatial
  503. extent (B) in three dimensions.
  504. Usage:
  505. .. code-block:: python
  506. >>> A = SpatialExtent(north=80, south=20, east=60, west=10,
  507. ... bottom=-50, top=50)
  508. >>> B = SpatialExtent(north=79, south=21, east=59, west=11,
  509. ... bottom=-49, top=49)
  510. >>> A.contain(B)
  511. True
  512. >>> B.contain(A)
  513. False
  514. :param extent: The spatial extent
  515. :return: True or False
  516. """
  517. return extent.is_in(self)
  518. def equivalent_2d(self, extent):
  519. """Return True if this extent (A) is equal to the provided spatial
  520. extent (B) in two dimensions.
  521. Usage:
  522. .. code-block:: python
  523. >>> A = SpatialExtent(north=80, south=20, east=60, west=10)
  524. >>> B = SpatialExtent(north=80, south=20, east=60, west=10)
  525. >>> A.equivalent_2d(B)
  526. True
  527. >>> B.equivalent_2d(A)
  528. True
  529. :param extent: The spatial extent
  530. :return: True or False
  531. """
  532. if self.get_projection() != extent.get_projection():
  533. self.msgr.error(_("Projections are different. Unable to compute "
  534. "equivalent_2d for spatial extents"))
  535. return False
  536. eN = extent.get_north()
  537. eS = extent.get_south()
  538. eE = extent.get_east()
  539. eW = extent.get_west()
  540. N = self.get_north()
  541. S = self.get_south()
  542. E = self.get_east()
  543. W = self.get_west()
  544. # Adjust the east and west in case of LL projection
  545. if self.get_projection() == "LL":
  546. while eE < W:
  547. eE += 360.0
  548. eW += 360.0
  549. while eW > E:
  550. eE -= 360.0
  551. eW -= 360.0
  552. if W != eW:
  553. return False
  554. if E != eE:
  555. return False
  556. if N != eN:
  557. return False
  558. if S != eS:
  559. return False
  560. return True
  561. def equivalent(self, extent):
  562. """Return True if this extent (A) is equal to the provided spatial
  563. extent (B) in three dimensions.
  564. Usage:
  565. .. code-block:: python
  566. >>> A = SpatialExtent(north=80, south=20, east=60, west=10,
  567. ... bottom=-50, top=50)
  568. >>> B = SpatialExtent(north=80, south=20, east=60, west=10,
  569. ... bottom=-50, top=50)
  570. >>> A.equivalent(B)
  571. True
  572. >>> B.equivalent(A)
  573. True
  574. :param extent: The spatial extent
  575. :return: True or False
  576. """
  577. if not self.equivalent_2d(extent):
  578. return False
  579. eT = extent.get_top()
  580. eB = extent.get_bottom()
  581. T = self.get_top()
  582. B = self.get_bottom()
  583. if B != eB:
  584. return False
  585. if T != eT:
  586. return False
  587. return True
  588. def cover_2d(self, extent):
  589. """Return True if this extent (A) covers the provided spatial
  590. extent (B) in two dimensions.
  591. ::
  592. _____ _____ _____ _____
  593. |A __| |__ A| |A | B| |B | A|
  594. | |B | | B| | | |__| |__| |
  595. |__|__| |__|__| |_____| |_____|
  596. _____ _____ _____ _____
  597. |A|B| | |A __| |A _ | |__ A|
  598. | |_| | | |__|B | |B| | B|__| |
  599. |_____| |_____| |_|_|_| |_____|
  600. _____ _____ _____ _____
  601. |A|B | |_____|A |A|B|A| |_____|A
  602. | | | |B | | | | | |_____|B
  603. |_|___| |_____| |_|_|_| |_____|A
  604. The following cases are excluded:
  605. - contain
  606. - in
  607. - equivalent
  608. :param extent: The spatial extent
  609. :return: True or False
  610. """
  611. if self.get_projection() != extent.get_projection():
  612. self.msgr.error(_("Projections are different. Unable to compute"
  613. " cover_2d for spatial extents"))
  614. return False
  615. # Exclude equivalent_2d
  616. if self.equivalent_2d(extent):
  617. return False
  618. eN = extent.get_north()
  619. eS = extent.get_south()
  620. eE = extent.get_east()
  621. eW = extent.get_west()
  622. N = self.get_north()
  623. S = self.get_south()
  624. E = self.get_east()
  625. W = self.get_west()
  626. # Adjust the east and west in case of LL projection
  627. if self.get_projection() == "LL":
  628. while eE < W:
  629. eE += 360.0
  630. eW += 360.0
  631. while eW > E:
  632. eE -= 360.0
  633. eW -= 360.0
  634. # Edges of extent located outside of self are not allowed
  635. if E <= eW:
  636. return False
  637. if W >= eE:
  638. return False
  639. if N <= eS:
  640. return False
  641. if S >= eN:
  642. return False
  643. # First we check that at least one edge of extent meets an edge of self
  644. if W != eW and E != eE and N != eN and S != eS:
  645. return False
  646. # We check that at least one edge of extent is located in self
  647. edge_count = 0
  648. if W < eW and E > eW:
  649. edge_count += 1
  650. if E > eE and W < eE:
  651. edge_count += 1
  652. if N > eN and S < eN:
  653. edge_count += 1
  654. if S < eS and N > eS:
  655. edge_count += 1
  656. if edge_count == 0:
  657. return False
  658. return True
  659. def cover(self, extent):
  660. """Return True if this extent covers the provided spatial
  661. extent in three dimensions.
  662. The following cases are excluded:
  663. - contain
  664. - in
  665. - equivalent
  666. :param extent: The spatial extent
  667. :return: True or False
  668. """
  669. if self.get_projection() != extent.get_projection():
  670. self.msgr.error(_("Projections are different. Unable to compute "
  671. "cover for spatial extents"))
  672. return False
  673. # Exclude equivalent_2d
  674. if self.equivalent_2d(extent):
  675. return False
  676. eN = extent.get_north()
  677. eS = extent.get_south()
  678. eE = extent.get_east()
  679. eW = extent.get_west()
  680. eT = extent.get_top()
  681. eB = extent.get_bottom()
  682. N = self.get_north()
  683. S = self.get_south()
  684. E = self.get_east()
  685. W = self.get_west()
  686. T = self.get_top()
  687. B = self.get_bottom()
  688. # Adjust the east and west in case of LL projection
  689. if self.get_projection() == "LL":
  690. while eE < W:
  691. eE += 360.0
  692. eW += 360.0
  693. while eW > E:
  694. eE -= 360.0
  695. eW -= 360.0
  696. # Edges of extent located outside of self are not allowed
  697. if E <= eW:
  698. return False
  699. if W >= eE:
  700. return False
  701. if N <= eS:
  702. return False
  703. if S >= eN:
  704. return False
  705. if T <= eB:
  706. return False
  707. if B >= eT:
  708. return False
  709. # First we check that at least one edge of extent meets an edge of self
  710. if W != eW and E != eE and N != eN and S != eS and B != eB and T != eT:
  711. return False
  712. # We check that at least one edge of extent is located in self
  713. edge_count = 0
  714. if W < eW and E > eW:
  715. edge_count += 1
  716. if E > eE and W < eE:
  717. edge_count += 1
  718. if N > eN and S < eN:
  719. edge_count += 1
  720. if S < eS and N > eS:
  721. edge_count += 1
  722. if N > eN and S < eN:
  723. edge_count += 1
  724. if S < eS and N > eS:
  725. edge_count += 1
  726. if T > eT and B < eT:
  727. edge_count += 1
  728. if B < eB and T > eB:
  729. edge_count += 1
  730. if edge_count == 0:
  731. return False
  732. return True
  733. def covered_2d(self, extent):
  734. """Return True if this extent is covered by the provided spatial
  735. extent in two dimensions.
  736. The following cases are excluded:
  737. - contain
  738. - in
  739. - equivalent
  740. :param extent: The spatial extent
  741. :return: True or False
  742. """
  743. return extent.cover_2d(self)
  744. def covered(self, extent):
  745. """Return True if this extent is covered by the provided spatial
  746. extent in three dimensions.
  747. The following cases are excluded:
  748. - contain
  749. - in
  750. - equivalent
  751. :param extent: The spatial extent
  752. :return: True or False
  753. """
  754. return extent.cover(self)
  755. def overlap_2d(self, extent):
  756. """Return True if this extent (A) overlaps with the provided spatial
  757. extent (B) in two dimensions.
  758. Code is lend from wind_overlap.c in lib/gis
  759. ::
  760. _____
  761. |A __|__
  762. | | | B|
  763. |__|__| |
  764. |_____|
  765. The following cases are excluded:
  766. - contain
  767. - in
  768. - cover
  769. - covered
  770. - equivalent
  771. :param extent: The spatial extent
  772. :return: True or False
  773. """
  774. if self.contain_2d(extent):
  775. return False
  776. if self.is_in_2d(extent):
  777. return False
  778. if self.cover_2d(extent):
  779. return False
  780. if self.covered_2d(extent):
  781. return False
  782. if self.equivalent_2d(extent):
  783. return False
  784. N = extent.get_north()
  785. S = extent.get_south()
  786. E = extent.get_east()
  787. W = extent.get_west()
  788. # Adjust the east and west in case of LL projection
  789. if self.get_projection() == "LL":
  790. while E < self.get_west():
  791. E += 360.0
  792. W += 360.0
  793. while W > self.get_east():
  794. E -= 360.0
  795. W -= 360.0
  796. if(self.get_north() <= S):
  797. return False
  798. if(self.get_south() >= N):
  799. return False
  800. if self.get_east() <= W:
  801. return False
  802. if self.get_west() >= E:
  803. return False
  804. return True
  805. def overlap(self, extent):
  806. """Return True if this extent overlaps with the provided spatial
  807. extent in three dimensions.
  808. The following cases are excluded:
  809. - contain
  810. - in
  811. - cover
  812. - covered
  813. - equivalent
  814. :param extent: The spatial extent
  815. :return: True or False
  816. """
  817. if self.is_in(extent):
  818. return False
  819. if self.contain(extent):
  820. return False
  821. if self.cover(extent):
  822. return False
  823. if self.covered(extent):
  824. return False
  825. if self.equivalent(extent):
  826. return False
  827. N = extent.get_north()
  828. S = extent.get_south()
  829. E = extent.get_east()
  830. W = extent.get_west()
  831. T = extent.get_top()
  832. B = extent.get_bottom()
  833. # Adjust the east and west in case of LL projection
  834. if self.get_projection() == "LL":
  835. while E < self.get_west():
  836. E += 360.0
  837. W += 360.0
  838. while W > self.get_east():
  839. E -= 360.0
  840. W -= 360.0
  841. if(self.get_north() <= S):
  842. return False
  843. if(self.get_south() >= N):
  844. return False
  845. if self.get_east() <= W:
  846. return False
  847. if self.get_west() >= E:
  848. return False
  849. if self.get_top() <= B:
  850. return False
  851. if self.get_bottom() >= T:
  852. return False
  853. return True
  854. def meet_2d(self, extent):
  855. """Return True if this extent (A) meets with the provided spatial
  856. extent (B) in two dimensions.
  857. ::
  858. _____ _____
  859. | A | B |
  860. |_____| |
  861. |_____|
  862. _____ _____
  863. | B | A |
  864. | | |
  865. |_____|_____|
  866. ___
  867. | A |
  868. | |
  869. |___|
  870. | B |
  871. | |
  872. |_____|
  873. _____
  874. | B |
  875. | |
  876. |_____|_
  877. | A |
  878. | |
  879. |_____|
  880. :param extent: The spatial extent
  881. :return: True or False
  882. """
  883. eN = extent.get_north()
  884. eS = extent.get_south()
  885. eE = extent.get_east()
  886. eW = extent.get_west()
  887. N = self.get_north()
  888. S = self.get_south()
  889. E = self.get_east()
  890. W = self.get_west()
  891. # Adjust the east and west in case of LL projection
  892. if self.get_projection() == "LL":
  893. while eE < W:
  894. eE += 360.0
  895. eW += 360.0
  896. while eW > E:
  897. eE -= 360.0
  898. eW -= 360.0
  899. edge = None
  900. edge_count = 0
  901. if E == eW:
  902. edge = "E"
  903. edge_count += 1
  904. if W == eE:
  905. edge = "W"
  906. edge_count += 1
  907. if N == eS:
  908. edge = "N"
  909. edge_count += 1
  910. if S == eN:
  911. edge = "S"
  912. edge_count += 1
  913. # Meet a a single edge only
  914. if edge_count != 1:
  915. return False
  916. # Check boundaries of the faces
  917. if edge == "E" or edge == "W":
  918. if N < eS or S > eN:
  919. return False
  920. if edge == "N" or edge == "S":
  921. if E < eW or W > eE:
  922. return False
  923. return True
  924. def meet(self, extent):
  925. """Return True if this extent meets with the provided spatial
  926. extent in three dimensions.
  927. :param extent: The spatial extent
  928. :return: True or False
  929. """
  930. eN = extent.get_north()
  931. eS = extent.get_south()
  932. eE = extent.get_east()
  933. eW = extent.get_west()
  934. eT = extent.get_top()
  935. eB = extent.get_bottom()
  936. N = self.get_north()
  937. S = self.get_south()
  938. E = self.get_east()
  939. W = self.get_west()
  940. T = self.get_top()
  941. B = self.get_bottom()
  942. # Adjust the east and west in case of LL projection
  943. if self.get_projection() == "LL":
  944. while eE < W:
  945. eE += 360.0
  946. eW += 360.0
  947. while eW > E:
  948. eE -= 360.0
  949. eW -= 360.0
  950. edge = None
  951. edge_count = 0
  952. if E == eW:
  953. edge = "E"
  954. edge_count += 1
  955. if W == eE:
  956. edge = "W"
  957. edge_count += 1
  958. if N == eS:
  959. edge = "N"
  960. edge_count += 1
  961. if S == eN:
  962. edge = "S"
  963. edge_count += 1
  964. if T == eB:
  965. edge = "T"
  966. edge_count += 1
  967. if B == eT:
  968. edge = "B"
  969. edge_count += 1
  970. # Meet a single edge only
  971. if edge_count != 1:
  972. return False
  973. # Check boundaries of the faces
  974. if edge == "E" or edge == "W":
  975. if N < eS or S > eN:
  976. return False
  977. if T < eB or B > eT:
  978. return False
  979. if edge == "N" or edge == "S":
  980. if E < eW or W > eE:
  981. return False
  982. if T < eB or B > eT:
  983. return False
  984. if edge == "T" or edge == "B":
  985. if E < eW or W > eE:
  986. return False
  987. if N < eS or S > eN:
  988. return False
  989. return True
  990. def disjoint_2d(self, extent):
  991. """Return True if this extent (A) is disjoint with the provided spatial
  992. extent (B) in three dimensions.
  993. ::
  994. _____
  995. | A |
  996. |_____|
  997. _______
  998. | B |
  999. |_______|
  1000. :param extent: The spatial extent
  1001. :return: True or False
  1002. """
  1003. if self.is_in_2d(extent):
  1004. return False
  1005. if self.contain_2d(extent):
  1006. return False
  1007. if self.cover_2d(extent):
  1008. return False
  1009. if self.covered_2d(extent):
  1010. return False
  1011. if self.equivalent_2d(extent):
  1012. return False
  1013. if self.overlapping_2d(extent):
  1014. return False
  1015. if self.meet_2d(extent):
  1016. return False
  1017. return True
  1018. def disjoint(self, extent):
  1019. """Return True if this extent is disjoint with the provided spatial
  1020. extent in three dimensions.
  1021. :param extent: The spatial extent
  1022. :return: True or False
  1023. """
  1024. if self.is_in(extent):
  1025. return False
  1026. if self.contain(extent):
  1027. return False
  1028. if self.cover(extent):
  1029. return False
  1030. if self.covered(extent):
  1031. return False
  1032. if self.equivalent(extent):
  1033. return False
  1034. if self.overlapping(extent):
  1035. return False
  1036. if self.meet(extent):
  1037. return False
  1038. return True
  1039. def spatial_relation_2d(self, extent):
  1040. """Returns the two dimensional spatial relation between this
  1041. extent and the provided spatial extent in two dimensions.
  1042. Spatial relations are:
  1043. - disjoint
  1044. - meet
  1045. - overlap
  1046. - cover
  1047. - covered
  1048. - in
  1049. - contain
  1050. - equivalent
  1051. Usage: see self.spatial_relation()
  1052. """
  1053. if self.equivalent_2d(extent):
  1054. return "equivalent"
  1055. if self.contain_2d(extent):
  1056. return "contain"
  1057. if self.is_in_2d(extent):
  1058. return "in"
  1059. if self.cover_2d(extent):
  1060. return "cover"
  1061. if self.covered_2d(extent):
  1062. return "covered"
  1063. if self.overlap_2d(extent):
  1064. return "overlap"
  1065. if self.meet_2d(extent):
  1066. return "meet"
  1067. if self.disjoint_2d(extent):
  1068. return "disjoint"
  1069. return "unknown"
  1070. def spatial_relation(self, extent):
  1071. """Returns the two dimensional spatial relation between this
  1072. extent and the provided spatial extent in three dimensions.
  1073. Spatial relations are:
  1074. - disjoint
  1075. - meet
  1076. - overlap
  1077. - cover
  1078. - covered
  1079. - in
  1080. - contain
  1081. - equivalent
  1082. Usage:
  1083. .. code-block:: python
  1084. >>> A = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  1085. >>> B = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  1086. >>> A.spatial_relation(B)
  1087. 'equivalent'
  1088. >>> B.spatial_relation(A)
  1089. 'equivalent'
  1090. >>> B = SpatialExtent(north=70, south=20, east=60, west=10, bottom=-50, top=50)
  1091. >>> A.spatial_relation_2d(B)
  1092. 'cover'
  1093. >>> A.spatial_relation(B)
  1094. 'cover'
  1095. >>> B = SpatialExtent(north=70, south=30, east=60, west=10, bottom=-50, top=50)
  1096. >>> A.spatial_relation_2d(B)
  1097. 'cover'
  1098. >>> A.spatial_relation(B)
  1099. 'cover'
  1100. >>> B.spatial_relation_2d(A)
  1101. 'covered'
  1102. >>> B.spatial_relation(A)
  1103. 'covered'
  1104. >>> B = SpatialExtent(north=70, south=30, east=50, west=10, bottom=-50, top=50)
  1105. >>> A.spatial_relation_2d(B)
  1106. 'cover'
  1107. >>> B.spatial_relation_2d(A)
  1108. 'covered'
  1109. >>> A.spatial_relation(B)
  1110. 'cover'
  1111. >>> B = SpatialExtent(north=70, south=30, east=50, west=20, bottom=-50, top=50)
  1112. >>> B.spatial_relation(A)
  1113. 'covered'
  1114. >>> B = SpatialExtent(north=70, south=30, east=50, west=20, bottom=-50, top=50)
  1115. >>> A.spatial_relation_2d(B)
  1116. 'contain'
  1117. >>> A.spatial_relation(B)
  1118. 'cover'
  1119. >>> B = SpatialExtent(north=70, south=30, east=50, west=20, bottom=-40, top=50)
  1120. >>> A.spatial_relation(B)
  1121. 'cover'
  1122. >>> B = SpatialExtent(north=70, south=30, east=50, west=20, bottom=-40, top=40)
  1123. >>> A.spatial_relation(B)
  1124. 'contain'
  1125. >>> B.spatial_relation(A)
  1126. 'in'
  1127. >>> B = SpatialExtent(north=90, south=30, east=50, west=20, bottom=-40, top=40)
  1128. >>> A.spatial_relation_2d(B)
  1129. 'overlap'
  1130. >>> A.spatial_relation(B)
  1131. 'overlap'
  1132. >>> B = SpatialExtent(north=90, south=5, east=70, west=5, bottom=-40, top=40)
  1133. >>> A.spatial_relation_2d(B)
  1134. 'in'
  1135. >>> A.spatial_relation(B)
  1136. 'overlap'
  1137. >>> B = SpatialExtent(north=90, south=5, east=70, west=5, bottom=-40, top=60)
  1138. >>> A.spatial_relation(B)
  1139. 'overlap'
  1140. >>> B = SpatialExtent(north=90, south=5, east=70, west=5, bottom=-60, top=60)
  1141. >>> A.spatial_relation(B)
  1142. 'in'
  1143. >>> A = SpatialExtent(north=80, south=60, east=60, west=10, bottom=-50, top=50)
  1144. >>> B = SpatialExtent(north=60, south=20, east=60, west=10, bottom=-50, top=50)
  1145. >>> A.spatial_relation_2d(B)
  1146. 'meet'
  1147. >>> A.spatial_relation(B)
  1148. 'meet'
  1149. >>> A = SpatialExtent(north=60, south=40, east=60, west=10, bottom=-50, top=50)
  1150. >>> B = SpatialExtent(north=80, south=60, east=60, west=10, bottom=-50, top=50)
  1151. >>> A.spatial_relation_2d(B)
  1152. 'meet'
  1153. >>> A.spatial_relation(B)
  1154. 'meet'
  1155. >>> A = SpatialExtent(north=80, south=40, east=60, west=40, bottom=-50, top=50)
  1156. >>> B = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1157. >>> A.spatial_relation_2d(B)
  1158. 'meet'
  1159. >>> A.spatial_relation(B)
  1160. 'meet'
  1161. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1162. >>> B = SpatialExtent(north=90, south=30, east=60, west=40, bottom=-50, top=50)
  1163. >>> A.spatial_relation_2d(B)
  1164. 'meet'
  1165. >>> A.spatial_relation(B)
  1166. 'meet'
  1167. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1168. >>> B = SpatialExtent(north=70, south=50, east=60, west=40, bottom=-50, top=50)
  1169. >>> A.spatial_relation_2d(B)
  1170. 'meet'
  1171. >>> A.spatial_relation(B)
  1172. 'meet'
  1173. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1174. >>> B = SpatialExtent(north=60, south=20, east=60, west=40, bottom=-50, top=50)
  1175. >>> A.spatial_relation_2d(B)
  1176. 'meet'
  1177. >>> A.spatial_relation(B)
  1178. 'meet'
  1179. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1180. >>> B = SpatialExtent(north=40, south=20, east=60, west=40, bottom=-50, top=50)
  1181. >>> A.spatial_relation_2d(B)
  1182. 'disjoint'
  1183. >>> A.spatial_relation(B)
  1184. 'disjoint'
  1185. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1186. >>> B = SpatialExtent(north=60, south=20, east=60, west=40, bottom=-60, top=60)
  1187. >>> A.spatial_relation(B)
  1188. 'meet'
  1189. >>> A = SpatialExtent(north=80, south=40, east=40, west=20, bottom=-50, top=50)
  1190. >>> B = SpatialExtent(north=90, south=30, east=60, west=40, bottom=-40, top=40)
  1191. >>> A.spatial_relation(B)
  1192. 'meet'
  1193. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1194. >>> B = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1195. >>> A.spatial_relation(B)
  1196. 'meet'
  1197. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1198. >>> B = SpatialExtent(north=80, south=50, east=60, west=30, bottom=-50, top=0)
  1199. >>> A.spatial_relation(B)
  1200. 'meet'
  1201. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1202. >>> B = SpatialExtent(north=70, south=50, east=50, west=30, bottom=-50, top=0)
  1203. >>> A.spatial_relation(B)
  1204. 'meet'
  1205. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1206. >>> B = SpatialExtent(north=90, south=30, east=70, west=10, bottom=-50, top=0)
  1207. >>> A.spatial_relation(B)
  1208. 'meet'
  1209. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1210. >>> B = SpatialExtent(north=70, south=30, east=50, west=10, bottom=-50, top=0)
  1211. >>> A.spatial_relation(B)
  1212. 'meet'
  1213. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1214. >>> B = SpatialExtent(north=80, south=40, east=60, west=20, bottom=0, top=50)
  1215. >>> A.spatial_relation(B)
  1216. 'meet'
  1217. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1218. >>> B = SpatialExtent(north=80, south=50, east=60, west=30, bottom=0, top=50)
  1219. >>> A.spatial_relation(B)
  1220. 'meet'
  1221. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1222. >>> B = SpatialExtent(north=70, south=50, east=50, west=30, bottom=0, top=50)
  1223. >>> A.spatial_relation(B)
  1224. 'meet'
  1225. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1226. >>> B = SpatialExtent(north=90, south=30, east=70, west=10, bottom=0, top=50)
  1227. >>> A.spatial_relation(B)
  1228. 'meet'
  1229. >>> A = SpatialExtent(north=80, south=40, east=60, west=20, bottom=-50, top=0)
  1230. >>> B = SpatialExtent(north=70, south=30, east=50, west=10, bottom=0, top=50)
  1231. >>> A.spatial_relation(B)
  1232. 'meet'
  1233. >>> A = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  1234. >>> B = SpatialExtent(north=90, south=81, east=60, west=10, bottom=-50, top=50)
  1235. >>> A.spatial_relation(B)
  1236. 'disjoint'
  1237. >>> A = SpatialExtent(north=80, south=20, east=60, west=10, bottom=-50, top=50)
  1238. >>> B = SpatialExtent(north=90, south=80, east=60, west=10, bottom=-50, top=50)
  1239. >>> A.spatial_relation(B)
  1240. 'meet'
  1241. """
  1242. if self.equivalent(extent):
  1243. return "equivalent"
  1244. if self.contain(extent):
  1245. return "contain"
  1246. if self.is_in(extent):
  1247. return "in"
  1248. if self.cover(extent):
  1249. return "cover"
  1250. if self.covered(extent):
  1251. return "covered"
  1252. if self.overlap(extent):
  1253. return "overlap"
  1254. if self.meet(extent):
  1255. return "meet"
  1256. if self.disjoint(extent):
  1257. return "disjoint"
  1258. return "unknown"
  1259. def set_spatial_extent_from_values(self, north, south, east, west, top,
  1260. bottom):
  1261. """Set the three dimensional spatial extent
  1262. :param north: The northern edge
  1263. :param south: The southern edge
  1264. :param east: The eastern edge
  1265. :param west: The western edge
  1266. :param top: The top edge
  1267. :param bottom: The bottom edge
  1268. """
  1269. self.set_north(north)
  1270. self.set_south(south)
  1271. self.set_east(east)
  1272. self.set_west(west)
  1273. self.set_top(top)
  1274. self.set_bottom(bottom)
  1275. def set_spatial_extent(self, spatial_extent):
  1276. """Set the three dimensional spatial extent
  1277. :param spatial_extent: An object of type SpatialExtent or its
  1278. subclasses
  1279. """
  1280. self.set_north(spatial_extent.get_north())
  1281. self.set_south(spatial_extent.get_south())
  1282. self.set_east(spatial_extent.get_east())
  1283. self.set_west(spatial_extent.get_west())
  1284. self.set_top(spatial_extent.get_top())
  1285. self.set_bottom(spatial_extent.get_bottom())
  1286. def set_projection(self, proj):
  1287. """Set the projection of the spatial extent it should be XY or LL.
  1288. As default the projection is XY
  1289. """
  1290. if proj is None or (proj != "XY" and proj != "LL"):
  1291. self.D["proj"] = "XY"
  1292. else:
  1293. self.D["proj"] = proj
  1294. def set_spatial_extent_from_values_2d(self, north, south, east, west):
  1295. """Set the two dimensional spatial extent from values
  1296. :param north: The northern edge
  1297. :param south: The southern edge
  1298. :param east: The eastern edge
  1299. :param west: The western edge
  1300. """
  1301. self.set_north(north)
  1302. self.set_south(south)
  1303. self.set_east(east)
  1304. self.set_west(west)
  1305. def set_spatial_extent_2d(self, spatial_extent):
  1306. """Set the three dimensional spatial extent
  1307. :param spatial_extent: An object of type SpatialExtent or its
  1308. subclasses
  1309. """
  1310. self.set_north(spatial_extent.north)
  1311. self.set_south(spatial_extent.south)
  1312. self.set_east(spatial_extent.east)
  1313. self.set_west(spatial_extent.west)
  1314. def set_id(self, ident):
  1315. """Convenient method to set the unique identifier (primary key)"""
  1316. self.ident = ident
  1317. self.D["id"] = ident
  1318. def set_north(self, north):
  1319. """Set the northern edge of the map"""
  1320. if north is not None:
  1321. self.D["north"] = float(north)
  1322. else:
  1323. self.D["north"] = None
  1324. def set_south(self, south):
  1325. """Set the southern edge of the map"""
  1326. if south is not None:
  1327. self.D["south"] = float(south)
  1328. else:
  1329. self.D["south"] = None
  1330. def set_west(self, west):
  1331. """Set the western edge of the map"""
  1332. if west is not None:
  1333. self.D["west"] = float(west)
  1334. else:
  1335. self.D["west"] = None
  1336. def set_east(self, east):
  1337. """Set the eastern edge of the map"""
  1338. if east is not None:
  1339. self.D["east"] = float(east)
  1340. else:
  1341. self.D["east"] = None
  1342. def set_top(self, top):
  1343. """Set the top edge of the map"""
  1344. if top is not None:
  1345. self.D["top"] = float(top)
  1346. else:
  1347. self.D["top"] = None
  1348. def set_bottom(self, bottom):
  1349. """Set the bottom edge of the map"""
  1350. if bottom is not None:
  1351. self.D["bottom"] = float(bottom)
  1352. else:
  1353. self.D["bottom"] = None
  1354. def get_id(self):
  1355. """Convenient method to get the unique identifier (primary key)
  1356. :return: None if not found
  1357. """
  1358. if "id" in self.D:
  1359. return self.D["id"]
  1360. else:
  1361. return None
  1362. def get_projection(self):
  1363. """Get the projection of the spatial extent"""
  1364. return self.D["proj"]
  1365. def get_volume(self):
  1366. """Compute the volume of the extent, in case z is zero
  1367. (top == bottom or top - bottom = 1) the area is returned"""
  1368. if self.get_projection() == "LL":
  1369. self.msgr.error(_("Volume computation is not supported "
  1370. "for LL projections"))
  1371. area = self.get_area()
  1372. bbox = self.get_spatial_extent_as_tuple()
  1373. z = abs(bbox[4] - bbox[5])
  1374. if z == 0:
  1375. z = 1.0
  1376. return area * z
  1377. def get_area(self):
  1378. """Compute the area of the extent, extent in z direction is ignored"""
  1379. if self.get_projection() == "LL":
  1380. self.msgr.error(_("Area computation is not supported "
  1381. "for LL projections"))
  1382. bbox = self.get_spatial_extent_as_tuple()
  1383. y = abs(bbox[0] - bbox[1])
  1384. x = abs(bbox[2] - bbox[3])
  1385. return x * y
  1386. def get_spatial_extent_as_tuple(self):
  1387. """Return a tuple (north, south, east, west, top, bottom)
  1388. of the spatial extent"""
  1389. return (
  1390. self.north, self.south, self.east, self.west,
  1391. self.top, self.bottom)
  1392. def get_spatial_extent_as_tuple_2d(self):
  1393. """Return a tuple (north, south, east, west,) of the 2d spatial extent
  1394. """
  1395. return (self.north, self.south, self.east, self.west)
  1396. def get_north(self):
  1397. """Get the northern edge of the map
  1398. :return: None if not found"""
  1399. if "north" in self.D:
  1400. return self.D["north"]
  1401. else:
  1402. return None
  1403. def get_south(self):
  1404. """Get the southern edge of the map
  1405. :return: None if not found"""
  1406. if "south" in self.D:
  1407. return self.D["south"]
  1408. else:
  1409. return None
  1410. def get_east(self):
  1411. """Get the eastern edge of the map
  1412. :return: None if not found"""
  1413. if "east" in self.D:
  1414. return self.D["east"]
  1415. else:
  1416. return None
  1417. def get_west(self):
  1418. """Get the western edge of the map
  1419. :return: None if not found"""
  1420. if "west" in self.D:
  1421. return self.D["west"]
  1422. else:
  1423. return None
  1424. def get_top(self):
  1425. """Get the top edge of the map
  1426. :return: None if not found"""
  1427. if "top" in self.D:
  1428. return self.D["top"]
  1429. else:
  1430. return None
  1431. def get_bottom(self):
  1432. """Get the bottom edge of the map
  1433. :return: None if not found"""
  1434. if "bottom" in self.D:
  1435. return self.D["bottom"]
  1436. else:
  1437. return None
  1438. id = property(fget=get_id, fset=set_id)
  1439. north = property(fget=get_north, fset=set_north)
  1440. south = property(fget=get_south, fset=set_south)
  1441. east = property(fget=get_east, fset=set_east)
  1442. west = property(fget=get_west, fset=set_west)
  1443. top = property(fget=get_top, fset=set_top)
  1444. bottom = property(fget=get_bottom, fset=set_bottom)
  1445. def print_info(self):
  1446. """Print information about this class in human readable style"""
  1447. # 0123456789012345678901234567890
  1448. print(" +-------------------- Spatial extent ----------------------------------------+")
  1449. print(" | North:...................... " + str(self.get_north()))
  1450. print(" | South:...................... " + str(self.get_south()))
  1451. print(" | East:.. .................... " + str(self.get_east()))
  1452. print(" | West:....................... " + str(self.get_west()))
  1453. print(" | Top:........................ " + str(self.get_top()))
  1454. print(" | Bottom:..................... " + str(self.get_bottom()))
  1455. def print_shell_info(self):
  1456. """Print information about this class in shell style"""
  1457. print("north=" + str(self.get_north()))
  1458. print("south=" + str(self.get_south()))
  1459. print("east=" + str(self.get_east()))
  1460. print("west=" + str(self.get_west()))
  1461. print("top=" + str(self.get_top()))
  1462. print("bottom=" + str(self.get_bottom()))
  1463. ###############################################################################
  1464. class RasterSpatialExtent(SpatialExtent):
  1465. def __init__(self, ident=None, north=None, south=None, east=None,
  1466. west=None, top=None, bottom=None):
  1467. SpatialExtent.__init__(self, "raster_spatial_extent",
  1468. ident, north, south, east, west, top, bottom)
  1469. class Raster3DSpatialExtent(SpatialExtent):
  1470. def __init__(self, ident=None, north=None, south=None, east=None,
  1471. west=None, top=None, bottom=None):
  1472. SpatialExtent.__init__(self, "raster3d_spatial_extent",
  1473. ident, north, south, east, west, top, bottom)
  1474. class VectorSpatialExtent(SpatialExtent):
  1475. def __init__(self, ident=None, north=None, south=None, east=None,
  1476. west=None, top=None, bottom=None):
  1477. SpatialExtent.__init__(self, "vector_spatial_extent",
  1478. ident, north, south, east, west, top, bottom)
  1479. class STRDSSpatialExtent(SpatialExtent):
  1480. def __init__(self, ident=None, north=None, south=None, east=None,
  1481. west=None, top=None, bottom=None):
  1482. SpatialExtent.__init__(self, "strds_spatial_extent",
  1483. ident, north, south, east, west, top, bottom)
  1484. class STR3DSSpatialExtent(SpatialExtent):
  1485. def __init__(self, ident=None, north=None, south=None, east=None,
  1486. west=None, top=None, bottom=None):
  1487. SpatialExtent.__init__(self, "str3ds_spatial_extent",
  1488. ident, north, south, east, west, top, bottom)
  1489. class STVDSSpatialExtent(SpatialExtent):
  1490. def __init__(self, ident=None, north=None, south=None, east=None,
  1491. west=None, top=None, bottom=None):
  1492. SpatialExtent.__init__(self, "stvds_spatial_extent",
  1493. ident, north, south, east, west, top, bottom)
  1494. ###############################################################################
  1495. if __name__ == "__main__":
  1496. import doctest
  1497. doctest.testmod()