source: BOL/Multi_atlas/atlas/atlas.py @ 3942

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