source: trunk/LMDZ.PLUTO/util/script_figures/mapmeanwinds.py @ 3833

Last change on this file since 3833 was 3833, checked in by afalco, 17 hours ago

Pluto: updated plots scripts.
Fixed some issues with reading XIOS, etc.
Included display_netcdf.py tool from Mars PCM.
AF

  • Property svn:executable set to *
File size: 2.6 KB
Line 
1#! /usr/bin/env python
2from    netCDF4               import    Dataset
3from    numpy                 import    *
4import  numpy                 as        np
5import  matplotlib.pyplot     as        mpl
6from matplotlib.cm import get_cmap
7import pylab
8from matplotlib import ticker
9import matplotlib.colors as colors
10from mpl_toolkits.basemap import Basemap, shiftgrid
11from FV3_utils import *
12from input import *
13
14fa='sans-serif'
15hfont = {'fontname':'Arial'}
16mpl.rc('font',family=fa)
17mpl.rc('pdf',fonttype=42)
18font=30
19cc=['k']
20pal=get_cmap(name="rainbow")
21norm=colors.LogNorm()
22#lvls=np.logspace(-6,-4,21)
23norm=None #colors.LogNorm()
24lvls=np.linspace(37,57,21)
25
26altitude=1 # in km
27
28### Data
29filename1=name+'_A.nc'
30# name=name+'.nc'
31filename2=name+'.nc'
32print("Plot "+filename1)
33nc1=Dataset(filename1)
34nc2=Dataset(filename2)
35ts=nc2.variables["temperature"][:,0,:,:]
36# ts=nc2.variables["tsurf"][:,:,:]
37u=nc1.variables["u"][:,:,:,:]
38v=nc1.variables["v"][:,:,:,:]
39lat=getvar(nc1,"latitude")
40alt=getvar(nc1,"altitude")
41lon=getvar(nc1,"longitude")
42ps=nc2.variables["ps"][:,:]
43# ps=nc2.variables["phisinit"][:,:]/0.6169/1000.  # altitude km
44
45numalt=getind(altitude,alt)
46print('numalt =',numalt,'altitude=',alt[numalt])
47u=u[:,numalt,:,:]
48v=v[:,numalt,:,:]
49u=np.mean(u,axis=0)
50v=np.mean(v,axis=0)
51ts=np.mean(ts,axis=0)
52ps=np.mean(ps,axis=0)
53
54if lon[0]<0:
55    ts=switchlon(ts)
56    u=switchlon(u)
57    v=switchlon(v)
58    ps=switchlon(ps)
59    # topo=switchlon(topo)
60    lon=lon+180.
61
62### Figure
63fig=mpl.figure(figsize=(20, 10))
64# CF=mpl.contourf(lon, lat, ts,lvls,cmap=pal,norm=norm,extend='both')
65# cbar=mpl.colorbar(CF, shrink=1, ticks=lvls[::2]) #,format="%1.0f")
66# cbar.ax.set_title("[K]",y=1.04,fontsize=font)
67# for t in cbar.ax.get_yticklabels():
68#   t.set_fontsize(font)
69# mpl.title('Surface temperatures and winds',fontsize=font)
70
71#vect=lvls
72# CS=mpl.contour(lon,lat,topo,np.linspace(-6,6,21),colors='k',linewidths=0.5)
73#lab=mpl.clabel(CS, inline=1, fontsize=20, fmt='%1.2e',inline_spacing=1)
74#for l in lab:
75#  l.set_rotation(0)
76
77CF=mpl.contourf(lon, lat, ps,cmap=pal,norm=norm,extend='both')
78cbar=mpl.colorbar(CF, label="Pa")
79cbar.set_label("Pa", size=font)
80for t in cbar.ax.get_yticklabels():
81  t.set_fontsize(font)
82mpl.title('Surface pressures and winds',fontsize=font)
83
84getwinds(lon,lat,u,v,1,1,100,0.002,5)
85
86mpl.grid()
87mpl.ylabel(r'Latitude',labelpad=10,fontsize=font)
88mpl.xlabel('Longitude',labelpad=10, fontsize=font)
89pylab.ylim([-90,90])
90yticks=np.linspace(-90,90,13)
91pylab.xlim([0,360])
92xticks=np.linspace(0,360,7)
93mpl.yticks(yticks,fontsize=font)
94mpl.xticks(xticks,fontsize=font)
95mpl.tight_layout()
96mpl.savefig('mapwinds1km.png',bbox_inches='tight',dpi=70)
97mpl.show()
98
Note: See TracBrowser for help on using the repository browser.