0

I am trying to calculate the vertical integral of the horizontal advection of Moist Static Energy uisng ERA5 data over a particular region. I am using xarray hence I know how to interpolate and integrate:

ds6 = ds5.sel(year=2000)
m = Cp*ds6.t + Lv*ds6.q + ds6.z
hadv = mpcalc.advection(m, u=ds6.u, v=ds6.v)
inter = hadv.interp(level=xs, method="cubic")
inter.integrate('level')

The above integrates from 1000 mb to 100 mb. I have surface pressure data for each coordinate. How do I do the integral at each coordinate only from the surface pressure to 100 mb instead?

Here is what I have tried:

Ps = ds10.sp.sel(time='2021-12-01')
selected = inter.where(inter.level < Ps) #selecting data points above Ps only
selected.dropna(dim='level') #removing Nan so that we can integrate

This does not work because it removes all data points below the maximum Ps of any data point, and spatial variability because of Ps is gone. What do I do?

2
  • hi, doesn't answer your question but I think ERA5 has vertical integrals of energy fluxes and moisture fluxes in case that helps (see under section "vertical integrals") Commented Jun 5, 2023 at 11:07
  • Thanks for the reply! The available vertical integrals do not have what I want unfortunately, they only provide vertical integrals of basic quantities. I wish to calculate the vertical integral of MSE,
    – Rakshith
    Commented Jun 6, 2023 at 4:50

1 Answer 1

1

Say you have the following dataset ds:

<xarray.Dataset>
Dimensions:  (lon: 10, lat: 20, level: 50)
Coordinates:
  * lon      (lon) int64 0 1 2 3 4 5 6 7 8 9
  * lat      (lat) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
  * level    (level) float64 1.1e+03 1.078e+03 1.056e+03 ... 54.49 32.24 10.0
Data variables:
    sfc_p    (lon, lat) int64 1031 980 1047 1028 988 ... 990 978 1026 1055 1003
    adv      (lon, lat, level) float64 -0.9586 -0.04061 1.627 ... -1.318 0.8371

Then you can apply a function to every individual vertical profile using groupby. Check this out:

all_dims_except_level = [d for d in ds.dims if (d != "level")]
adv_gridpoint = ds.stack(gridpoint=all_dims_except_level)


def integral_from_sfc_p_to_p(data, target_pressure):
    sfc_p_at_gridpoint = data.sfc_p.values.item()

    levels_between_sfc_p_and_upper_limit = data.level.where(
        (data.level > target_pressure) & (data.level < data.sfc_p), drop=True
    ).values

    target_levels = np.concatenate(
        [[sfc_p_at_gridpoint], levels_between_sfc_p_and_upper_limit, [target_pressure]]
    )

    data_interp_to_sfc_p = data.interp(level=target_levels)
    data_integrated = data_interp_to_sfc_p.integrate("level")
    return data_integrated


print(
    adv_gridpoint.groupby("gridpoint")
    .apply(integral_from_sfc_p_to_p, target_pressure=100)
    .unstack()
)

yields:

<xarray.Dataset>
Dimensions:  (lon: 10, lat: 20)
Coordinates:
  * lon      (lon) int64 0 1 2 3 4 5 6 7 8 9
  * lat      (lat) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
Data variables:
    sfc_p    (lon, lat) int64 1031 980 1047 1028 988 ... 990 978 1026 1055 1003
    adv      (lon, lat) float64 -172.2 61.96 -109.3 22.05 ... -23.87 182.0 235.1

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.