TANGERANG SELATAN WEATHER

Sabtu, 23 Mei 2026

Menghitung Rasio Pencampuran Jenuh dari ERA5

Latar gradien atmosfer biru-hijau melambangkan kapasitas uap air di troposfer tropis

Mengapa Rasio Pencampuran Jenuh Penting untuk Prakiraan Cuaca

Berapa banyak uap air yang bisa ditampung udara tropis sebelum turun sebagai hujan? Jawabannya ditentukan oleh saturation mixing ratio \(w_s\) — batas atas kandungan uap air yang bisa dipegang oleh satu kilogram udara kering pada suhu dan tekanan tertentu. Ketika \(w\) aktual mendekati \(w_s\), hujan lebat menjadi sangat mungkin terjadi.

Di Indonesia, SST warm-pool yang secara konsisten berada di kisaran 28–30 °C mendorong nilai \(w_s\) permukaan hingga mendekati 27–30 g/kg — jauh di atas rata-rata global. Kaitan antara suhu dan kapasitas uap air ini berakar pada hubungan Clausius–Clapeyron: setiap kenaikan 1 K suhu meningkatkan tekanan uap jenuh sekitar 6–7%, sehingga \(w_s\) ikut naik. Artinya di dunia yang lebih hangat, atmosfer dapat menyimpan lebih banyak kelembapan sebelum mencapai saturasi.

Monthly average total column water vapour (TCWV) from NASA AIRS, July 2025, showing high values over tropical oceans Sumber: NASA GES DISC / Giovanni; Mitchell, B. (2025), NASA Earthdata Blog (tautan)

Tutorial ini menunjukkan cara menghitung \(w_s\) dari data ERA5 Indonesia secara end-to-end: mulai dari download via cdsapi, konversi unit, sampai peta kontur dengan cartopy, dan ringkasan musiman DJF vs JJA.

Definisi Rasio Pencampuran dan Konstanta 0,622

Mixing ratio \(w\) adalah massa uap air per satuan massa udara kering, dinyatakan dalam g/kg. Ini berbeda dari specific humidity \(q\), yang menggunakan massa udara lembap (basah) sebagai penyebut:

$$q = \frac{w}{1 + w}$$

Untuk nilai \(w\) yang umum di atmosfer (\(w < 40\) g/kg), selisih antara \(w\) dan \(q\) hanya beberapa persen, sehingga keduanya bisa dianggap setara untuk keperluan praktis.

Saturation mixing ratio \(w_s\) adalah nilai maksimum \(w\) yang mungkin ada pada suhu dan tekanan tertentu — yaitu mixing ratio ketika udara sepenuhnya jenuh. Rumusnya:

$$w_s = 0{,}622 \cdot \frac{e_s}{p - e_s}$$

Konstanta \(0{,}622\) (disebut \(\varepsilon\) atau epsilon) berasal dari perbandingan massa molar:

$$\varepsilon = \frac{M_{\text{uap air}}}{M_{\text{udara kering}}} = \frac{18{,}02}{28{,}97} \approx 0{,}622$$

Konstanta yang sama juga merupakan rasio konstanta gas \(R_d / R_v\). Karena tekanan uap biasanya kurang dari 7% tekanan total atmosfer, dalam praktek \(p - e_s \approx p\), sehingga \(w_s \approx q_s\) — perbedaannya di bawah ambang batas signifikan untuk analisis ERA5.

Formula Bolton untuk Tekanan Uap Jenuh

Untuk menghitung \(e_s\) dari suhu, kita gunakan Formula Bolton (1980) Persamaan 10 — formula empiris dengan akurasi \(\pm 0{,}3\%\) untuk rentang suhu \(-35\ \text{°C}\) hingga \(35\ \text{°C}\):

$$e_s(T_C) = 6{,}112 \exp\!\left(\frac{17{,}67 \cdot T_C}{T_C + 243{,}5}\right) \text{ hPa}$$

