Changeset 776 in lmdz_wrf


Ignore:
Timestamp:
May 25, 2016, 5:10:25 PM (9 years ago)
Author:
lfita
Message:

Adding pressure on `draw_vertical_levels'

Location:
trunk/tools
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/tools/drawing.py

    r774 r776  
    36683668    draw_topo_geogrid(ncfile, values, varn)
    36693669      ncfile= file to use
    3670       values= [zlog]:[dzlog]:[title]:[graphic_kind]:[legloc]
     3670      values= [zlogs]:[plogs]:[title]:[graphic_kind]:[legloc]
     3671        zlogs= zlog,dzlog
    36713672        zlog: to use logarithmic scale on the height axis ('true/false')
    36723673        dzlog: to use logarithmic scale on the difference of height between levels axis ('true/false')
     3674        plogs= plog,dplog
     3675        plog: to use logarithmic scale on the height axis ('true/false')
     3676        dplog: to use logarithmic scale on the difference of height between levels axis ('true/false')
    36733677        title: title of the graph ('!' for spaces)
    36743678        graphic_kind: kind of figure (jpg, pdf, png)
     
    36773681          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
    36783682          9: 'upper center', 10: 'center'      kfig= kind of figure
    3679       varn= name of the variable with the vertical levels
    3680         'WRFz': for WRF z-levels (computed as (PH + PHB)/g, from a PHB(0,i,j) = 0)
     3683      varn= [varnheight],[varnpres]
     3684        varnheight: name of the variable with the height of the vertical levels
     3685          'WRFz': for WRF z-levels (computed as (PH + PHB)/g, from a PHB(0,i,j) = 0)
     3686        varnpres: name of the variable with the pressure of the vertical levels ('None', for no pressure plot)
     3687          'WRFp': for WRF p-levels (computed as P + PB, from a PHB(0,i,j) = 0)
    36813688    """
    36823689    fname = 'draw_vertical_levels'
     
    36873694        quit()
    36883695
    3689     expectargs = '[zlog]:[dzlog]:[title]:[graphic_kind]:[legloc]'
     3696    expectargs = '[zlogs]:[plogs]:[title]:[graphic_kind]:[legloc]'
    36903697 
    36913698    drw.check_arguments(fname,values,expectargs,':')
    36923699
    3693     zlog = values.split(':')[0]
    3694     dzlog = values.split(':')[1]
     3700    zlog = values.split(':')[0].split(',')[0]
     3701    dzlog = values.split(':')[0].split(',')[1]
     3702    plog = values.split(':')[1].split(',')[0]
     3703    dplog = values.split(':')[1].split(',')[1]
    36953704    title = values.split(':')[2].replace('!',' ')
    36963705    kindfig = values.split(':')[3]
    36973706    legloc = int(values.split(':')[4])
    36983707
     3708    if varn.find(',') == -1:
     3709        varnheight = varn
     3710        varnpres = None
     3711        pvals = None
     3712        print warnmsg
     3713        print '  ' + fname + ': assuming no pressure variable!!'
     3714    else:
     3715        varnheight = varn.split(',')[0]
     3716        varnpres = varn.split(',')[1]
     3717        if varnpres == 'None':
     3718            varnpres = None
     3719            pvals = None
     3720
    36993721    if not os.path.isfile(ncfile):
    37003722        print errormsg
     
    37043726    objf = NetCDFFile(ncfile, 'r')
    37053727
    3706     if varn == 'WRFz':
     3728    if varnheight == 'WRFz':
    37073729        if not gen.searchInlist(objf.variables,'PH'):
    37083730            print errormsg
     
    37233745        zvals = geop[0, :, ijz0[0], ijz0[1]] / 9.8
    37243746    else:
    3725         if not gen.searchInlist(objf.variables, varn):
     3747        if not gen.searchInlist(objf.variables, varnheight):
    37263748            print errormsg
    3727             print '  ' + fname + ": file '" + ncfile + "' does not have variable '" +\
    3728               varn + "' !!"
     3749            print '  ' + fname + ": file '" + ncfile + "' does not have height " +   \
     3750              " variable '" + varnheight + "' !!"
    37293751            quit(-1)
    37303752        objvar = objf.variables[varn]
    37313753        if len(objvar.shape) == 4:
    37323754            print warnmsg
    3733             print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
    3734               "shape: dt,dz,dy,dx. Tacking first time-step"
     3755            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
     3756              "' with shape: dt,dz,dy,dx. Tacking first time-step"
    37353757
    37363758            ijz0 = gen.index_mat(objvar[0,0,], 0.)
     
    37383760        elif len(objvar.shape) == 3:
    37393761            print warnmsg
    3740             print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
    3741               "shape: dz,dy,dx"
     3762            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
     3763              "' with shape: dz,dy,dx"
    37423764
    37433765            ijz0 = gen.index_mat(objvar[0,], 0.)
     
    37463768        elif len(objvar.shape) == 2:
    37473769            print warnmsg
    3748             print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
    3749               "shape: dz,dyx"
     3770            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
     3771              "' with shape: dz,dyx"
    37503772
    37513773            ijz0 = gen.index_mat(objvar[0,], 0.)
     
    37543776            zvals = objvar[:]
    37553777
     3778# Pressure
     3779    if varnpres is not None:
     3780        if varnpres == 'WRFp':
     3781            if not gen.searchInlist(objf.variables,'P'):
     3782                print errormsg
     3783                print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
     3784                  "variable 'P' !!"
     3785                quit(-1)
     3786            if not gen.searchInlist(objf.variables,'PB'):
     3787                print errormsg
     3788                print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
     3789                  "variable 'PB' !!"
     3790                quit(-1)
     3791
     3792            objph = objf.variables['P']
     3793            objphb = objf.variables['PB']
     3794            pres = objph[:] + objphb[:]
     3795
     3796            pvals = pres[0, :, ijz0[0], ijz0[1]]
     3797        else:
     3798            if not gen.searchInlist(objf.variables, varnpres):
     3799                print errormsg
     3800                print '  ' + fname + ": file '" + ncfile + "' does not have pressure " + \
     3801                  " variable '" + varnpres + "' !!"
     3802                quit(-1)
     3803            objvar = objf.variables[varnpres]
     3804            if len(objvar.shape) == 4:
     3805                print warnmsg
     3806                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
     3807                  "' with shape: dt,dz,dy,dx. Tacking first time-step"
     3808   
     3809                pvals = objvar[0, :, ijz0[0], ijz0[1]]
     3810            elif len(objvar.shape) == 3:
     3811                print warnmsg
     3812                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
     3813                  "' with shape: dz,dy,dx"
     3814   
     3815                pvals = objvar[:, ijz0[0], ijz0[1]]
     3816           
     3817            elif len(objvar.shape) == 2:
     3818                print warnmsg
     3819                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
     3820                  "' with shape: dz,dyx"
     3821   
     3822                pvals = objvar[:, ijz0[0]]
     3823            else:
     3824                pvals = objvar[:]
     3825
     3826# Logarithmic axes
    37563827    if zlog == 'true':
    37573828        zlogv = True
     
    37633834        print "    must be either: 'true' or 'false'"
    37643835        quit(-1)
     3836
    37653837    if dzlog == 'true':
    37663838        dzlogv = True
     
    37733845        quit(-1)
    37743846
    3775     drw.plot_vertical_lev(zvals, zlogv, dzlogv, title, kindfig, legloc)
     3847    if pvals is not None:
     3848        if plog == 'true':
     3849            plogv = True
     3850        elif plog == 'false':
     3851            plogv = False
     3852        else:
     3853            print errormsg
     3854            print '  ' + fname + ": wrong value for plog: '" + plog + "' !!"
     3855            print "    must be either: 'true' or 'false'"
     3856            quit(-1)
     3857        if dplog == 'true':
     3858            dplogv = True
     3859        elif dplog == 'false':
     3860            dplogv = False
     3861        else:
     3862            print errormsg
     3863            print '  ' + fname + ": wrong value for dplog: '" + dplog + "' !!"
     3864            print "    must be either: 'true' or 'false'"
     3865            quit(-1)
     3866
     3867    drw.plot_vertical_lev(zvals, pvals, zlogv, dzlogv, plogv, dplogv, title, kindfig, legloc)
    37763868
    37773869    objf.close()
  • trunk/tools/drawing_tools.py

    r773 r776  
    14841484        elif u == 'C': lu='$^{\circ}C$'
    14851485        elif u == 'days': lu='$day$'
     1486        elif u == 'degrees_East': lu='$degrees\ East$'
    14861487        elif u == 'degrees_east': lu='$degrees\ East$'
    14871488        elif u == 'degree_east': lu='$degrees\ East$'
    14881489        elif u == 'degrees longitude': lu='$degrees\ East$'
    14891490        elif u == 'degrees latitude': lu='$degrees\ North$'
     1491        elif u == 'degrees_North': lu='$degrees\ North$'
    14901492        elif u == 'degrees_north': lu='$degrees\ North$'
    14911493        elif u == 'degree_north': lu='$degrees\ North$'
     
    65626564
    65636565
    6564 def plot_vertical_lev(vertz, zlog, dzlog, gtit, kfig, lloc):
     6566def plot_vertical_lev(vertz, vertp, zlog, dzlog, plog, dplog, gtit, kfig, lloc):
    65656567    """ plotting vertical levels distribution
    65666568    plot_vertical_lev(vertz, gtit, kfig, lloc)
    65676569      vertz= distribution of vertical heights [z]
     6570      vertp= distribution of vertical pressures [Pa]
    65686571      zlog: to use logarithmic scale on the height axis ('true/false')
    65696572      dzlog: to use logarithmic scale on the difference of height between levels axis ('true/false')
     6573      plog: to use logarithmic scale on the pressure axis ('true/false')
     6574      pzlog: to use logarithmic scale on the difference of pressure between levels axis ('true/false')
    65706575      gtit= title of the graph ('!' for spaces)
    65716576      kfig= kind of figure (jpg, pdf, png)
     
    65906595
    65916596    plt.rc('text', usetex=True)
    6592 
    6593     fig, ax1 = plt.subplots()
     6597# Height plot
     6598##
     6599    if vertp is not None:
     6600        print plt.subplots.__doc__
     6601        fig, (ax1, ax3) = plt.subplots(nrows=2, ncols=1, sharex=True)
     6602    else:
     6603        fig, ax1 = plt.subplots()
    65946604    ax2 = ax1.twinx()
    65956605
    65966606    plt.xlim(0,Nlev)
    6597 
     6607 
    65986608    l1 = ax1.plot(range(Nlev), vertz, 'r-x', label='height')
    65996609    l2 = ax2.plot(range(1,Nlev), dvertz, 'b-x', label='dheight')
    66006610
    6601     ax1.set_xlabel('level (\#)')
     6611#    ax1.set_xlabel('level (\#)')
    66026612    ax1.set_ylabel('height (m)', color='r')
    66036613    ax1.set_ylim(1,np.max(vertz))
     
    66056615    ax2.set_ylabel('difference between levels (m)', color='b')
    66066616
    6607     plt.title(gtit)
     6617    if vertp is not None: ax1.set_title('height')
    66086618    if zlog: ax1.set_yscale('log')
    66096619    if dzlog: ax2.set_yscale('log')
    66106620
     6621# Pressure plot
     6622##
     6623    if vertp is not None:
     6624
     6625        Nlev = len(vertp)
     6626        dvertp = vertp[0:Nlev-1] - vertp[1:Nlev]
     6627
     6628        ax4 = ax3.twinx()
     6629
     6630        l3 = ax3.plot(range(Nlev), vertp, 'r-o', label='pressure')
     6631        l4 = ax4.plot(range(1,Nlev), dvertp, 'b-o', label='dpressure')
     6632
     6633        ax3.set_xlim(0, Nlev)
     6634        ax3.set_xlabel('level (\#)')
     6635        ax3.set_ylabel('pressure (Pa)', color='r')
     6636        ax3.set_ylim(np.min(vertp), np.max(vertp))
     6637        ax3.grid()
     6638        ax4.set_ylabel('difference between levels (Pa)', color='b')
     6639
     6640        if plog: ax3.set_yscale('log')
     6641        if dplog: ax4.set_yscale('log')
     6642
     6643    if vertp is not None:
     6644        fig.suptitle(gtit)
     6645        ax3.set_title('pressure')
     6646    else:
     6647        plt.title(gtit)
     6648
    66116649    output_kind(kfig, figname, True)
    66126650
Note: See TracChangeset for help on using the changeset viewer.