1  """ Map Series  
  2   
  3  """ 
  4   
  5   
  6   
  7  from GP_Util import *  
  8   
  9 -def Import(study, varname, inshape, datefield, datetype, statefield, extract=True, exporttable=False, addtsf=False, adddates=True, ): 
  10          """Import a Map Series (such as vegetation maps) 
 11          This function imports a map series. A map series in this sense  
 12          is a polygon shapefile that represents discrete observations of  
 13          a value or state made at a specific places and times. There are  
 14          two output variables created for each map series which represent  
 15          the next and previous state observed at each sample point and  
 16          date. they are named with respect the the input varname as Pname  
 17          and Nname corresponding to the previous and next observations of  
 18          a series called name.  
 19           
 20      Inputs               
 21          varname         variable name for source data and output data 
 22          inshape         input polygon shapefile 
 23          datefield       name of date field 
 24          datetype        format of date field ( YYYY | YYYYMMDD ) 
 25          statefield      name of the state (mapped quantity) field 
 26   
 27      Options 
 28          adddates        add new unique dates to the study sample dates 
 29          addtsf          compute the time since fire at that mapdate 
 30          exporttable     export the variable specific data table """ 
 31          vardict = PrepMapSeries(study, varname, inshape, datefield, datetype, statefield, adddates) 
 32           
 33          PopulateIntervals(vardict['datagrid'] , 
 34                                            sortexpr = 'POINTID A; EVTDATE A' ) 
 35           
 36          if addtsf == True: 
 37                  AddFireAttributes(basegrid = vardict['datagrid'], 
 38                          addgrid         = study.infodict['FIRE']['datagrid'],  
 39                          cellsize        = study.cellsize) 
 40                   
 41          if exporttable == True: 
 42                  fields = vardict['fields'] 
 43                  if addtsf == True:  
 44                          fields.append('TSF') 
 45                           
 46                  ExportShapeTable(inshape=vardict['datagrid'],  
 47                          fieldlist       = fields,  
 48                          outfile         = varname, 
 49                          sortexpr        = 'POINTID A; EVTDATE A', ) 
  50   
 52          """ """ 
 53          vardict = study.infodict[varname] 
 54          (pdict,ndict) = ExtractRasters( 
 55                  varname = varname, 
 56                  datagrid = vardict['datagrid'], 
 57                  dates = study.enddates, 
 58                  cellsize = study.cellsize, ) 
 59   
 60          study.AddSeries(pdict) 
 61          study.AddSeries(ndict) 
  62           
 63 -def PrepMapSeries(study, varname, inshape, datefield, datetype, statefield,adddates=True): 
  64          """ """  
 65          log.info("Data Prep ("+varname+")") 
 66           
 67          workingdir = os.path.join(study.workingdir,varname) 
 68          MakeFolder(workingdir) 
 69          os.chdir(workingdir) 
 70          gp.Workspace = workingdir 
 71           
 72          datashp = os.path.abspath('data.shp') 
 73          datagrid = os.path.abspath('datagrid.shp') 
 74   
 75           
 76          gp.copy(inshape,datashp) 
 77          CopyField(datashp,'EVTDATE',datefield) 
 78          CopyField(datashp,'STATE',statefield) 
 79          ConvertDateField(datashp,'EVTDATE',datetype) 
 80           
 81           
 82          gp.Append_management(study.studyshp, datashp, "NO_TEST") 
 83          gp.Append_management(study.finstudyshp, datashp, "NO_TEST") 
 84   
 85           
 86          gp.Intersect(study.pointgrid+';'+ datashp, datagrid) 
 87           
 88           
 89          gp.AddField_management(datagrid,'PREVDATE',"DOUBLE") 
 90          gp.AddField_management(datagrid,'PREVSTATE',"DOUBLE") 
 91           
 92           
 93          fields = ['POINTID', 'EVTDATE', 'STATE', 'PREVDATE', 'PREVSTATE']  
 94                   
 95          uniquedates = num.unique(GetShpField(datashp,'EVTDATE')) 
 96          if adddates == True: 
 97                  study.UpdateDates(uniquedates) 
 98           
 99           
100          vardict = { 
101                  "type"                  : 'map', 
102                  "varname"               : varname , 
103                  "workingdir"    : workingdir , 
104                  "datashape"     : datashp , 
105                  "datagrid"              : datagrid , 
106                  "uniquedates"   : uniquedates, 
107                  "fields"                : fields, 
108                  } 
109                   
110          study.AddSource(varname,vardict) 
111          return vardict 
 112   