di mana \(T_C\) adalah suhu dalam Celsius. Beberapa nilai referensi yang berguna untuk sanity-check: pada \(T_C = 25\ \text{°C}\) dan \(p = 1013\ \text{hPa}\), \(e_s \approx 31{,}67\ \text{hPa}\) dan \(w_s \approx 20\ \text{g/kg}\); pada \(T_C = 30\ \text{°C}\) dan \(p = 1010\ \text{hPa}\), \(e_s \approx 42{,}5\ \text{hPa}\) dan \(w_s \approx 27\ \text{g/kg}\).

Snippet berikut mengimplementasikan kedua formula sebagai fungsi Python yang bisa dipakai ulang di seluruh artikel:

import numpy as np

def bolton_es(t_c):
    """Saturation vapour pressure (hPa) via Bolton (1980) Eq. 10.
    t_c : temperature in Celsius (scalar or array)
    Returns e_s in hPa. Accuracy +/- 0.3% for -35 <= t_c <= 35 degC.
    """
    return 6.112 * np.exp(17.67 * t_c / (t_c + 243.5))

def saturation_mixing_ratio(t_c, p_hpa):
    """Saturation mixing ratio w_s in g/kg.
    t_c   : temperature in Celsius
    p_hpa : total pressure in hPa
    """
    es = bolton_es(t_c)
    ws_kgkg = 0.622 * es / (p_hpa - es)   # kg/kg
    return ws_kgkg * 1000.0                # g/kg

# --- sanity checks ---
tests = [
    (25.0, 1013.0, "Jakarta coast ~sea level"),
    (30.0, 1010.0, "Tropical warm-pool surface"),
    (15.0, 1013.0, "Temperate reference (Bolton: e_s ~ 17.04 hPa)"),
]
print(f"{'T (°C)':>8}  {'p (hPa)':>8}  {'e_s (hPa)':>10}  {'w_s (g/kg)':>11}  note")
print("-" * 65)
for t, p, note in tests:
    es = bolton_es(t)
    ws = saturation_mixing_ratio(t, p)
    print(f"{t:8.1f}  {p:8.1f}  {es:10.3f}  {ws:11.3f}  {note}")
  T (°C)   p (hPa)   e_s (hPa)   w_s (g/kg)  note
-----------------------------------------------------------------
    25.0    1013.0      31.674       20.076  Jakarta coast ~sea level
    30.0    1010.0      42.456       27.293  Tropical warm-pool surface
    15.0    1013.0      17.040       10.642  Temperate reference (Bolton: e_s ~ 17.04 hPa)

Output di atas mengonfirmasi nilai referensi dari dossier penelitian: \(w_s \approx 20\ \text{g/kg}\) pada 25 °C dan \(\approx 27\ \text{g/kg}\) pada 30 °C — angka yang konsisten dengan kondisi permukaan warm-pool Indonesia.

Mengunduh dan Membuka Data ERA5 Indonesia

Kita download dua variabel ERA5 dari cds.climate.copernicus.eu: suhu 2 m (t2m, dalam K) dan mean sea-level pressure (msl, dalam Pa). Keduanya tersedia di product reanalysis-era5-single-levels. Snippet ini memakai cache guard if not os.path.exists(OUT) sehingga download hanya berjalan sekali — iterasi berikutnya langsung baca file lokal.

import os, cdsapi, xarray as xr

T2M = "era5_t2m_indonesia_2024_6h.nc"
MSL = "era5_msl_indonesia_2024_6h.nc"

bbox   = [6, 95, -11, 141]   # N, W, S, E — Indonesia
months = [f"{m:02d}" for m in range(1, 13)]
days   = [f"{d:02d}" for d in range(1, 32)]
times  = ["00:00", "06:00", "12:00", "18:00"]

c = cdsapi.Client(quiet=True)

if not os.path.exists(T2M):
    c.retrieve(
        "reanalysis-era5-single-levels",
        {
            "product_type": "reanalysis",
            "variable": "2m_temperature",
            "year": "2024",
            "month": months,
            "day": days,
            "time": times,
            "area": bbox,
            "format": "netcdf",
        },
        T2M,
    )

if not os.path.exists(MSL):
    c.retrieve(
        "reanalysis-era5-single-levels",
        {
            "product_type": "reanalysis",
            "variable": "mean_sea_level_pressure",
            "year": "2024",
            "month": months,
            "day": days,
            "time": times,
            "area": bbox,
            "format": "netcdf",
        },
        MSL,
    )

