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