kmz_creation.py

   1 import sys, os
   2 import datetime as dt
   3 import numpy as np
   4 import saveload as sl
   5 from utils import Dict
   6 from JPSSD import retrieve_af_data, time_iso2num
   7 from interpolation import sort_dates
   8 from mpl_toolkits.basemap import Basemap
   9 from plot_pixels import basemap_scatter_mercator, create_kml
  10
  11 only_fire = True
  12 dlon = .01
  13 dlat = .005
  14
  15 def exist(path):
  16         return (os.path.isfile(path) and os.access(path,os.R_OK))
  17
  18 print ''
  19 print '>> Reading the arguments <<'
  20 dti = dt.datetime.strptime(sys.argv[1],'%Y%m%d%H%M%S')
  21 time_start_iso = '%d-%02d-%02dT%02d:%02d:%02dZ' % (dti.year,dti.month,dti.day,dti.hour,dti.minute,dti.second)
  22 dtf = dti+dt.timedelta(days=float(sys.argv[2]))
  23 time_final_iso = '%d-%02d-%02dT%02d:%02d:%02dZ' % (dtf.year,dtf.month,dtf.day,dtf.hour,dtf.minute,dtf.second)
  24 time_iso = (time_start_iso,time_final_iso)
  25
  26 bounds = (float(sys.argv[3]),float(sys.argv[4]),float(sys.argv[5]),float(sys.argv[6]))
  27 fxlon,fxlat = np.meshgrid(np.arange(bounds[0],bounds[1],dlon),
  28             np.arange(bounds[2],bounds[3],dlat))
  29 maxsize = 500
  30 coarsening=np.int(1+np.max(fxlon.shape)/maxsize)
  31 fxlon = fxlon[0::coarsening,0::coarsening]
  32 fxlat = fxlat[0::coarsening,0::coarsening]
  33 bbox = (fxlon.min(),fxlon.max(),fxlat.min(),fxlat.max())
  34 time_num = map(time_iso2num,time_iso)
  35
  36 print ''
  37 print '>> Retrieving satellite data <<'
  38 sys.stdout.flush()
  39 data = retrieve_af_data(bbox,time_iso)
  40
  41 print ''
  42 print '>> Saving data file <<'
  43 sl.save((data,fxlon,fxlat,time_num),'data')
  44 # sort the granules by dates
  45 sdata=sort_dates(data)
  46
  47 print ''
  48 print '>> Creating KMZ file <<'
  49 # creating KMZ overlay of each information
  50 # create the Basemap to plot into
  51 bmap = Basemap(projection='merc',llcrnrlat=bbox[2], urcrnrlat=bbox[3], llcrnrlon=bbox[0], urcrnrlon=bbox[1])
  52 # initialize array
  53 kmld = []
  54 # for each observed information
  55 for idx, g in enumerate(sdata):
  56     # create png name
  57     pngfile = g[0]+'.png'
  58     # create timestamp for KML
  59     timestamp = g[1].acq_date + 'T' + g[1].acq_time[0:2] + ':' + g[1].acq_time[2:4] + 'Z'
  60     if not exist(pngfile):
  61         # plot a scatter basemap
  62         raster_png_data,corner_coords = basemap_scatter_mercator(g[1],bbox,bmap,only_fire)
  63         # compute bounds
  64         bounds = (corner_coords[0][0],corner_coords[1][0],corner_coords[0][1],corner_coords[2][1])
  65         # write PNG file
  66         with open(pngfile, 'w') as f:
  67             f.write(raster_png_data)
  68         print '> File %s saved.' % g[0]
  69     else:
  70         print '> File %s already created.' % g[0]
  71     # append dictionary information for the KML creation
  72     kmld.append(Dict({'name': g[0], 'png_file': pngfile, 'bounds': bbox, 'time': timestamp}))
  73 # create KML
  74 create_kml(kmld,'./doc.kml')
  75 # create KMZ with all the PNGs included
  76 os.system('zip -r %s doc.kml *_A*_*.png' % 'googlearth.kmz')
  77 print 'Created file googlearth.kmz'
  78
  79 print ''
  80 print '>> DONE <<'
  81 sys.exit()