ds_t = xr.open_dataset(T2M)
ds_p = xr.open_dataset(MSL)

print("=== T2M ===")
print(ds_t)
print("\n=== MSL ===")
print(ds_p)
=== T2M ===
<xarray.Dataset> Size: 75MB
Dimensions:     (valid_time: 1464, latitude: 69, longitude: 185)
Coordinates:
  * valid_time  (valid_time) datetime64[ns] 12kB 2024-01-01 ... 2024-12-31T18...
    expver      (valid_time) <U4 23kB ...
  * latitude    (latitude) float64 552B 6.0 5.75 5.5 5.25 ... -10.5 -10.75 -11.0
  * longitude   (longitude) float64 1kB 95.0 95.25 95.5 ... 140.5 140.8 141.0
    number      int64 8B ...
Data variables:
    t2m         (valid_time, latitude, longitude) float32 75MB ...
Attributes:
    GRIB_centre:             ecmf
    GRIB_centreDescription:  European Centre for Medium-Range Weather Forecasts
    GRIB_subCentre:          0
    Conventions:             CF-1.7
    institution:             European Centre for Medium-Range Weather Forecasts
    history:                 2026-05-10T03:51 GRIB to CDM+CF via cfgrib-0.9.1...

=== MSL ===
<xarray.Dataset> Size: 75MB
Dimensions:     (valid_time: 1464, latitude: 69, longitude: 185)
Coordinates:
  * valid_time  (valid_time) datetime64[ns] 12kB 2024-01-01 ... 2024-12-31T18...
    expver      (valid_time) <U4 23kB ...
  * latitude    (latitude) float64 552B 6.0 5.75 5.5 5.25 ... -10.5 -10.75 -11.0
  * longitude   (longitude) float64 1kB 95.0 95.25 95.5 ... 140.5 140.8 141.0
    number      int64 8B ...
Data variables:
    msl         (valid_time, latitude, longitude) float32 75MB ...
Attributes:
    GRIB_centre:             ecmf
    GRIB_centreDescription:  European Centre for Medium-Range Weather Forecasts
    GRIB_subCentre:          0
    Conventions:             CF-1.7
    institution:             European Centre for Medium-Range Weather Forecasts
    history:                 2026-05-10T03:54 GRIB to CDM+CF via cfgrib-0.9.1...

Dari output di atas kita lihat dimensi dataset: valid_time (atau time) untuk sumbu waktu, latitude menurun dari 6 N ke 11 S, dan longitude dari 95 E ke 141 E. Variabel t2m dalam K dan msl dalam Pa — keduanya perlu dikonversi sebelum masuk ke formula \(w_s\).

Menerapkan Formula w_s ke Satu Timestep ERA5

Kita pilih timestep 2024-01-01 06 UTC sebagai contoh. Konversi unit: suhu K → °C dengan mengurangi 273,15; tekanan Pa → hPa dengan membagi 100. Setelah itu formula Bolton dan \(w_s\) diterapkan element-wise ke seluruh grid.

import numpy as np

# --- pilih satu timestep ---
t_step = "2024-01-01T06:00"

# ambil slice dan konversi unit
t_c   = ds_t["t2m"].sel(valid_time=t_step, method="nearest") - 273.15   # K -> °C
p_hpa = ds_p["msl"].sel(valid_time=t_step, method="nearest") / 100.0    # Pa -> hPa

# hitung w_s menggunakan fungsi dari snippet-1
ws = saturation_mixing_ratio(t_c.values, p_hpa.values)   # g/kg, numpy array

print(f"Timestep  : {t_step}")
print(f"Grid shape: {ws.shape}")
print(f"\nSpatial stats w_s (g/kg) — Indonesia 2024-01-01 06Z:")
print(f"  min  = {ws.min():.2f}")
print(f"  mean = {ws.mean():.2f}")
print(f"  max  = {ws.max():.2f}")

