spatial_extent.py 55 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853
  1. """
  2. Spatial extents classes for map layer and space time datasets
  3. Usage:
  4. .. code-block:: python
  5. >>> import grass.temporal as tgis
  6. >>> tgis.init()
  7. >>> extent = tgis.RasterSpatialExtent(
  8. ... ident="raster@PERMANENT", north=90, south=90, east=180, west=180,
  9. ... top=100, bottom=-20)
  10. >>> extent = tgis.Raster3DSpatialExtent(
  11. ... ident="raster3d@PERMANENT", north=90, south=90, east=180, west=180,
  12. ... top=100, bottom=-20)
  13. >>> extent = tgis.VectorSpatialExtent(
  14. ... ident="vector@PERMANENT", north=90, south=90, east=180, west=180,
  15. ... top=100, bottom=-20)
  16. >>> extent = tgis.STRDSSpatialExtent(
  17. ... ident="strds@PERMANENT", north=90, south=90, east=180, west=180,
  18. ... top=100, bottom=-20)
  19. >>> extent = tgis.STR3DSSpatialExtent(
  20. ... ident="str3ds@PERMANENT", north=90, south=90, east=180, west=180,
  21. ... top=100, bottom=-20)
  22. >>> extent = tgis.STVDSSpatialExtent(
  23. ... ident="stvds@PERMANENT", north=90, south=90, east=180, west=180,
  24. ... top=100, bottom=-20)
  25. (C) 2012-2013 by the GRASS Development Team
  26. This program is free software under the GNU General Public
  27. License (>=v2). Read the file COPYING that comes with GRASS
  28. for details.
  29. :authors: Soeren Gebbert
  30. """
  31. from __future__ import print_function
  32. from .base import SQLDatabaseInterface
  33. 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()