| 1 | ## Python script to check model outputs within a given set of values from an external ASII file |
|---|
| 2 | # L. Fita, LMD. April 2015 |
|---|
| 3 | ## e.g. # checking_output.py -f histins.nc -v t2m -l 'limitsVariables.inf' |
|---|
| 4 | ## e.g. # python checking_output.py -f limit.nc -v SST,ALB -l 'limitsVariables.inf' |
|---|
| 5 | # NOTE: |
|---|
| 6 | # `nc_var_tools.py' set of python scripts is required |
|---|
| 7 | # can be found in subversion `http://svn.lmd.jussieu.fr/LMDZ_WRF/trunk/tools/nc_var_tools.py' |
|---|
| 8 | import numpy as np |
|---|
| 9 | from netCDF4 import Dataset as NetCDFFile |
|---|
| 10 | import os |
|---|
| 11 | import re |
|---|
| 12 | from optparse import OptionParser |
|---|
| 13 | import nc_var_tools as ncvar |
|---|
| 14 | |
|---|
| 15 | main = 'checking_output.py' |
|---|
| 16 | version = 0.1 |
|---|
| 17 | author = 'L. Fita' |
|---|
| 18 | institution = 'Laboratoire de Meteorologie Dynamique' |
|---|
| 19 | city = ' Paris' |
|---|
| 20 | country = 'France' |
|---|
| 21 | |
|---|
| 22 | errormsg = 'ERROR -- error -- ERROR -- error' |
|---|
| 23 | warnmsg = 'WARNING -- warning -- WARNING -- warning' |
|---|
| 24 | |
|---|
| 25 | fillValue = 1.e20 |
|---|
| 26 | |
|---|
| 27 | def compute_stddev(vals): |
|---|
| 28 | """ Function to compute standard deviation |
|---|
| 29 | vals= values |
|---|
| 30 | """ |
|---|
| 31 | fname = 'compute_stddev' |
|---|
| 32 | |
|---|
| 33 | mvals = np.mean(vals) |
|---|
| 34 | m2vals = np.mean(vals*vals) |
|---|
| 35 | |
|---|
| 36 | stddev = np.sqrt(m2vals - mvals*mvals) |
|---|
| 37 | |
|---|
| 38 | return stddev |
|---|
| 39 | |
|---|
| 40 | def checking(cn,valc,vern,lval,Vn): |
|---|
| 41 | """ Function to check a variable with |
|---|
| 42 | cn= checking name |
|---|
| 43 | valc= value to check |
|---|
| 44 | vern= verification value |
|---|
| 45 | lval= value to use for checking |
|---|
| 46 | Vn= variable name |
|---|
| 47 | """ |
|---|
| 48 | fname = 'checking' |
|---|
| 49 | |
|---|
| 50 | VerN = vern[0:1] |
|---|
| 51 | |
|---|
| 52 | good = True |
|---|
| 53 | if VerN == '<' and valc < lval: |
|---|
| 54 | print warnmsg |
|---|
| 55 | print ' ' + fname + " variable '" + vn + "' has a '" +cn+ "'",valc,vern,lval |
|---|
| 56 | good = False |
|---|
| 57 | elif VerN == '>' and valc > lval: |
|---|
| 58 | print warnmsg |
|---|
| 59 | print ' ' + fname + " variable '" + vn + "' has a '" +cn+ "'",valc,vern,lval |
|---|
| 60 | good = False |
|---|
| 61 | elif VerN == '%' and np.abs((valc-lval)*100./lval)>np.float(vern[1:len(vern)+1]): |
|---|
| 62 | print warnmsg |
|---|
| 63 | print ' ' + fname + " variable '" + vn + "' has a '" +cn+ "'", valc, '>', \ |
|---|
| 64 | VerN, np.float(vern[1:len(vern)+1]), 'than',lval |
|---|
| 65 | good = False |
|---|
| 66 | |
|---|
| 67 | return good |
|---|
| 68 | |
|---|
| 69 | |
|---|
| 70 | ####### ###### ##### #### ### ## # |
|---|
| 71 | |
|---|
| 72 | parser = OptionParser() |
|---|
| 73 | parser.add_option("-f", "--file", dest="ncfile", |
|---|
| 74 | help="file to check", metavar="FILE") |
|---|
| 75 | parser.add_option("-v", "--variable", dest="varns", |
|---|
| 76 | help="var to check (',' list of variables; 'all', for all variables)", metavar="VALUE") |
|---|
| 77 | parser.add_option("-l", "--ASCIIvarlimits", dest="lfile", |
|---|
| 78 | help="ASCII file with the limits for the variables (long explanation with -l h)", metavar="FILE") |
|---|
| 79 | |
|---|
| 80 | (opts, args) = parser.parse_args() |
|---|
| 81 | |
|---|
| 82 | ####### ####### |
|---|
| 83 | ## MAIN |
|---|
| 84 | ####### |
|---|
| 85 | |
|---|
| 86 | # Available statistics |
|---|
| 87 | Availablestats = ['min', 'max', 'mean', 'stddev'] |
|---|
| 88 | LAvailablestats = {'min':'minimum', 'max':'maximum', 'mean':'mean', \ |
|---|
| 89 | 'stddev':'standard deviation'} |
|---|
| 90 | |
|---|
| 91 | Availableverifs = ['<', '>', '%'] |
|---|
| 92 | SAvailableverifs = {'<':'lt', '>':'gt', '%':'percen'} |
|---|
| 93 | LAvailableverifs = {'<':'smaller', '>':'bigger', '%':'percentage'} |
|---|
| 94 | |
|---|
| 95 | ####### |
|---|
| 96 | ofile = 'wrongvariables.nc' |
|---|
| 97 | |
|---|
| 98 | if not os.path.isfile(opts.ncfile): |
|---|
| 99 | print errormsg |
|---|
| 100 | print ' ' + main + ": file '" + opts.ncfile + "' does not exist !!" |
|---|
| 101 | quit(-1) |
|---|
| 102 | |
|---|
| 103 | if opts.varns is None: |
|---|
| 104 | print warnmsg |
|---|
| 105 | print ' ' + main + ": no variables given!!" |
|---|
| 106 | print ' checking all' |
|---|
| 107 | varns = 'all' |
|---|
| 108 | else: |
|---|
| 109 | varns = opts.varns |
|---|
| 110 | |
|---|
| 111 | if opts.lfile is None: |
|---|
| 112 | print warnmsg |
|---|
| 113 | print ' ' + main + ": no ASCII file with limits for the variables is given!!" |
|---|
| 114 | quit(-1) |
|---|
| 115 | elif opts.lfile == 'h': |
|---|
| 116 | print 'HELP: ASCII limits of the variables_____________' |
|---|
| 117 | print 'First line: varname [stat1] [stat2] .... [statn]' |
|---|
| 118 | print 'Second line: - [verif1] [verif2] ... [verifn]' |
|---|
| 119 | print 'afterwards: [varname] [value1] [value2] ... [valuen]' |
|---|
| 120 | print ' [stat]: statistics/parameter to compute. Has to exist wihin the script' |
|---|
| 121 | print ' directly with an if or as function compute[stat]. Available ones are:' |
|---|
| 122 | print ' ',Availablestats |
|---|
| 123 | print ' [verif]: verification to make with the stat' |
|---|
| 124 | print ' <: Error if variable has values smaller than the limit' |
|---|
| 125 | print ' >: Error if variable has values bigger than the limit' |
|---|
| 126 | print ' %N: Error if variable has values +/- N % than the limit' |
|---|
| 127 | quit() |
|---|
| 128 | |
|---|
| 129 | ncobj = NetCDFFile(opts.ncfile, 'r') |
|---|
| 130 | filevars = ncobj.variables.keys() |
|---|
| 131 | |
|---|
| 132 | # Getting variables names |
|---|
| 133 | if varns == 'all': |
|---|
| 134 | vns = filevars |
|---|
| 135 | else: |
|---|
| 136 | vns = opts.varns.split(',') |
|---|
| 137 | |
|---|
| 138 | Nvars = len(vns) |
|---|
| 139 | print "From file '" + opts.ncfile + "' checking",Nvars,'variables...' |
|---|
| 140 | |
|---|
| 141 | # Getting ASCII file and its values |
|---|
| 142 | ## |
|---|
| 143 | # Assuming first line with something like varname min max mean stddev |
|---|
| 144 | # Continuted with a line for each variable and value like: |
|---|
| 145 | # t2m -70. 50. 20. 40. |
|---|
| 146 | # (...) |
|---|
| 147 | # It is required than this script some function called compute_[value] must be placed |
|---|
| 148 | # '#' for comment |
|---|
| 149 | olval = open(opts.lfile, 'r') |
|---|
| 150 | ili = 0 |
|---|
| 151 | statsverif = {} |
|---|
| 152 | limitvals = {} |
|---|
| 153 | |
|---|
| 154 | for line in olval: |
|---|
| 155 | linevals = ncvar.reduce_spaces(line) |
|---|
| 156 | if len(linevals) != 0 and linevals[0][0:0] != '#': |
|---|
| 157 | # print ili,linevals |
|---|
| 158 | if ili == 0: |
|---|
| 159 | NVals = len(linevals) - 1 |
|---|
| 160 | stats = linevals[1:NVals + 1] |
|---|
| 161 | elif ili == 1: |
|---|
| 162 | for ist in range(NVals): |
|---|
| 163 | statsverif[stats[ist]] = linevals[1+ist] |
|---|
| 164 | else: |
|---|
| 165 | limitvals[linevals[0]] = np.zeros((NVals), dtype=np.float) |
|---|
| 166 | for ist in range(NVals): |
|---|
| 167 | limitvals[linevals[0]][ist] = np.float(linevals[1+ist]) |
|---|
| 168 | |
|---|
| 169 | ili = ili + 1 |
|---|
| 170 | |
|---|
| 171 | # Checking stats to compute |
|---|
| 172 | for st in stats: |
|---|
| 173 | if not ncvar.searchInlist(Availablestats,st): |
|---|
| 174 | print errormsg |
|---|
| 175 | print ' ' + main + ": statistics/value '" + st + "' not ready!!" |
|---|
| 176 | print ' available parameters:', Availablestats |
|---|
| 177 | quit(-1) |
|---|
| 178 | |
|---|
| 179 | # Checking verifiaction to make |
|---|
| 180 | for st in stats: |
|---|
| 181 | verin = statsverif[st] |
|---|
| 182 | if not ncvar.searchInlist(Availableverifs,verin[0:1]): |
|---|
| 183 | print errormsg |
|---|
| 184 | print ' ' + main + ": verification '" + verin + "' not ready!!" |
|---|
| 185 | print ' available verifications:', Availableverifs |
|---|
| 186 | quit(-1) |
|---|
| 187 | |
|---|
| 188 | print 'Checking',NVals,'parameters:',stats,'verifyed as:',statsverif.values() |
|---|
| 189 | |
|---|
| 190 | # Creating output file |
|---|
| 191 | oobj = NetCDFFile(ofile, 'w') |
|---|
| 192 | |
|---|
| 193 | # Dimensinos |
|---|
| 194 | newdim = oobj.createDimension('stats', NVals) |
|---|
| 195 | newdim = oobj.createDimension('verif', NVals) |
|---|
| 196 | newdim = oobj.createDimension('checks', 2) |
|---|
| 197 | newdim = oobj.createDimension('Lstring', 50) |
|---|
| 198 | |
|---|
| 199 | # Variables with Dimensions values |
|---|
| 200 | newvar = oobj.createVariable('stats', 'c', ('stats', 'Lstring')) |
|---|
| 201 | newattr = ncvar.basicvardef(newvar, 'stats', 'statistics/parameters used to check', \ |
|---|
| 202 | '-') |
|---|
| 203 | newvals = ncvar.writing_str_nc(newvar, stats, 50) |
|---|
| 204 | |
|---|
| 205 | newvar = oobj.createVariable('verifs', 'c', ('stats', 'Lstring')) |
|---|
| 206 | newattr = ncvar.basicvardef(newvar, 'verifs', 'verifications used to check', '-') |
|---|
| 207 | newvals = ncvar.writing_str_nc(newvar, stats, 50) |
|---|
| 208 | |
|---|
| 209 | newvar = oobj.createVariable('checks', 'c', ('checks', 'Lstring')) |
|---|
| 210 | newattr = ncvar.basicvardef(newvar, 'checks', 'values used to check', '-') |
|---|
| 211 | newvals = ncvar.writing_str_nc(newvar, ['computed', 'limits'], 50) |
|---|
| 212 | |
|---|
| 213 | # Global attributes |
|---|
| 214 | newattr = ncvar.set_attribute(oobj, 'script', main) |
|---|
| 215 | newattr = ncvar.set_attribute(oobj, 'version', version) |
|---|
| 216 | newattr = ncvar.set_attribute(oobj, 'author', author) |
|---|
| 217 | newattr = ncvar.set_attribute(oobj, 'institution', institution) |
|---|
| 218 | newattr = ncvar.set_attribute(oobj, 'city', city) |
|---|
| 219 | newattr = ncvar.set_attribute(oobj, 'country', country) |
|---|
| 220 | |
|---|
| 221 | newattr = ncvar.set_attribute(oobj, 'description', 'Variables from ' + opts.ncfile + \ |
|---|
| 222 | ' with a wrong value according to ' + opts.lfile) |
|---|
| 223 | |
|---|
| 224 | oobj.sync() |
|---|
| 225 | |
|---|
| 226 | # Getting netCDF file and its values |
|---|
| 227 | ## |
|---|
| 228 | valsstats = np.zeros((NVals,Nvars), dtype = np.float) |
|---|
| 229 | wrongVars = {} |
|---|
| 230 | iv = 0 |
|---|
| 231 | for iv in range(Nvars): |
|---|
| 232 | vn = vns[iv] |
|---|
| 233 | print vn + '...' |
|---|
| 234 | if not ncvar.searchInlist(filevars, vn): |
|---|
| 235 | print errormsg |
|---|
| 236 | print ' ' + main + ": file '" + opts.ncfile + "' does not have variable '" +\ |
|---|
| 237 | vn + "' !!" |
|---|
| 238 | quit(-1) |
|---|
| 239 | |
|---|
| 240 | vals = ncobj.variables[vn][:] |
|---|
| 241 | wrongstats = [] |
|---|
| 242 | for ist in range(NVals): |
|---|
| 243 | stn = stats[ist] |
|---|
| 244 | Lstn = LAvailablestats[stn] |
|---|
| 245 | statpos = stats.index(stn) |
|---|
| 246 | limitval = limitvals[vn][statpos] |
|---|
| 247 | verif = statsverif[stn] |
|---|
| 248 | # print ' ', ist, stn, ':',limitval,verif |
|---|
| 249 | if stn == 'max': |
|---|
| 250 | valsstats[ist,iv] = np.max(vals) |
|---|
| 251 | elif stn == 'mean': |
|---|
| 252 | valsstats[ist,iv] = np.mean(vals) |
|---|
| 253 | elif stn == 'min': |
|---|
| 254 | valsstats[ist,iv] = np.min(vals) |
|---|
| 255 | elif stn == 'stddev': |
|---|
| 256 | valsstats[ist,iv] = compute_stddev(vals) |
|---|
| 257 | |
|---|
| 258 | # checks |
|---|
| 259 | if not checking(stn,valsstats[ist,iv],verif,limitval,vn): |
|---|
| 260 | vval = valsstats[ist,iv] |
|---|
| 261 | Sverif = SAvailableverifs[verif[0:1]] |
|---|
| 262 | Lverif = LAvailableverifs[verif[0:1]] |
|---|
| 263 | if verif[0:1] != '<' and verif[0:1] != '>': |
|---|
| 264 | Vverif = np.float(verif[1:len(verif)]) |
|---|
| 265 | VverifS = str(Vverif) |
|---|
| 266 | else: |
|---|
| 267 | VverifS = '' |
|---|
| 268 | |
|---|
| 269 | wrongstats.append(stn) |
|---|
| 270 | if not ncvar.searchInlist(oobj.variables.keys(), vn): |
|---|
| 271 | oobj.close() |
|---|
| 272 | ncvar.fvaradd(opts.ncfile + ':' + vn, ofile) |
|---|
| 273 | |
|---|
| 274 | oobj = NetCDFFile(ofile ,'a') |
|---|
| 275 | vobj = oobj.variables[vn] |
|---|
| 276 | if ncvar.searchInlist(vobj.ncattrs(),('units')): |
|---|
| 277 | vunits = vobj.getncattr('units') |
|---|
| 278 | else: |
|---|
| 279 | vunits = '-' |
|---|
| 280 | |
|---|
| 281 | dims = vobj.dimensions |
|---|
| 282 | # 1/0 Matrix with wrong values |
|---|
| 283 | if stn == 'max': |
|---|
| 284 | wrongvals = np.where(vals > limitval, 1, 0) |
|---|
| 285 | newvar = oobj.createVariable(vn+'Wrong'+stn+Sverif, 'i', dims) |
|---|
| 286 | newattr = ncvar.basicvardef(newvar, vn+'_wrong_'+stn+'_'+Sverif, \ |
|---|
| 287 | 'wrong ' + vn + ' ' + str(vval) + ' ' + Lstn + ' ' + VverifS + ' '+\ |
|---|
| 288 | Lverif + ' ' + str(limitval), vunits) |
|---|
| 289 | newvar[:] = wrongvals |
|---|
| 290 | elif stn == 'mean': |
|---|
| 291 | newattr = ncvar.set_attribute(vobj, 'wrong_'+stn+'_'+Sverif+VverifS, \ |
|---|
| 292 | 'wrong ' + vn + ' ' + str(vval) + ' ' + Lstn + ' ' + VverifS + ' '+\ |
|---|
| 293 | Lverif + ' ' + str(limitval)) |
|---|
| 294 | elif stn == 'min': |
|---|
| 295 | wrongvals = np.where(vals < limitval, 1, 0) |
|---|
| 296 | newvar = oobj.createVariable(vn+'Wrong'+stn+Sverif, 'i', dims) |
|---|
| 297 | newattr = ncvar.basicvardef(newvar, vn+'_wrong_'+stn+'_'+Sverif, \ |
|---|
| 298 | 'wrong ' + vn + ' ' + str(vval) + ' ' + Lstn + ' ' + VverifS + ' '+\ |
|---|
| 299 | Lverif + ' ' + str(limitval), vunits) |
|---|
| 300 | newvar[:] = wrongvals |
|---|
| 301 | elif stn == 'stddev': |
|---|
| 302 | newattr = ncvar.set_attribute(vobj, 'wrong_'+stn+'_'+Sverif+VverifS, \ |
|---|
| 303 | 'wrong ' + vn + ' ' + str(vval) + ' ' + Lstn + ' ' + VverifS + ' '+\ |
|---|
| 304 | Lverif + ' ' + str(limitval)) |
|---|
| 305 | |
|---|
| 306 | if len(wrongstats) != 0: wrongVars[vn] = wrongstats |
|---|
| 307 | newattr = ncvar.set_attribute(vobj, 'wrong', \ |
|---|
| 308 | ncvar.numVector_String(wrongVars[vn], ' ')) |
|---|
| 309 | |
|---|
| 310 | newvar = oobj.createVariable(vn + 'checks', 'f4', ('checks', 'stats')) |
|---|
| 311 | newattr = ncvar.basicvardef(newvar, vn + '_checks', vn + \ |
|---|
| 312 | ' computed/checked', vunits) |
|---|
| 313 | for ist in range(NVals): |
|---|
| 314 | newvar[:,ist] = [valsstats[ist,iv], limitvals[vn][ist]] |
|---|
| 315 | |
|---|
| 316 | oobj.sync() |
|---|
| 317 | oobj.close() |
|---|
| 318 | |
|---|
| 319 | # Results |
|---|
| 320 | ## |
|---|
| 321 | print ' Variables with Wrong results _______________________________________________' |
|---|
| 322 | print wrongVars |
|---|
| 323 | |
|---|
| 324 | print "succesfull writting of checking output file '" + ofile + "' !!" |
|---|
![(please configure the [header_logo] section in trac.ini)](/LMDZ_WRF/chrome/site/your_project_logo.png){kind=logo}
