source: BOL/Multi_atlas/atlas/zmean_multi.py @ 3902

Last change on this file since 3902 was 3684, 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.
This version is still under development.

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