spatial_extent.py 56 KB

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