temporal_raster3d_algebra.py 5.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164
  1. """!@package grass.temporal
  2. Temporal raster algebra
  3. (C) 2013 by the GRASS Development Team
  4. This program is free software under the GNU General Public
  5. License (>=v2). Read the file COPYING that comes with GRASS
  6. for details.
  7. :authors: Thomas Leppelt and Soeren Gebbert
  8. """
  9. from __future__ import print_function
  10. try:
  11. import ply.yacc as yacc
  12. except ImportError:
  13. pass
  14. from .temporal_raster_base_algebra import (
  15. TemporalRasterBaseAlgebraParser,
  16. TemporalRasterAlgebraLexer,
  17. )
  18. import grass.pygrass.modules as pymod
  19. from .space_time_datasets import Raster3DDataset
  20. class TemporalRaster3DAlgebraParser(TemporalRasterBaseAlgebraParser):
  21. """The temporal raster algebra class"""
  22. def __init__(
  23. self,
  24. pid=None,
  25. run=False,
  26. debug=True,
  27. spatial=False,
  28. register_null=False,
  29. dry_run=False,
  30. nprocs=1,
  31. ):
  32. TemporalRasterBaseAlgebraParser.__init__(
  33. self,
  34. pid=pid,
  35. run=run,
  36. debug=debug,
  37. spatial=spatial,
  38. register_null=register_null,
  39. dry_run=dry_run,
  40. nprocs=nprocs,
  41. )
  42. self.m_mapcalc = pymod.Module("r3.mapcalc")
  43. self.m_mremove = pymod.Module("g.remove")
  44. def parse(self, expression, basename=None, overwrite=False):
  45. # Check for space time dataset type definitions from temporal algebra
  46. l = TemporalRasterAlgebraLexer()
  47. l.build()
  48. l.lexer.input(expression)
  49. while True:
  50. tok = l.lexer.token()
  51. if not tok:
  52. break
  53. if tok.type == "STVDS" or tok.type == "STRDS" or tok.type == "STR3DS":
  54. raise SyntaxError("Syntax error near '%s'" % (tok.type))
  55. self.lexer = TemporalRasterAlgebraLexer()
  56. self.lexer.build()
  57. self.parser = yacc.yacc(module=self, debug=self.debug, write_tables=False)
  58. self.overwrite = overwrite
  59. self.count = 0
  60. self.stdstype = "str3ds"
  61. self.maptype = "raster_3d"
  62. self.mapclass = Raster3DDataset
  63. self.basename = basename
  64. self.expression = expression
  65. self.parser.parse(expression)
  66. return self.process_chain_dict
  67. def p_statement_assign(self, t):
  68. # The expression should always return a list of maps.
  69. """
  70. statement : stds EQUALS expr
  71. """
  72. TemporalRasterBaseAlgebraParser.p_statement_assign(self, t)
  73. def p_ts_neighbor_operation(self, t):
  74. # Examples:
  75. # A[1,0,-1]
  76. # B[-2]
  77. # C[1,-2,1,3]
  78. """
  79. expr : stds L_SPAREN number COMMA number COMMA number R_SPAREN
  80. | stds L_SPAREN number R_SPAREN
  81. | stds L_SPAREN number COMMA number COMMA number COMMA number R_SPAREN
  82. """
  83. # Check input stds.
  84. maplist = self.check_stds(t[1])
  85. t_neighbor = 0
  86. row_neighbor = None
  87. col_neighbor = None
  88. depth_neighbor = None
  89. if len(t) == 5:
  90. t_neighbor = t[3]
  91. elif len(t) == 9:
  92. row_neighbor = t[3]
  93. col_neighbor = t[5]
  94. depth_neighbor = t[7]
  95. elif len(t) == 11:
  96. t_neighbor = t[9]
  97. row_neighbor = t[3]
  98. col_neighbor = t[5]
  99. depth_neighbor = t[7]
  100. if self.run:
  101. resultlist = []
  102. max_index = len(maplist)
  103. for map_i in maplist:
  104. # Get map index and temporal extent.
  105. map_index = maplist.index(map_i)
  106. new_index = map_index + t_neighbor
  107. if new_index < max_index and new_index >= 0:
  108. map_i_t_extent = map_i.get_temporal_extent()
  109. # Get neighboring map and set temporal extent.
  110. map_n = maplist[new_index]
  111. # Generate an intermediate map for the result map list.
  112. map_new = self.generate_new_map(map_n, bool_op="and", copy=True)
  113. map_new.set_temporal_extent(map_i_t_extent)
  114. # Create r.mapcalc expression string for the operation.
  115. if "cmd_list" in dir(map_new) and len(t) == 5:
  116. cmdstring = "%s" % (map_new.cmd_list)
  117. elif "cmd_list" not in dir(map_new) and len(t) == 5:
  118. cmdstring = "%s" % (map_n.get_id())
  119. elif "cmd_list" in dir(map_new) and len(t) in (9, 11):
  120. cmdstring = "%s[%s,%s,%s]" % (
  121. map_new.cmd_list,
  122. row_neighbor,
  123. col_neighbor,
  124. depth_neighbor,
  125. )
  126. elif "cmd_list" not in dir(map_new) and len(t) in (9, 11):
  127. cmdstring = "%s[%s,%s,%s]" % (
  128. map_n.get_id(),
  129. row_neighbor,
  130. col_neighbor,
  131. depth_neighbor,
  132. )
  133. # Set new command list for map.
  134. map_new.cmd_list = cmdstring
  135. # Append map with updated command list to result list.
  136. resultlist.append(map_new)
  137. t[0] = resultlist
  138. if __name__ == "__main__":
  139. import doctest
  140. doctest.testmod()