source: BOL/Multi_atlas/atlas/atlas_SH.py @ 4284

Last change on this file since 4284 was 3743, checked in by idelkadi, 4 years ago

Repository under svn of a first version of Multiatlas diagnostics for LMDZ. This version is adapted to be able to run a LMDZ multiatlas on an individual account on the ciclad machine of the IPSL. In this version, the parts to be modified are identified so as to subsequently adapt it to other machines.
(cleaning)

File size: 11.1 KB
Line 
1# -*- coding: iso-8859-1 -*-
2# Created : S.Sénési - nov 2015
3
4#
5desc="\nCreation d'un atlas pour une simu, une grille et une liste de variables et de saisons \n"+\
6"  Exemples : \n"+\
7"  >>> python ./atlas.py -v tas,hfls -s NPv3.1ada_1982_1991\n"+\
8""
9# Avec CliMAF, cette etape est loin d'etre necessaire; on la réalise pour 'exposer' ces fichiers
10# dans une arborescence à laquelle sont habitues certains utilisateurs
11
12
13# Répertoire de base pour les entrées et les résultats
14#dir_default='/data/hourdin/LMDZ6/SE/ORIG'
15#dir_default='/prodigfs/fabric/LMDZ6/SE/ORIG'
16#dir_default='/prodigfs/ipslfs/dods/fabric/lmdz/SE/ORIG'
17
18# Gestion des options et arguments d'appel
19from optparse import OptionParser
20parser = OptionParser(desc) ; parser.set_usage("%%prog [-h]\n%s" % desc)
21parser.add_option("-i", "--input", help="repertoire des donnes d'entree(defaut : %s)"%dir_default, 
22                  action="store",default=dir_default)
23parser.add_option("-g", "--grid", help="nom de grille (default: VLR)", action="store",default='VLR')
24parser.add_option("-P", "--projection", help="Projection NH, SH ...", action="store",default='GLOB')
25parser.add_option("-p", "--season", help="saison a traiter " "(eg : JJA, DJF, YEAR, defaut=%YEAR)", 
26                  action="store", default='YEAR')
27parser.add_option("-s", "--simulation", help="simulation+annees a traiter (sim_YYY1_YYY2) - laisser vide pour lister",
28                  action="store",default=None)
29parser.add_option("-t", "--reference", help="simulation de reference (sim_YYY1_YYY2, default=OBS) ",
30                  action="store",default='OBS')
31parser.add_option("-v", "--variables", help="liste des variables (separees par des virgules)", action="store",default=None)
32parser.add_option("--root", help="Path to the root directory", action="store",default=None)
33parser.add_option("-f", "--force", help="force le recalcul de champs existants", 
34                  action="store_true",default=None)
35parser.add_option("-o", "--pdf", help="nom du pdf de sortie (default: atlas_<SIMU>_<SAISON>.pdf)", action="store")
36(opts, args) = parser.parse_args()
37
38#---------------------------------------------------------------------------------
39import math
40from climaf.api import *
41from climaf.html import * 
42# La description de l'organisation des données SE et des alias et rescalings
43# est partagée dans une micro-librairie :
44from lmdz_SE import * # svsg, all_SE_simulations
45from plot_params import plot_params
46#---------------------------------------------------------------------------------
47#
48def apply_scale_offset(dat,scale,offset):
49    return ccdo(ccdo(dat,operator='mulc,'+str(float(scale))),operator='addc,'+str(float(offset)))
50#
51#craz()
52if opts.simulation is None:
53    print "Available simulations at %s are : "%opts.input,
54    for s in all_SE_simulations() : print s,
55    exit(0)
56#
57lvars=opts.variables
58if lvars is not None : lvars=lvars.split(',')
59else : lvars=variables_list
60#
61# Preparons une commande pour assembler les sorties Pdf
62if opts.pdf : pdffile=opts.pdf
63else: pdffile="atlas_"+opts.simulation+"_"+opts.season+".pdf"
64pdfargs=["pdfjam","--landscape","-o ",pdffile]
65#
66# Initialisation de l'index html
67index= header("LMDZ Atlas for "+opts.simulation+ " versus "+opts.reference+" ("+opts.season+")") 
68index += cell('PDF',pdffile)
69index += section("2d vars", level=4)
70index += open_table()
71#
72# Titres de colonnes
73ref=opts.reference ; 
74if (ref == 'OBS' ) : text_diff='bias'
75else:                text_diff='diff'
76index+=open_line('VARIABLE')+cell('bias')+cell('rmse')+cell('mean')+cell(ref)+cell(text_diff)+\
77        cell('zonal')+cell('all')+cell('pdf')+close_line()
78#
79# -- Declare the script ml2pl for vertical interpolation
80cscript("ml2pl", "/home/jservon/Evaluation/CliMAF/Atlas_LMDz/ml2pl.sh -p ${var_2} -v ${var_1} ${in_1} ${out} ${in_2}",
81    commuteWithTimeConcatenation=True, commuteWithSpaceConcatenation=True)
82# -- Vertical levels for the vertical interpolation
83fixed_fields("ml2pl",("press_levels.txt","/home/jservon/Evaluation/CliMAF/press_levels.txt"))
84#
85for variable  in lvars :
86    # Get the model and the reference
87    if opts.root:
88       simu=svsg(opts.simulation,variable,opts.season,opts.grid, root=opts.root)
89    else:
90       simu=svsg(opts.simulation,variable,opts.season,opts.grid)
91    print 'variable = ',variable
92    reff=svsg(opts.reference,variable,opts.season,opts.grid)
93    #
94    # If the variable is a 3D field:
95    #  - interpolate the variable on the standard pressure levels with ml2pl (L. Guez)
96    #  - Compute the difference model-ref with diff_zonmean (computes the zonal mean lat/pressure fields,
97    #    interpolates the model on the ref, both vertically and horizontally, and returns the difference)
98    if is3d(variable) :
99       simu_pres = svsg(opts.simulation,'pres',opts.season,opts.grid)
100       simu = ml2pl(simu,simu_pres)
101       simu = zonmean(simu)
102       reff = zonmean(reff)
103       diff = diff_zonmean(simu,reff)
104    else:
105        if (opts.grid == '' ) : reff=regrid(reff,simu)
106        diff=minus(simu,reff)
107
108    pparams = plot_params(variable,'full_field')
109    vertical_interval = 'trYMaxF=1000|trYMinF=1'
110    stringFontHeight=0.018
111    if is3d(variable):
112        pparams.update({'options':vertical_interval})
113        stringFontHeight=0.03
114    # Map for simulation
115    simu_fig=plot(simu,title="",
116                  proj=opts.projection,
117                  gsnLeftString=variable,
118                  gsnCenterString=opts.simulation,
119                  gsnRightString=opts.season,
120                  gsnStringFontHeightF=stringFontHeight,
121                  mpCenterLonF=0,
122                  **pparams)
123    simu_avg=cvalue(space_average(simu))
124    #
125    # Map for reference
126    ref_fig=plot(reff,title="",
127                  proj=opts.projection,
128                 gsnLeftString=variable,
129                 gsnCenterString=ref,
130                 gsnRightString=opts.season,
131                 gsnStringFontHeightF=stringFontHeight,
132                 mpCenterLonF=0,
133                 **pparams)
134    ref_avg=cvalue(space_average(reff))
135    #
136    # Bias (or difference between simulations) map
137    if (ref == 'OBS' ) : p=plot_params(variable,'bias')
138    else:                p=plot_params(variable,'model_model')
139    tmp_aux_params = plot_params(variable,'full_field')
140    scale = 1.0 ; offset = 0.0
141    if 'offset' in tmp_aux_params or 'scale' in tmp_aux_params:
142       if 'offset' in tmp_aux_params:
143          offset = tmp_aux_params['offset']
144       else:
145          offset=0.0
146       if 'scale' in tmp_aux_params:
147          scale = tmp_aux_params['scale']
148       else:
149          scale=1.0
150       wreff = apply_scale_offset(reff,scale,offset)
151       wsimu = apply_scale_offset(simu,scale,offset)
152    else:
153       wreff = reff
154       wsimu = simu
155    #
156    if is3d(variable):
157        p.update({'options':vertical_interval})
158    if variable in ['ua','va','ta','hus']:
159        tmp_levs = tmp_aux_params['colors']
160        p.update({'contours':tmp_levs})
161        diff_fig=plot(diff,wreff,title="", format='png', mpCenterLonF=0,
162                  gsnLeftString=variable,
163                  gsnCenterString=opts.simulation+' - '+ref,
164                  gsnRightString=opts.season,
165                  gsnStringFontHeightF=stringFontHeight,
166                  aux_options='cnLineThicknessF=2|cnLineLabelsOn=True', **p)
167    else:
168        p.update({'contours':1})
169        diff_fig=plot(diff,title="", format='png', mpCenterLonF=0,
170                  proj=opts.projection,
171                  gsnLeftString=variable,
172                  gsnCenterString=opts.simulation+' - '+ref,
173                  gsnRightString=opts.season,
174                  gsnStringFontHeightF=stringFontHeight,
175                  **p)
176
177    #
178    # Bias mean value, and RMSD/RMSE
179    diff_avg=cvalue(space_average(diff))
180    rmsd=math.sqrt(cvalue(space_average(ccdo(diff,operator='-b F64 sqr'))))
181    #
182    # Zonal means
183    if not is3d(variable):
184        # -- apply a mask corresponding to the reference
185        mask = div(reff,reff)
186        msimu = mul(wsimu,mask)
187        # -- Compute the zonal mean
188        zmean=ccdo(msimu, operator='zonmean')
189        ref_zmean=ccdo(wreff, operator='zonmean')
190        #
191        sim=opts.simulation
192        #if variable in ['zg500']:
193        #   ref_zmean = ccdo(ref_zmean,operator='-b F32 mulc,1')
194        #   zmean = ccdo(zmean,operator='-b F32 mulc,1')
195        zmean_fig=curves(cens([sim,ref],zmean,ref_zmean),
196                         title="",
197                         lgcols=3,
198                         options=#'tiYAxisString=""|'+\
199                                 #'+\'+\
200                                 'tmYROn=True|'+\
201                                 'tmYRBorderOn=True|'+\
202                                 'tmYLOn=False|'+\
203                                 'tmYUseRight=True|'+\
204                                 'vpXF=0|'+\
205                                 'vpWidthF=0.66|'+\
206                                 'vpHeightF=0.33|'+\
207                                 'tmYRLabelsOn=True|'+\
208                                 'tmXBLabelFontHeightF=0.018|'+\
209                                 'tmYLLabelFontHeightF=0.016|'+\
210                                 'lgLabelFontHeightF=0.018|'+\
211                                 #'pmLegendSide=Bottom|'+\
212                                 'pmLegendOrthogonalPosF=-0.32|'+\
213                                 'pmLegendParallelPosF=1.0|'+\
214                                 'tmXMajorGrid=True|'+\
215                                 'tmYMajorGrid=True|'+\
216                                 'tmXMajorGridLineDashPattern=2|'+\
217                                 'tmYMajorGridLineDashPattern=2|'+\
218                                 'xyLineThicknessF=8|'+\
219                                 'gsnLeftString='+variable+'|'+\
220                                 'gsnCenterString='+opts.simulation+' vs '+ref+'|'+\
221                                 'gsnRightString='+opts.season+'|'+\
222                                 'gsnStringFontHeightF='+str(stringFontHeight))
223    #    # Composite figure
224    if is3d(variable):
225        page=cpage([[simu_fig,ref_fig],[diff_fig,None]],orientation='landscape', page_trim=True, fig_trim=True)
226        pdf_page=cpage([[simu_fig,ref_fig],[diff_fig,None]],orientation='landscape',
227                   page_trim=True, fig_trim=True, format='pdf', title=variable+" "+opts.simulation)
228    else:
229        page=cpage([[simu_fig,ref_fig],[diff_fig,zmean_fig]],orientation='landscape', page_trim=True, fig_trim=True)
230        pdf_page=cpage([[simu_fig,ref_fig],[diff_fig,zmean_fig]],orientation='landscape', 
231                   page_trim=True, fig_trim=True, format='pdf', title=variable+" "+opts.simulation)
232    pdfargs.append(cfile(pdf_page))
233    #
234    thumbnail_size = 200
235    if is3d(variable):
236            index+=open_line(varlongname(variable)+' ('+variable+')')+\
237                    cell("%.2g"%diff_avg,cfile(diff_fig))+\
238                    cell("%.2g"%rmsd,cfile(diff_fig))+\
239                    cell(simu,cfile(simu_fig),thumbnail=thumbnail_size,hover=False)+\
240                    cell(ref,cfile(ref_fig),thumbnail=thumbnail_size,hover=False)+\
241                    cell(text_diff,cfile(diff_fig),thumbnail=thumbnail_size,hover=False)+\
242                    ' '+\
243                    cell('page',cfile(page),thumbnail=thumbnail_size,hover=False)+\
244                    cell('Pdf',cfile(pdf_page))
245            close_line()
246    else:
247            index+=open_line(varlongname(variable)+' ('+variable+')')+\
248                    cell("%.2g"%diff_avg,cfile(diff_fig))+\
249                    cell("%.2g"%rmsd,cfile(diff_fig))+\
250                    cell(simu,cfile(simu_fig),thumbnail=thumbnail_size,hover=False)+\
251                   cell(ref,cfile(ref_fig),thumbnail=thumbnail_size,hover=False)+\
252                    cell(text_diff,cfile(diff_fig),thumbnail=thumbnail_size,hover=False)+\
253                    cell('zonal mean',cfile(zmean_fig),thumbnail=thumbnail_size,hover=False)+\
254                    cell('page',cfile(page),thumbnail=thumbnail_size,hover=False)+\
255                    cell('Pdf',cfile(pdf_page))
256            close_line()
257#
258# Finalisons l'index html
259index += close_table()
260index += trailer()
261#out="index_example.html"
262out="index_example_"+opts.season+"_"+opts.simulation+".html"
263with open(out,"w") as filout : filout.write(index)
264#
265# Creation du Pdf multi-pages
266comm=subprocess.Popen(pdfargs, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
267#
268import os,os.path ; 
269# print("Attendez un bon peu : lancemement de firefox sur Ciclad....")
270# os.system("firefox file://"+os.path.abspath(os.path.curdir)+"/"+out+"&")
Note: See TracBrowser for help on using the repository browser.