Usuage request: Sigma to pressure conversion using ERA-I

[avatar101]

Usage request

I'm trying to interpolate ECMWF ERA-I data which is on hybrid sigma levels to isobaric levels using metpy.log_interpolate_idfollowing this example. Below is my try and the description of dataset.

ds = xr.open_dataset('P20001201_00').load()
print(ds)
>><xarray.Dataset>
Dimensions:  (lat: 181, lev: 60, lev_2: 1, lon: 360, nhyi: 61, nhym: 60, time: 1)
Coordinates:
  * lon      (lon) float64 -180.0 -179.0 -178.0 -177.0 ... 177.0 178.0 179.0
  * lat      (lat) float64 -90.0 -89.0 -88.0 -87.0 -86.0 ... 87.0 88.0 89.0 90.0
  * lev      (lev) float64 1.0 2.0 3.0 4.0 5.0 6.0 ... 56.0 57.0 58.0 59.0 60.0
  * lev_2    (lev_2) float64 1.0
  * time     (time) datetime64[ns] 2000-12-01
Dimensions without coordinates: nhyi, nhym
Data variables:
    hyai     (nhyi) float64 0.0 20.0 38.43 63.65 95.64 ... 65.89 7.368 0.0 0.0
    hybi     (nhyi) float64 0.0 0.0 0.0 0.0 0.0 ... 0.9883 0.994 0.9976 1.0
    hyam     (nhym) float64 10.0 29.21 51.04 79.64 ... 138.1 36.63 3.684 0.0
    hybm     (nhym) float64 0.0 0.0 0.0 0.0 0.0 ... 0.984 0.9911 0.9958 0.9988
    T        (time, lev, lat, lon) float32 -59.238907 -59.238907 ... -20.415268
    OMEGA    (time, lev, lat, lon) float32 -1.1673663e-05 ... -0.008483887
    U        (time, lev, lat, lon) float32 -3.511917 -3.402542 ... 4.682083
    V        (time, lev, lat, lon) float32 -6.4184265 -6.4770203 ... -5.4629974
    Q        (time, lev, lat, lon) float32 4.044671e-06 ... 0.00056174514
    PS       (time, lev_2, lat, lon) float32 696.69763 696.69763 ... 1015.3609
    LNSP     (time, lev_2, lat, lon) float32 11.151522 11.151522 ... 11.52817
Attributes:
    CDI:                       Climate Data Interface version 1.6.4 (http://c...
    Conventions:               CF-1.4
    history:                   Fri Sep 30 04:11:57 2016: cdo -f nc merge P200...
    institution:               European Centre for Medium-Range Weather Forec...
    constants_file_name:       P20001201_00
    NCO:                       20160930
    nco_openmp_thread_number:  1
    CDO:                       Climate Data Operators version 1.6.4 (http://c...

# converting to pressure to model levels
# [Description on ECMWF model levels](https://rda.ucar.edu/datasets/ds627.2/docs/Eta_coordinate/)
p_sigma = ds.hyai + ds.hybi * ds.PS 
p_sigma_units = p_sigma.values * units.hPa 
plevs = [250.] * units.hPa  # desired level
U, V = ds.U.values * units.meter/units.second, ds.V.values * units.meter/units.second
# begin interpolation
u_inter, v_inter = log_interpolate_1d(plevs, p_sigma_units, U, V, axis=1)

• Problem description: Maybe I'm doing some fundamental error in pressure conversion? I'm getting a runtime error of dividing by zero and also getting the following type error:

TypeError                                 Traceback (most recent call last)
<ipython-input-21-3326c5c4c9e3> in <module>
----> 1 u_inter, v_inter = log_interpolate_1d(plevs, p_sigma_units, U, V, axis=1)

/usr/local/anaconda3/lib/python3.6/site-packages/metpy/xarray.py in wrapper(*args, **kwargs)
    434         kwargs = {name: (v.metpy.unit_array if isinstance(v, xr.DataArray) else v)
    435                   for name, v in kwargs.items()}
--> 436         return func(*args, **kwargs)
    437     return wrapper
    438 

/usr/local/anaconda3/lib/python3.6/site-packages/pint/registry_helpers.py in wrapper(*values, **kw)
    192             new_values, values_by_name = converter(ureg, values, strict)
    193 
--> 194             result = func(*new_values, **kw)
    195 
    196             if container:

/usr/local/anaconda3/lib/python3.6/site-packages/metpy/interpolate/one_dimension.py in log_interpolate_1d(x, xp, *args, **kwargs)
    226     log_x = np.log(x)
    227     log_xp = np.log(xp)
--> 228     return interpolate_1d(log_x, log_xp, *args, axis=axis, fill_value=fill_value)

/usr/local/anaconda3/lib/python3.6/site-packages/metpy/xarray.py in wrapper(*args, **kwargs)
    434         kwargs = {name: (v.metpy.unit_array if isinstance(v, xr.DataArray) else v)
    435                   for name, v in kwargs.items()}
--> 436         return func(*args, **kwargs)
    437     return wrapper
    438 

/usr/local/anaconda3/lib/python3.6/site-packages/pint/registry_helpers.py in wrapper(*values, **kw)
    192             new_values, values_by_name = converter(ureg, values, strict)
    193 
--> 194             result = func(*new_values, **kw)
    195 
    196             if container:

/usr/local/anaconda3/lib/python3.6/site-packages/metpy/interpolate/one_dimension.py in interpolate_1d(x, xp, *args, **kwargs)
    117     xp = xp[sorter]
    118     # Ensure pressure in increasing order
--> 119     variables = [arr[sorter] for arr in args]
    120 
    121     # Make x broadcast with xp

/usr/local/anaconda3/lib/python3.6/site-packages/metpy/interpolate/one_dimension.py in <listcomp>(.0)
    117     xp = xp[sorter]
    118     # Ensure pressure in increasing order
--> 119     variables = [arr[sorter] for arr in args]
    120 
    121     # Make x broadcast with xp

TypeError: tuple indices must be integers or slices, not tuple

• Versions.
  • python --version 3.6.7
  • `python -c 'import metpy; print(metpy.version)' 0.10.0

[dopplershift]

I'm pretty sure this is a duplicate of #998 and is caused by pint 0.9. Can you try installing pint=0.8.1 and see if that fixes your problem?

[deepwindlee]

hello i meet the similar problem,
the code is from help(log_interpolate_1d), but it raise error
x_log = np.array([1e3, 1e4, 1e5, 1e6])
y_log = np.log(x_log) * 2 + 3
x_interp = np.array([5e3, 5e4, 5e5])
print(x_log.shape)
print(y_log.shape)
#metpy.calc.log_interp(x_interp, x_log, y_log)
type(x_log)
metpy.calc.log_interp(x_interp, x_log, y_log)

---

TypeError: tuple indices must be integers or slices, not tuple

python --version 3.6.4
`python -c 'import metpy; print(metpy.version)' 0.10.0

[dopplershift]

@deepwindlee Can you install pint 0.8.1? Either:
pip install pint==0.8.1
or
conda install -c conda-forge pint=0.8.1
depending on whether you've installed metpy using pip or conda.

[zbruick]

This was a duplicate of #997 and fixed by #1087. Going to close, but this can be reopened if an issue still exists.