Handle dates with read_DYNAMICO and read_CHIMERE

metrics/half_plume.py

@@ -2,6 +2,8 @@ import numpy as np
 import netCDF4 as nc
 import matplotlib.pyplot as plt
 
+CHIMERE_days = ['04', '05', '06', '07', '08', '09', '10', '11']
+
 # flatten the conc arrays TAKE ONLY FINAL SIMULATION TIME, ONLY ONE SPECIES ON DYNAMICO SIDE
 
 def compute_pressure(ap, bp, p0, ps):
@@ -40,10 +42,11 @@ def max_volume_plot(variable_chim, variable_dyn, day_list):
     plt.savefig('metric_figures/max_concentrations_versus_time.png')
 
 
-def read_DYNAMICO(dynamico_nc, dynamico_hybrids, hour=427):
+def read_DYNAMICO(dynamico_nc, dynamico_hybrids, day=8, hour=22):
+    hour = 2*( (day-1)*24+hour ) # assume output every 30 min  
     g, p0 = 9.81, 1e5 # gravity, reference pressure for hybrid coordinate
 
-    print('Reading DYNAMICO data ...')
+    print('Reading DYNAMICO data from %s ...'%dynamico_nc)
 
     ds = nc.Dataset(dynamico_nc)
     print(ds)
@@ -91,7 +94,7 @@ def read_DYNAMICO(dynamico_nc, dynamico_hybrids, hour=427):
     return volume, SO2_density
     
 
-def read_CHIMERE(chimere_nc, hour=24):
+def read_CHIMERE(chimere_nc, day=8, hour=24):
     #	double airm(Time, bottom_top, south_north, west_east) ;
     #		airm:units = "molecule/cm**3" ;
     #		airm:long_name = "Air density" ;
@@ -103,9 +106,9 @@ def read_CHIMERE(chimere_nc, hour=24):
     # * we neglect horizontal variations of cell area
     # => all cells have the same air mass
     # * we work with molecule number rather than mass (constant factor = molar mass)
-    print('Reading CHIMERE data ...')
-
-    ds2 = nc.Dataset(chimere_nc)
+    filename = chimere_nc%CHIMERE_days[day-1]
+    print('Reading CHIMERE data from %s ...'%filename)
+    ds2 = nc.Dataset(filename)
     GTRC1 = ds2['GTRC1'][hour,:,:,:].flatten() # SO2 molecular mixing ratio (ppb)
     air_density = ds2['airm'][hour,:,:,:].flatten()   # air molecule number per unit volume
     SO2_density = GTRC1 * air_density                # SO2 molecule number per unit volume
@@ -160,18 +163,19 @@ def bisection_method(half_mass, cum_mass):
 
 # --------------------------------------------------------------
 
-def main():
-    max_volume_plot('GTRC1', 'q', ['04', '05', '06', '07', '08', '09', '10', '11'])
+import argparse
 
+def main():
     dynamico_nc = '/data/PLT8/pub/smailler/PUY60-dynamico/dynamico.nc'
     dynamico_hybrids = '/data/PLT8/pub/smailler/PUY60-dynamico/apbp.nc'
     data2 = read_DYNAMICO(dynamico_nc, dynamico_hybrids)
     volume_half_mass_DYNAMICO= volume_half_mass(*data2)
 
-    chimere_nc = '/data/PLT8/pub/smailler/CHIMOUT-PUY60-tier2/out.2011061100_24_PUY60-tier2.nc'
+    chimere_nc = '/data/PLT8/pub/smailler/CHIMOUT-PUY60-tier2/out.201106%s00_24_PUY60-tier2.nc'
     data = read_CHIMERE(chimere_nc)
     volume_half_mass_CHIMERE = volume_half_mass(*data)
 
+    max_volume_plot('GTRC1', 'q', CHIMERE_days)
 
     #plots(metrics_CHIMERE, metrics_DYNAMICO)