...
Spherical harmonic transform available in NCL and pyNGL for lat/lon data.
Interpolate EAM native grid output to a (N+1)x2N lat/lon cap grid with TempestRemap’s “highorder” (native SE shape function) algorithm.
For cubed sphere grid with resolution NE, take N>=NE*6 (degrees of freedom, pole to pole)
Need instantaneous vorticity, divergence and PS on the GLL grid (not the PG2) grid. PS is needed if interpolating to pressure levels
in EAM, use: ‘VOR:I','DIV:I','DYN_PS:I'
Need ‘DYN_PS’ instead of ‘PS’ to get PS on the GLL grid instead of the PG2 grid.
For smooth results, average over at least 1 months month (post spinup) with 12h snapshots. Running longer rather than with higher frequency I/O results in smoother spectra in Aqua planet. Are there seasonal effects in Real planet?
At NE256, interpolating and computing the spectra for each snapshot takes about 5min and 30GB of memory. At NE1024 this will probably require ~500GB of memory and could take 320min per snapshot
At NE1024 needed to use homme_tool for interpolation (running on 1024 KNL nodes: 20min init time + 5min per 3D snapshot). NCL spherical harmonic transform code required 100GB of memory, taking about 3h per transform.
KE Spectra
KE spectra is useful for understanding and tuning the dissipation mechanisms in the model. For a global model, it can be computed via a vector spherical harmonic transform of the instantaneous velocity field (or spherical harmonic transforms of the instantaneous vorticity and divergence scalars) on a spherical surface. For all the data below, we interpolate to 250mb (I think because this is the most energetic part of the atmosphere - check this?). Spherical harmonics can be thought of as polynomials in Cartesian coordinates (x,y,z) restricted to the sphere. For each degree k there will be 2k+1 polynomials (spherical harmonics) of total degree k. To compute the KE spectra, E(k), we sum the power (coefficient squared) over all degree k spherical harmonics.
...
Example KE spectra. Plotted is the “compensated” spectra, E(k)*k^5/3, to better illustrate the possible transition to ^-5/3. E(k) was computed from many flow snapshots and then E(k) is averaged over all the results to get the nice curves shown below. Unfortunately, one cannot compute E(k) from time averaged flow fields.
...
In CAM4 at 1/8 degree (14km) resolution the model had a nice Nastrom-Gage transition. This has slowly deteriorated, perhaps due to increased dissipation being added to the model as it the physics has evolved. It may also have been that CAM4 had this transition for the wrong reasons. The breakdown of the spectra into divergence and vorticity components (right) shows that most of the differences are due to small scale divergence - associated with wave activity. CAM5.3 is EAM V0. It will be interesting to see the specta spectra for higher resolutions now that we are able to run at 6km and 3km resolutions.
...
Computing the KE spectra requires performing spherical harmonic transforms. This can be done via the long standing SPHEREPACK Fortran package ( https://www2.cisl.ucar.edu/resources/legacy/spherepack ), which has an easy to use interface in NCL (and pyNGL?as of 2020/7, not available in pyNGL).
SPHEREPACK requires gridded data on either a Gaussian grid or an equal angle lat/lon “cap” grid which includes points at the poles. It relies on a combination of associated Legendre transforms and FFTs. For efficiency the grid should have dimensions N x 2N (Gaussian), or (N+1) x 2N (equal angle cap) with N the product of powers of small primes. The most common lat/lon data is an equal angle offset grid, which avoids pole points (see. CAM-FV Grid Overview ). Data on this offset grid would have to first be interpolated to a cap grid before computing spherical harmonic transforms.
...
Computing E(k) at lower resolution or with certain algorithms introduces a roll off (reduced energy) and high frequencies. Since the behavior of the spectra, especially the rolloff starting around wave number 200 is of most interest, this should be avoided.
Using “PG2” output (cyan curve below) is not recommended - this is most likely due to the downscaling downsampling of the velocity when transforming from the dynamics to the PG2 physics grid, removing energy from small scales.
“intbilin” algorithm is not recommended (brow and green curves below). This algorithm is very nice for analysis, but does remove energy from small scales. “bilin” is a little better (orange, blue)
Best results are obtained by the highorder algorithm (GLL 2048HO and GLL 3072HO).
...