spatial_extent.py 55 KB

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