# --- titik stasiun ---
stations = {
    "Jakarta":   (-6.2,  106.8),
    "Pontianak": ( 0.0,  109.3),
    "Makassar":  (-5.1,  119.4),
}
print(f"\n{'Stasiun':<12}  {'T (°C)':>7}  {'p (hPa)':>8}  {'w_s (g/kg)':>11}")
print("-" * 44)
for name, (lat, lon) in stations.items():
    t_pt = float(ds_t["t2m"].sel(valid_time=t_step, latitude=lat,  longitude=lon, method="nearest") - 273.15)
    p_pt = float(ds_p["msl"].sel(valid_time=t_step, latitude=lat,  longitude=lon, method="nearest") / 100.0)
    ws_pt = saturation_mixing_ratio(t_pt, p_pt)
    print(f"{name:<12}  {t_pt:7.2f}  {p_pt:8.2f}  {ws_pt:11.3f}")
Timestep  : 2024-01-01T06:00
Grid shape: (69, 185)

Spatial stats w_s (g/kg) — Indonesia 2024-01-01 06Z:
  min  = 10.21
  mean = 24.90
  max  = 40.27

Stasiun        T (°C)   p (hPa)   w_s (g/kg)
--------------------------------------------
Jakarta         30.63   1009.94       28.352
Pontianak       30.18   1009.30       27.611
Makassar        28.06   1008.71       24.308

Output menunjukkan distribusi \(w_s\) di seluruh Indonesia pada timestep ini. Jakarta dan Pontianak — dua titik pesisir dengan \(T_{2m} > 30\ \text{°C}\) — memiliki \(w_s\) jauh di atas rata-rata grid (\(28{,}35\) dan \(27{,}61\) g/kg vs rata-rata \(24{,}90\) g/kg), sementara Makassar dengan suhu yang sedikit lebih rendah berada tepat di sekitar rata-rata (\(24{,}31\) g/kg). Pola ini mencerminkan kontrol Clausius–Clapeyron yang dominan: kenaikan beberapa derajat di \(T_{2m}\) langsung mendorong \(e_s\) naik secara eksponensial.

Visualisasi w_s sebagai Peta Kontur di Indonesia

Peta kontur memberikan gambaran spasial yang jauh lebih intuitif dari statistik numerik saja. Kita render \(w_s\) pada timestep yang sama dengan cartopy dan simpan hasilnya sebagai PNG.

import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import numpy as np

# --- data: reuse t_c dan p_hpa dari snippet-3 ---
t_step = "2024-01-01T06:00"
t_c_xr   = ds_t["t2m"].sel(valid_time=t_step, method="nearest") - 273.15
p_hpa_xr = ds_p["msl"].sel(valid_time=t_step, method="nearest") / 100.0

lats = t_c_xr["latitude"].values
lons = t_c_xr["longitude"].values
ws_grid = saturation_mixing_ratio(t_c_xr.values, p_hpa_xr.values)

# --- plot ---
fig, ax = plt.subplots(
    figsize=(12, 5),
    subplot_kw={"projection": ccrs.PlateCarree()},
)
ax.set_extent([95, 141, -11, 6], crs=ccrs.PlateCarree())
ax.add_feature(cfeature.COASTLINE, linewidth=0.7, edgecolor="black")
ax.add_feature(cfeature.BORDERS,   linewidth=0.4, edgecolor="gray")
ax.add_feature(cfeature.LAND,      facecolor="lightyellow", zorder=0)
ax.add_feature(cfeature.OCEAN,     facecolor="aliceblue",   zorder=0)

llon, llat = np.meshgrid(lons, lats)
cf = ax.contourf(
    llon, llat, ws_grid,
    levels=np.arange(14, 32, 1),
    cmap="YlGnBu",
    transform=ccrs.PlateCarree(),
    extend="both",
)
cbar = fig.colorbar(cf, ax=ax, orientation="vertical", pad=0.02, shrink=0.85)
cbar.set_label("$w_s$ (g/kg)", fontsize=11)

gl = ax.gridlines(draw_labels=True, linewidth=0.4, color="gray", alpha=0.6, linestyle="--")
gl.top_labels = False
gl.right_labels = False

ax.set_title(
    "Saturation Mixing Ratio ($w_s$) — Indonesia\n2024-01-01 06 UTC (ERA5)",
    fontsize=12,
    pad=10,
)