114          """Populate the point-events with  
115                  PREVSTATE, PREVDATE 
116          """ 
117          log.info("Populate Intervals") 
118           
119          rows = gp.Updatecursor(datafile,"","","",sortexpr) 
120          row = rows.Next() 
121          prevpid = -1 
122          while row: 
123                   
124                  if prevpid != row.POINTID: 
125                          prevstate = -9999 
126                          prevdate = -9999 
127                           
128                   
129                  row.PREVSTATE = prevstate        
130                  row.PREVDATE = prevdate 
131                   
132                   
133                  prevstate = row.STATE 
134                  prevdate = row.EVTDATE 
135                  prevpid = row.POINTID 
136                   
137                  rows.UpdateRow(row) 
138                  row = rows.Next() 
 139   
141          """Add Time Since Fire to each map in a map series""" 
142          log.info("Compute Map Series TSF") 
143          os.chdir(os.path.dirname(basegrid)) 
144          gp.Workspace = os.getcwd() 
145           
146          vegdates = num.unique(GetShpField(basegrid,'EVTDATE')) 
147          gp.MakeFeatureLayer(basegrid,"map") 
148          gp.MakeFeatureLayer(addgrid,"fire") 
149           
150          for date in vegdates: 
151                   
152                  gp.SelectLayerByAttribute("map", "NEW_SELECTION",'"EVTDATE" = '+str(date) ) 
153                  gp.CopyFeatures("map","submap.shp") 
154                           
155                   
156                  gp.SelectLayerByAttribute("fire", "NEW_SELECTION",'"EVTDATE" < '+str(date) ) 
157                   
158                   
159                  gp.PointToRaster_conversion('fire',"EVTDATE", 'daterast', 'MAXIMUM', "NONE", cellsize) 
160                  gp.Minus_sa(date,'daterast','tsfrast') 
161                   
162                  gp.ExtractValuesToPoints_sa("submap.shp",  'tsfrast', 'tsfsub.shp', "NONE", "VALUE_ONLY") 
163                   
164                  if date == vegdates[0]: 
165                          gp.copy('tsfsub.shp','tsfmap.shp') 
166                  else: 
167                          gp.append('tsfsub.shp','tsfmap.shp') 
168           
169           
170          gp.Delete('daterast') 
171          gp.Delete('tsfrast') 
172          gp.Delete('tsfsub.shp') 
173          gp.Delete('submap.shp') 
174          gp.Delete('map')   
175          gp.Delete('fire') 
176           
177          CopyField('tsfmap.shp',"TSF","RASTERVALU") 
178          gp.copy('tsfmap.shp', basegrid) 
 179   
181          os.chdir(os.path.dirname(datagrid)) 
182          gp.Workspace = MakeFolder(os.path.abspath('BURN')) 
183           
184          gp.MakeFeatureLayer(datagrid,"lyr") 
185          pfiles = [os.path.abspath('prev'+str(d)) for d in range(len(dates))] 
186          nfiles = [os.path.abspath('next'+str(d)) for d in range(len(dates))] 
187          prog = ProgDots(len(dates),title='Extract Map Series') 
188          for (date,pfile,nfile) in zip(dates,pfiles,nfiles): 
189                   
190                  gp.SelectLayerByAttribute("lyr", "NEW_SELECTION",'"PREVDATE" < ' + str(date) ) 
191                  gp.SelectLayerByAttribute("lyr", "SUBSET_SELECTION", '"EVTDATE" >= ' + str(date) ) 
192                   
193                  gp.PointToRaster_conversion('lyr',"PREVSTATE", pfile, 'MAXIMUM', "NONE", cellsize) 
194                   
195                  gp.PointToRaster_conversion('lyr',"STATE", nfile, 'MAXIMUM', "NONE", cellsize) 
196                  prog.step() 
197   
198           
199          prevNpy = RasterToNpySeries(pfiles,cleanup=False) 
200          nextNpy = RasterToNpySeries(nfiles,cleanup=False) 
201           
202          pfile = JoinDateRasters(prevNpy,'P'+varname) 
203          nfile = JoinDateRasters(nextNpy,'N'+varname) 
204           
205          pdict = {'name':'P'+varname,'file':pfile,'list':prevNpy} 
206          ndict = {'name':'N'+varname,'file':nfile,'list':nextNpy} 
207           
208          return (pdict, ndict) 
209