Changeset 503 in lmdz_wrf for trunk/tools

Timestamp:

Jun 15, 2015, 6:55:30 PM (10 years ago)

Author:

lfita

Message:

Working version with the new statistics

File:

• trunk/tools/validation_sim.py (modified) (10 diffs)

trunk/tools/validation_sim.py

@@ -1243 +1243 @@
 statsn = ['minimum', 'maximum', 'mean', 'mean2', 'standard deviation']
 gstatsn = ['bias', 'simobs_mean', 'sim_obsmin', 'sim_obsmax', 'sim_obsmean', 'mae',  \
-  'rmse', 'r_pearsoncorr', 'p_pearsoncorr']
+  'rmse', 'r_pearsoncorr', 'p_pearsoncorr', 'deviation_of_residuals_SDR',            \
+  'indef_of_efficiency_IOE', 'index_of_agreement_IOA', 'fractional_mean_bias_FMB']
 ostatsn = ['number of points', 'minimum', 'maximum', 'mean', 'mean2',                \
   'standard deviation']
@@ -1600 +1601 @@
 newdim = onewnc.createDimension('xaround',Ngrid*2+1)
 newdim = onewnc.createDimension('yaround',Ngrid*2+1)
-newdim = onewnc.createDimension('gstats',9)
+newdim = onewnc.createDimension('gstats',13)
 newdim = onewnc.createDimension('stats',5)
 newdim = onewnc.createDimension('tstats',6)
-newdim = onewnc.createDimension('Nstsim, 2)
+newdim = onewnc.createDimension('Nstsim', 2)
 
 # Variable dimensions
@@ -1742 +1743 @@
 # statisics obs
 # Masking 'nan'
+    Eobsmask0 = np.where(Earrayvals[:,1] != Earrayvals[:,1], fillValueF, 
+      Earrayvals[:,1]) 
     obsmask0 = np.where(arrayvals[:,1] != arrayvals[:,1], fillValueF, arrayvals[:,1]) 
+
+    Eobsmask = ma.masked_equal(Eobsmask0, fillValueF)
+    Eobsmask2 = Eobsmask*Eobsmask
     obsmask = ma.masked_equal(obsmask0, fillValueF)
     obsmask2 = obsmask*obsmask
 
-    obsstats = np.zeros((5), dtype=np.float)
-    obsstats[0] = obsmask.min()
-    obsstats[1] = obsmask.max()
-    obsstats[2] = obsmask.mean()
-    obsstats[3] = obsmask2.mean()
-    obsstats[4] = np.sqrt(obsstats[3] - obsstats[2]*obsstats[2])
+    obsstats = np.zeros((2,5), dtype=np.float)
+    obsstats[0,0] = obsmask.min()
+    obsstats[0,1] = obsmask.max()
+    obsstats[0,2] = obsmask.mean()
+    obsstats[0,3] = obsmask2.mean()
+    obsstats[0,4] = np.sqrt(obsstats[0,3] - obsstats[0,2]*obsstats[0,2])
+
+    obsstats[1,0] = obsmask.min()
+    obsstats[1,1] = obsmask.max()
+    obsstats[1,2] = obsmask.mean()
+    obsstats[1,3] = obsmask2.mean()
+    obsstats[1,4] = np.sqrt(obsstats[1,3] - obsstats[1,2]*obsstats[1,2])
  
 # Statistics sim-obs
@@ -1758 +1770 @@
     diffvals = np.zeros((Ncoindt), dtype=np.float)
 
-    Ediffvals = Earrayvals[:,0] - obsmask
+    Ediffvals = Earrayvals[:,0] - Eobsmask
     diffvals = arrayvals[:,0] - obsmask
 
@@ -1769 +1781 @@
     sumdiff2 = diff2vals.sum()
 
-    simobsstats[0,0] = simstats[0,0] - obsstats[0]
-    simobsstats[0,1] = np.mean(Earrayvals[:,0]*obsmask)
+    simobsstats[0,0] = simstats[0,0] - obsstats[0,0]
+    simobsstats[0,1] = np.mean(Earrayvals[:,0]*Eobsmask)
     simobsstats[0,2] = Ediffvals.min()
     simobsstats[0,3] = Ediffvals.max()
@@ -1777 +1789 @@
     simobsstats[0,6] = np.sqrt(Ediff2vals.mean())
     simobsstats[0,7], simobsstats[0,8] = sts.mstats.pearsonr(Earrayvals[:,0],        \
-      obsmask)
-
-    simobsstats[1,0] = simstats[1,0] - obsstats[0]
+      Eobsmask)
+
+    simobsstats[1,0] = simstats[1,0] - obsstats[1,0]
     simobsstats[1,1] = np.mean(arrayvals[:,0]*obsmask)
     simobsstats[1,2] = diffvals.min()
@@ -1789 +1801 @@
       obsmask)
 # From:
