Skip to content

Instantly share code, notes, and snippets.

@mathomp4
Last active July 27, 2026 14:49
Show Gist options
  • Select an option

  • Save mathomp4/756d9d2a9aa3c393b41d4d5b891ac0af to your computer and use it in GitHub Desktop.

Select an option

Save mathomp4/756d9d2a9aa3c393b41d4d5b891ac0af to your computer and use it in GitHub Desktop.

MERRA-2 Zonal Mean File Conversion Analysis

Executive Summary

This document details the reverse-engineering analysis of how 3D MERRA-2 NetCDF files (on a $576 \times 361 \times 72$ grid) are processed to produce the corresponding zonal-mean files (MERRA-2.inst3_3d_asm_Nv.monthly.YYYYMM.0001x0091x0072.nc4 on a $1 \times 91 \times 72$ grid), specifically focusing on ozone ($O_3$/$OX$).

It includes comparison between two input sources:

  1. Pre-made Monthly Means (MERRA-2.inst3_3d_asm_Nv.monthly.YYYYMM.nc4)
  2. Daily Instantaneous Files (MERRA2_400.inst3_3d_asm_Nv.YYYYMMDD.nc4, 3-hourly, 240 timesteps/month)

The June 2012 validation changes an important earlier conclusion: the dominant difference from the canonical ozone field is a unit conversion, not numerical precision or the use of monthly rather than daily input. The canonical field is approximately the zonal mean of MERRA-2 mass mixing ratio converted to volume (mole) mixing ratio:

$$ OX_{\mathrm{canonical}} \approx \frac{28.97}{48.0},\overline{O_3}_{\mathrm{mass}} \approx 0.60354,\overline{O_3}_{\mathrm{mass}}. $$

Here 28.97 and 48.0 are the molecular weights of dry air and ozone, respectively. The canonical per-month files retain O3:units = "kg kg-1", but the downstream pchem.species...MERRA2OX file correctly identifies the resulting OX field as mol mol-1.


1. Grid & Dimensional Comparison

Dimension Original Parent File Zonal Mean File (0001x0091x0072) Transformation Method
Longitude (lon) $576$ points ($-180^\circ$ to $179.375^\circ$, step $0.625^\circ$) $1$ point ($-180^\circ$) Unweighted arithmetic mean across all 576 longitudes
Latitude (lat) $361$ points ($-90^\circ$ to $90^\circ$, step $0.5^\circ$) $91$ points ($-90^\circ$ to $90^\circ$, step $2.0^\circ$) Direct sub-sampling (lat[::4])
Vertical Level (lev) $72$ model layer indices (1..72) $72$ nominal layer mid-point pressures in hPa ($\sim 0.015$ to $\sim 976.62\text{ hPa}$) Layer index coordinate populated with hPa pressure values

2. Input Data Sources: Daily Files vs Pre-Made Monthly Means

Approach 1: Pre-Made Monthly 3D File (MERRA-2.inst3_3d_asm_Nv.monthly.YYYYMM.nc4)

  • Pros: Fastest processing (under 2 seconds per month), small disk I/O (single 800 MB file).
  • Accuracy: For June 2012, the longitude mean alone differs from canonical O3 by up to $1.399 \times 10^{-5}$; this is expected because it is still in kg kg-1 while canonical O3/OX is effectively in mol mol-1.
  • After conversion: multiplying by $28.97/48.0$ reduces the June 2012 maximum absolute residual to $1.407 \times 10^{-8}$ and the mean absolute residual to $3.50 \times 10^{-10}$.

Approach 2: Aggregating Daily Instantaneous Files (MERRA2_400.inst3_3d_asm_Nv.YYYYMMDD.nc4)

  • Pros: Allows an independent check of the temporal averaging and is appropriate if a monthly product is unavailable or does not use the needed sampling/masking convention.
  • June 2012 result: mean over all 240 three-hourly samples and 576 longitudes agrees with the zonal mean of the pre-made monthly field to within $3.64 \times 10^{-12}$ before the molecular-weight conversion.
  • Conclusion for O3: daily input does not explain the canonical discrepancy for this month. Time and longitude means commute for these data; use the monthly field unless a future validation finds a non-linear mask or weighting difference.
  • I/O Overhead: Requires opening 30 daily files (~24 GB total input per month).

