temporal_raster3d_algebra.py 4.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131
  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. @author Thomas Leppelt and Soeren Gebbert
  8. """
  9. import grass.pygrass.modules as pygrass
  10. from temporal_raster_base_algebra import *
  11. ###############################################################################
  12. class TemporalRaster3DAlgebraParser(TemporalRasterBaseAlgebraParser):
  13. """The temporal raster algebra class"""
  14. def __init__(self, pid=None, run=False, debug=True, spatial = False, nprocs = 1, register_null = False):
  15. TemporalRasterBaseAlgebraParser.__init__(self, pid, run, debug, spatial, nprocs, register_null)
  16. self.m_mapcalc = pymod.Module('r3.mapcalc')
  17. self.m_remove = pymod.Module('g.remove')
  18. def parse(self, expression, basename = None, overwrite=False):
  19. self.lexer = TemporalRasterAlgebraLexer()
  20. self.lexer.build()
  21. self.parser = yacc.yacc(module=self, debug=self.debug)
  22. self.overwrite = overwrite
  23. self.count = 0
  24. self.stdstype = "str3ds"
  25. self.basename = basename
  26. self.expression = expression
  27. self.parser.parse(expression)
  28. def remove_empty_maps(self):
  29. """! Removes the intermediate vector maps.
  30. """
  31. if self.empty_maps:
  32. self.msgr.message(_("Removing empty 3D raster maps"))
  33. namelist = self.empty_maps.values()
  34. max = 100
  35. chunklist = [namelist[i:i + max] for i in range(0, len(namelist), max)]
  36. for chunk in chunklist:
  37. stringlist = ",".join(chunk)
  38. if self.run:
  39. m = copy.deepcopy(self.m_remove)
  40. m.inputs["type"].value = "rast3d"
  41. m.inputs["pattern"].value = stringlist
  42. m.flags["f"].value = True
  43. m.run()
  44. ######################### Temporal functions ##############################
  45. def p_statement_assign(self, t):
  46. # The expression should always return a list of maps.
  47. """
  48. statement : stds EQUALS expr
  49. """
  50. TemporalRasterBaseAlgebraParser.p_statement_assign(self, t)
  51. def p_ts_neighbor_operation(self, t):
  52. # Examples:
  53. # A[1,0,-1]
  54. # B[-2]
  55. # C[1,-2,1,3]
  56. """
  57. expr : stds L_SPAREN number COMMA number COMMA number R_SPAREN
  58. | stds L_SPAREN number R_SPAREN
  59. | stds L_SPAREN number COMMA number COMMA number COMMA number R_SPAREN
  60. """
  61. # Check input stds.
  62. maplist = self.check_stds(t[1])
  63. t_neighbor = 0
  64. row_neighbor = None
  65. col_neighbor = None
  66. depth_neighbor = None
  67. if len(t) == 5:
  68. t_neighbor = t[3]
  69. elif len(t) == 9:
  70. row_neighbor = t[3]
  71. col_neighbor = t[5]
  72. depth_neighbor = t[7]
  73. elif len(t) == 11:
  74. t_neighbor = t[9]
  75. row_neighbor = t[3]
  76. col_neighbor = t[5]
  77. depth_neighbor = t[7]
  78. if self.run:
  79. resultlist = []
  80. max_index = len(maplist)
  81. for map_i in maplist:
  82. # Get map index and temporal extent.
  83. map_index = maplist.index(map_i)
  84. new_index = map_index + t_neighbor
  85. if new_index < max_index and new_index >= 0:
  86. map_i_t_extent = map_i.get_temporal_extent()
  87. # Get neighboring map and set temporal extent.
  88. map_n = maplist[new_index]
  89. # Generate an intermediate map for the result map list.
  90. map_new = self.generate_new_map(map_n, bool_op = 'and', copy = True)
  91. map_new.set_temporal_extent(map_i_t_extent)
  92. # Create r.mapcalc expression string for the operation.
  93. if "cmd_list" in dir(map_new) and len(t) == 5:
  94. cmdstring = "%s" %(map_new.cmd_list)
  95. elif "cmd_list" not in dir(map_new) and len(t) == 5:
  96. cmdstring = "%s" %(map_n.get_id())
  97. elif "cmd_list" in dir(map_new) and len(t) in (9,11):
  98. cmdstring = "%s[%s,%s,%s]" %(map_new.cmd_list, row_neighbor, col_neighbor, depth_neighbor)
  99. elif "cmd_list" not in dir(map_new) and len(t) in (9,11):
  100. cmdstring = "%s[%s,%s,%s]" %(map_n.get_id(), row_neighbor, col_neighbor, depth_neighbor)
  101. # Set new command list for map.
  102. map_new.cmd_list = cmdstring
  103. # Append map with updated command list to result list.
  104. resultlist.append(map_new)
  105. t[0] = resultlist
  106. ###############################################################################
  107. if __name__ == "__main__":
  108. import doctest
  109. doctest.testmod()