-#Willmott, C. J. 1981. On the validation of models. Physical Geography, 2, 184–194
-#Willmott, C. J. (1984). On the evaluation of model performance in physical geography. Spatial Statistics and Models, G. L. Gaile and C. J. Willmott, eds., 443-460
-#Willmott, C. J., S. G. Ackleson, R. E. Davis, J. J. Feddema, K. M. Klink, D. R. Legates, J. O'Donnell, and C. M. Rowe (1985), Statistics for the Evaluation and Comparison of Models, J. Geophys. Res., 90(C5), 8995-9005
-#Legates, D. R., and G. J. McCabe Jr. (1999), Evaluating the Use of "Goodness-of-Fit" Measures in Hydrologic and Hydroclimatic Model Validation, Water Resour. Res., 35(1), 233–241 
+#Willmott, C. J. 1981. 'On the validation of models. Physical Geography', 2, 184-194
+#Willmott, C. J. (1984). 'On the evaluation of model performance in physical 
+#  geography'. Spatial Statistics and Models, G. L. Gaile and C. J. Willmott, eds.,
+#  443-460
+#Willmott, C. J., S. G. Ackleson, R. E. Davis, J. J. Feddema, K. M. Klink, D. R. 
+#  Legates, J. O'Donnell, and C. M. Rowe (1985), 'Statistics for the Evaluation and
+#  Comparison of Models', J. Geophys. Res., 90(C5), 8995-9005
+#Legates, D. R., and G. J. McCabe Jr. (1999), 'Evaluating the Use of "Goodness-of-Fit"
+#   Measures in Hydrologic and Hydroclimatic Model Validation', Water Resour. Res.,
+#   35(1), 233-241 
 #
 # Deviation of residuals (SDR)
     simobsstats[0,9] = np.mean(np.sqrt(np.abs((Ediffvals-simobsstats[0,0])*(Ediffvals-\
-      obsstats[0]))))
+      obsstats[0,0]))))
     simobsstats[1,9] = np.mean(np.sqrt(np.abs((diffvals-simobsstats[1,0])*(diffvals- \
-      obsstats[0]))))
+      obsstats[1,0]))))
 # Index of Efficiency (IOE)
-    obsbias2series = (obsmask - obsstats[0])*(obsmask - obsstats[0])
+    Eobsbias2series = (Eobsmask - obsstats[0,0])*(Eobsmask - obsstats[0,0])
+    Esumobsbias2series = Eobsbias2series.sum()
+    obsbias2series = (obsmask - obsstats[1,0])*(obsmask - obsstats[1,0])
     sumobsbias2series = obsbias2series.sum()
 
-    simobsstats[1,10] = 1. - Esumdiff2/(sumobsbias2series)
+    simobsstats[0,10] = 1. - Esumdiff2/(Esumobsbias2series)
     simobsstats[1,10] = 1. - sumdiff2/(sumobsbias2series)
 # Index of Agreement (IOA)
-    Esimbias2series = Earrayvals[:,0] - obsstats[0]
-    simbias2series = arrayvals[:,0] - obsstats[0]
-
-    obsbias2series = obsmask - obsstats[0]
+    Esimbias2series = Earrayvals[:,0] - obsstats[0,0]
+    simbias2series = arrayvals[:,0] - obsstats[1,0]
+
+    Eobsbias2series = Eobsmask - obsstats[0,0]
+    obsbias2series = obsmask - obsstats[1,0]
 
     Eabssimbias2series = np.abs(Esimbias2series)
     abssimbias2series = np.abs(simbias2series)
-    absobsbias2series = obbias2series.abs()
-
-    Eabssimobsbias2series = (Eabssimbias2series+absobsbias2series)*(Eabssimbias2series +\
-      absobsbias2series)
+    Eabsobsbias2series = np.where(Eobsbias2series<0, -Eobsbias2series, Eobsbias2series)
+    absobsbias2series = np.where(obsbias2series<0, -obsbias2series, obsbias2series)
+
+    Eabssimobsbias2series=(Eabssimbias2series+Eabsobsbias2series)*(Eabssimbias2series+\
+      Eabsobsbias2series)
     abssimobsbias2series = (abssimbias2series+absobsbias2series)*(abssimbias2series +\
       absobsbias2series)
@@ -1823 +1845 @@
     simobsstats[1,11] = 1. - sumdiff2/(abssimobsbias2series.sum())
 # Fractional Mean Bias (FMB)
-    simobsstats[0,12] = (simstats[0,0]-obsstats[0])/(0.5*(simstats[0,0]+obsstats[0]))
-    simobsstats[1,12] = (simstats[1,0]-obsstats[0])/(0.5*(simstats[1,0]+obsstats[0]))
+    simobsstats[0,12]=(simstats[0,0]-obsstats[0,0])/(0.5*(simstats[0,0]+obsstats[0,0]))
+    simobsstats[1,12]=(simstats[1,0]-obsstats[1,0])/(0.5*(simstats[1,0]+obsstats[1,0]))
 
 # Statistics around sim values
@@ -1881 +1903 @@
 
 # sim Statistics
-    for ist in range(Nstsim):
     if not searchInlist(onewnc.variables,valvars[ivar][0] + 'stsim'):
         newvar = onewnc.createVariable(valvars[ivar][0] + 'stsim', 'f', ('Nstsim',   \
@@ -1891 +1912 @@
 # obs Statistics
     if not searchInlist(onewnc.variables,valvars[ivar][1] + 'stobs'):
-        newvar = onewnc.createVariable(valvars[ivar][1] + 'stobs', 'f', ('stats'),   \
-          fill_value=fillValueF)
+        newvar = onewnc.createVariable(valvars[ivar][1] + 'stobs', 'f', ('Nstsim',   \
+          'stats'), fill_value=fillValueF)
         descvar = 'observed statisitcs: ' + valvars[ivar][1]
         basicvardef(newvar, valvars[ivar][1] + 'stobs', descvar,                     \