plt.tight_layout()
plt.savefig("snippet-4.png", dpi=120, bbox_inches="tight")
print("Saved: snippet-4.png")

snippet-4

Peta menunjukkan pola yang konsisten dengan distribusi SST Indonesia: nilai \(w_s\) tertinggi terkonsentrasi di wilayah maritim warm-pool (Laut Banda, Laut Sulawesi, perairan sekitar Maluku) sementara area pegunungan di Kalimantan dan Papua cenderung memiliki nilai lebih rendah karena tekanan permukaan yang lebih kecil di elevasi tinggi.

Variabilitas Musiman DJF dan JJA

\(w_s\) tidak konstan sepanjang tahun. Monsun boreal winter (DJF) membawa SST lebih tinggi dan angin baratan basah ke Indonesia bagian barat, sementara JJA ditandai angin tenggara kering dari Australia yang menekan kelembapan di Indonesia Timur. Mari kita kuantifikasi perbedaan ini dari data ERA5 2024.

import numpy as np

# Hitung w_s untuk seluruh time series
t_c_all   = ds_t["t2m"] - 273.15        # seluruh 2024, K -> °C
p_hpa_all = ds_p["msl"] / 100.0         # Pa -> hPa

# Vectorized — numpy broadcasting
ws_all = saturation_mixing_ratio(t_c_all.values, p_hpa_all.values)   # shape (time, lat, lon)

# Pasang kembali ke DataArray agar bisa groupby season
import xarray as xr
time_coord = ds_t["valid_time"] if "valid_time" in ds_t.coords else ds_t["time"]
ws_da = xr.DataArray(
    ws_all,
    dims=t_c_all.dims,
    coords=t_c_all.coords,
    name="ws",
    attrs={"units": "g/kg", "long_name": "Saturation mixing ratio"},
)

# Seasonal means
ws_seasonal = ws_da.groupby("valid_time.season").mean(dim="valid_time")

print(f"{'Season':<8}  {'mean (g/kg)':>12}  {'std':>7}  {'min':>7}  {'max':>7}")
print("-" * 48)
for season in ["DJF", "JJA"]:
    arr = ws_seasonal.sel(season=season).values
    print(
        f"{season:<8}  {arr.mean():12.3f}  {arr.std():7.3f}  {arr.min():7.3f}  {arr.max():7.3f}"
    )
Season     mean (g/kg)      std      min      max
------------------------------------------------
DJF             23.375    1.577    9.233   26.705
JJA             22.746    1.768    8.770   26.247

Output ringkasan di atas menggambarkan selisih rata-rata \(w_s\) antara DJF dan JJA di atas domain Indonesia. DJF biasanya menghasilkan nilai rata-rata yang lebih tinggi karena SST hangat dan konvergensi monsun aktif; JJA menunjukkan penurunan, terutama di Indonesia Timur akibat pengaruh massa udara dari kontinen Australia yang lebih kering.

Langkah Selanjutnya

Tutorial ini menutup loop dasar dari data ERA5 mentah hingga peta \(w_s\) yang siap diinterpretasi. Beberapa arah lanjutan yang bisa dieksplorasi:

  • Bandingkan \(w_s\) dengan \(q\) ERA5: ERA5 menyimpan specific humidity \(q\) asli (parameter ID 133, satuan kg/kg). Hitung \(w = q / (1 - q)\) dan lihat seberapa besar deviasinya dari \(w_s\) untuk sub-saturated air.
  • Korelasikan dengan IMERG precipitation extremes: Ambil data IMERG 2024 dan uji apakah event presipitasi ekstrem (\(> 50\) mm/6-jam) secara konsisten bertepatan dengan \(w / w_s > 0{,}85\).
  • Profil vertikal \(w_s\): Gunakan ERA5 pressure-level data (temperature dan geopotential di 850 dan 500 hPa, tersedia di era5_t_pl500-850_indonesia_2024_d.nc) untuk menghitung \(w_s\) pada dua level dan lihat bagaimana kapasitas uap air menurun dengan ketinggian.

Eksplorasi artikel meteorologi lainnya di meteo.my.id — kunjungi https://meteo.my.id untuk arsip lengkap.

Referensi

Tidak ada komentar:

Posting Komentar