3. Mathematical Details & Numerical Verification

Zonal Mean (Longitude)

An unweighted arithmetic mean is performed over all 576 longitudes:

$$ \text{ZonalMean}(\text{lev}, \text{lat}_{361}) = \frac{1}{576} \sum_{i=1}^{576} \text{Field}(\text{lev}, \text{lat}_{361}, \text{lon}_i) $$

Latitude Sub-sampling

Latitude values are sampled directly at a stride of 4 from the 361-point grid down to the 91-point grid:

$$ \text{Child}(\text{lev}, j) = \text{ZonalMean}(\text{lev}, 4 \times j) \quad \text{for } j = 0, 1, \dots, 90 $$

This selects latitudes $-90^\circ, -88^\circ, -86^\circ, \dots, +90^\circ$.

Ozone mass-to-volume conversion

The field to compare with the legacy OX forcing file is a volume mixing ratio, not the native MERRA-2 O3 mass mixing ratio:

$$ OX(\mathrm{mol\ mol^{-1}}) = O_3(\mathrm{kg\ kg^{-1}})\frac{M_{\mathrm{air}}}{M_{O_3}} \approx O_3\frac{28.97}{48.0}. $$

For a strict replacement workflow, apply this conversion after the longitude mean and before writing OX. Preserve the target file's pressure-level coordinate and latitude units (radians in pchem.species.CMIP-5.MERRA2OX.197902-201706.z_91x72.nc4).

June 2012 quantitative comparison

All figures below use the 72 by 91 O3/OX grid (6,552 values). daily is the mean over the 30 June daily files, all 240 three-hourly samples, and all 576 longitudes. monthly is the longitude mean of MERRA-2.inst3_3d_asm_Nv.monthly.201206.nc4; both are still mass mixing ratio until the indicated conversion.

Comparison Maximum absolute difference Mean absolute difference RMSE 99th-percentile absolute difference Maximum relative difference
daily O3 minus monthly O3 $2.97\times10^{-12}$ $1.22\times10^{-13}$ $2.50\times10^{-13}$ $9.58\times10^{-13}$ $1.65\times10^{-7}$
monthly O3 $\times 28.97/48$ minus canonical O3 $1.41\times10^{-8}$ $3.50\times10^{-10}$ $7.28\times10^{-10}$ $2.70\times10^{-9}$ $1.31\times10^{-2}$
daily O3 $\times 28.97/48$ minus canonical O3 $1.41\times10^{-8}$ $3.50\times10^{-10}$ $7.28\times10^{-10}$ $2.70\times10^{-9}$ $1.31\times10^{-2}$
canonical O3 versus legacy pchem OX (same new-level grid) 0 0 0 0 0

For daily minus monthly O3, 55 of 6,552 values exceed $10^{-12}$ and none exceed $10^{-11}$. Therefore daily input cannot materially improve the June 2012 replacement: its residual is about four orders of magnitude smaller than the residual remaining after the required mass-to-volume conversion.

Operational recommendation. Use the pre-made monthly 3D field for production zonal means. This is established for June 2012, not a mathematical guarantee for every future MERRA collection: retain a one-month daily-vs-monthly sentinel test when changing MERRA collection/version, variable, missing-value treatment, or temporal weighting.


4. How to Reproduce (Verified Code Examples)

Option A1: Python (xarray) - From Pre-Made Monthly File (Fastest)

import xarray as xr

# Open parent 3D monthly file
ds = xr.open_dataset('MERRA-2.inst3_3d_asm_Nv.monthly.201206.nc4', decode_times=False)

# 1. Select native O3 mass mixing ratio.
# 2. Arithmetic mean across longitude.
# 3. Convert mass mixing ratio to volume mixing ratio.
# 4. Sub-sample latitude (361 -> 91).
ox = (ds['O3'].mean(dim='lon') * (28.97 / 48.0)).isel(lat=slice(None, None, 4))
ox.name = 'OX'
ox.attrs.update(units='mol mol-1', long_name='odd_oxygen_volume_mixing_ratio')
ds_o3 = ox.to_dataset()

# 5. Expand lon dimension to length 1 (-180.0) & restore dimension order
ds_o3 = ds_o3.expand_dims(lon=[-180.0]).transpose('time', 'lev', 'lat', 'lon')

