spatial_extent.py 53 KB

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