spatial_extent.py 55 KB

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