Using air pressure as the vertical coordinate in ParaView

[huiwanpnnl]

In atmospheric science, it is useful and common to examine fields on isobaric surfaces (i.e., on a surface in the atmosphere where the air pressure is a constant value).

EAM uses a pressure-based terrain-following coordinate for the vertical discretization. The grid layers are not isobaric near the Earth's surface (although they do transition to isobaric higher up).

In EAM's output files, we have the 1D coordinate variables lev and ilev, which are sometimes referred to as the nominal pressure - they are the air pressure at layer midpoints and layer interfaces assuming the surface press (PS) is 1000 hPa. But because the actual PS can deviate substantially from 1000 hPa, especially over the high mountains, the actual air pressure in specific grid cells can deviate substantially from lev or ilev. In order to visualize EAM output on isobaric surfaces without having to vertically interpolate the output beforehand (i.e., before loading the data files in ParaView), Phil Rasch suggested adding the option of using the actual pressure field for the vertical coordinate inside ParaView.

I have checked some EAM output files and the source code and realized the code does not write the pressure field to output. I think the reason is that pressure can be calculated using PS and the coordinate parameters hyam, hybm, hyai, hybi. The formulas are given below:

For layer interfaces,

p_int(itime, klev, icol) =  hyai(klev) * P0 + hybi(klev) * PS(itime,icol)

For layer midpoints,

p_mid(itime, klev, icol) =  hyam(klev) * P0 + hybm(klev) * PS(itime,icol)

I've checked our sample output PD_1800_ad4fd8_ANN_climo_SE.nc shared at https://compy-dtn.pnl.gov/wanh895/ParaView_EAM/cubed_sphere/model_output_samples/ and confirmed that it has all the needed variables PS, hyam, hybm, hyai, hybi and the constant PS.

Can we implement the pressure-coordinate option using the pressure values calculated on the fly? (Tagging @berkgeveci, @cjh1, @bnmajor)

[cjh1]

Absolutely, we can do that.

[cjh1]

Do we a similar calculation when using Z3?

[huiwanpnnl]

Hmm, in principle we can, if Z3 is not an available variable in the same NetCDF file.

For now, since we do have Z3 in the NetCDF file, can we simply use its values? I think this case of "using one variable as the vertical coordinate for visualizing another variable" can have general value.

[cjh1]

Of course, I think I was just a little confused as the the relationskip of lev to Z3, but I see it is Z3(time, lev, ncol)