# Save to NetCDF4
ds_o3.to_netcdf('MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.OX.python.nc4', format='NETCDF4_CLASSIC')

Option A2: Python (xarray) - From 30 Daily Instantaneous Files (Full Time-Series Aggregation)

import xarray as xr
import glob

# Find all 30 daily files for the target month
daily_files = sorted(glob.glob('/discover/nobackup/projects/gmao/merra2/data/pub/products/MERRA2_400/Y2012/M06/MERRA2_400.inst3_3d_asm_Nv.201206*.nc4'))

# Open multi-file dataset across all 240 timesteps
ds_daily = xr.open_mfdataset(daily_files, combine='by_coords')

# Mean across both time and longitude, convert units, then sub-sample latitude
ox = (ds_daily['O3'].mean(dim=['time', 'lon']) * (28.97 / 48.0)).isel(lat=slice(None, None, 4))
ox.name = 'OX'
ox.attrs.update(units='mol mol-1', long_name='odd_oxygen_volume_mixing_ratio')
ds_o3 = ox.to_dataset()
ds_o3 = ds_o3.expand_dims(time=[ds_daily.time.values[0]], lon=[-180.0])
ds_o3 = ds_o3.transpose('time', 'lev', 'lat', 'lon')

# Save to NetCDF4
ds_o3.to_netcdf('MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.OX.from_dailies.nc4', format='NETCDF4_CLASSIC')

Option B1: Climate Data Operators (CDO) - From Monthly File

cdo -mulc,0.6035416666666667 -chname,O3,OX -select,name=O3 -zonmean -samplegrid,4 \
    MERRA-2.inst3_3d_asm_Nv.monthly.201206.nc4 \
    MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.OX.cdo.nc4

Option B2: Climate Data Operators (CDO) - From Daily Files

cdo -mulc,0.6035416666666667 -chname,O3,OX -select,name=O3 -zonmean -samplegrid,4 -ensmean \
    /discover/nobackup/projects/gmao/merra2/data/pub/products/MERRA2_400/Y2012/M06/MERRA2_400.inst3_3d_asm_Nv.201206*.nc4 \
    MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.OX.cdo_dailies.nc4

Option C: GrADS / lats4d

export PATH=/discover/nobackup/projects/gmao/share/dasilva/opengrads/Contents:$PATH

lats4d.sh -i MERRA-2.inst3_3d_asm_Nv.monthly.201206.nc4 \
          -o MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.O3.lats4d \
          -vars o3 \
          -func "(28.97/48.0)*ave(@,x=1,x=576)" \
          -ftype sdf

nccopy -d 5 MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.O3.lats4d.nc \
          MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.O3.lats4d.nc4
rm MERRA-2.inst3_3d_asm_Nv.monthly.201206.0001x0091x0072.O3.lats4d.nc

5. Software Tooling & Historical Origin

  • Global Attributes:
    • Title: "GEOS-5 GCM (DTOA Conversion)"
    • History: "File written by GFIO v1.0.8"
  • Software Chain: These attributes are automatically assigned by GrADS / lats4d (Arlindo da Silva's LATS-based NetCDF exporter tool built on GMAO's GFIO library) or the GMAO dtoa (Data Transfer / Output Assembly) post-processing tool.
  • The control file xdf.tabl was used by OpenGrADS to template and access these monthly zonal-mean files (%y4%m2).

6. Metadata & Tolerance Comparison Notes (nccmp)

When running nccmp -dmfgsB between generated files and canonical files:

  • Format & Metadata: nccmp flags global attribute differences (Source, Contact, Title, History) and NetCDF attribute conventions (_FillValue=nan vs 1.0e15f).
  • Data interpretation: an unconverted O3 comparison is invalid because it compares kg kg-1 to the legacy OX-like mol mol-1 field. The June 2012 unconverted maximum difference is $1.399 \times 10^{-5}$, not float32 roundoff.
  • Legacy forcing-file check: June 2012 OX in pchem.species.CMIP-5.MERRA2OX.197902-201706.z_91x72.nc4 is bit-for-bit identical to the corresponding legacy canonical OX/new-levels field (6,552 of 6,552 values). Thus replacing OX requires preserving its mol mol-1 convention.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment