Extract efficiently profiles from a FA file (whose fields are 2D horizontal)¶
[1]:
%matplotlib inline
import epygram
epygram.init_env()
[2]:
r = epygram.resources.meta_resource('../inputs/ICMSHAROM+0001', 'r', 'CL') # CL = CombineLevels
[3]:
# in a "CL" meta-resource, fields are to be grabbed through a GRIB2 nomenclature
t_as_grib2 = dict(discipline=0, parameterCategory=0, parameterNumber=0, # this is temperature
typeOfFirstFixedSurface=119) # hybrid-pressure levels
[4]:
t3d = r.readfield(t_as_grib2)
[5]:
t3d.sp2gp()
[6]:
profiles = []
for position in [(2,45), (5,50), (-3,47)]:
profgeo = t3d.geometry.make_profile_geometry(*position)
profiles.append(t3d.extract_subdomain(profgeo))
[7]:
fig = ax = None
for p in profiles:
fig, ax = p.plotfield(labels='Temperature @ ({}, {})'.format(p.geometry.grid['longitudes'][0],
p.geometry.grid['latitudes'][0]),
over=(fig, ax)) # superpose profiles
/usr/lib64/python2.7/site-packages/matplotlib/figure.py:1744: UserWarning: This figure includes Axes that are not compatible with tight_layout, so its results might be incorrect.
warnings.warn("This figure includes Axes that are not "
and now on an Altitude z-axis¶
[8]:
from epygram.geometries.VGeometry import hybridP2altitude
# for the conversions, we need R and T profiles, so (at least) q and T
ps = epygram.formats.resource(r.resource.container.abspath, 'r').readfield('SURFPRESSION')
q_as_grib2 = dict(discipline=0, parameterCategory=1, parameterNumber=0, # this is specific humidity
typeOfFirstFixedSurface=119) # hybrid-pressure levels
q3d = r.readfield(q_as_grib2)
[9]:
ps.sp2gp()
ps.operation('exp')
[10]:
from bronx.meteo.conversion import q2R
[11]:
for p in profiles:
q_profile = q3d.extract_subdomain(p.geometry).data # profile of q at the same place
# specific humidity (and optionally hydrometeors) to moist air specific constant R
R_profile = q2R(q_profile)
t_profile = p.data
ps_local = ps.getvalue_ll(p.geometry.grid['longitudes'][0], p.geometry.grid['latitudes'][0])
p.geometry.vcoordinate = hybridP2altitude(p.geometry.vcoordinate,
R_profile,
t_profile,
ps_local,
vertical_mean='geometric',
Pdep=None) # we could also provide a (NH) Pressure-departure profile for more accuracy
[12]:
fig = ax = None
for p in profiles:
fig, ax = p.plotfield(labels='Temperature @ ({}, {})'.format(p.geometry.grid['longitudes'][0],
p.geometry.grid['latitudes'][0]),
over=(fig, ax)) # superpose profiles
