advection.advect_detect

Functions for detecting horizontal and vertical advection.

advection.advect_detect.apply_advection_correction(main_data, H_adv, V_adv, HA_Q=None, rn_min=75.0)

Fold measured advective fluxes into the surface energy balance.

Energy-balance bookkeeping (Moderow et al. 2021; Wang et al. 2024)

The advection-augmented surface energy balance, written in the Moderow et al. (2021) OUT-positive convention (positive flux = energy out of the control volume), places the advective terms on the turbulent-sum side alongside H and LE – the available-energy side (Rn - G) is left untouched:

Rn - G = H + LE + HA_T + HA_Q + VAT

so the corrected turbulent + advective sum is:

(H + LE)_corrected = H + LE + HA_T + HA_Q + VAT

Each advective term carries its own sign from the upstream computation (compute_advection_fluxes()): a negative term is energy advected INTO the field (the oasis fingerprint), a positive term is energy advected OUT. This function neither re-signs nor re-derives the terms – it simply adds them on the turbulent side. The residual is reported in the CLAUDE.md convention Residual = Rn - G - H - LE (positive = available energy exceeds the turbulent sum), so the correction moves a gated step’s residual toward zero.

Conditional inclusion gate (Wang et al. 2024)

The advective fluxes are folded in only at timesteps where both Wang’s conditions hold (this is what lifted closure from 89 % to 97 % in the alfalfa study):

1. Rn > rn_min (default 75 W/m^2) – sufficient daytime radiative forcing, AND

2. (H + LE) < (Rn - G) – the (spectrally-corrected) turbulent sum is below the available energy, i.e. there is an under-closure gap for the advected energy to fill.

Where the gate fails (night, low Rn, or already-over-closed steps) the timestep is left exactly uncorrected: (H + LE)_corrected == H + LE. The boolean 'included' mask records which timesteps were corrected.

NaN handling

A NaN in any of Rn/G/H/LE makes the gate comparisons False for that step, so it is left uncorrected (never silently forced closed). A NaN in an advective term contributes 0 to the corrected sum (an unmeasured advective component adds nothing) rather than poisoning the whole corrected sum with NaN.

param main_data:

Main-tower energy balance. Required keys 'H', 'LE', 'Rn', 'G' [all W/m^2], each a length-n series.

type main_data:

dict

param H_adv:

Horizontal heat advection HA_T [W/m^2] (Wang Eq. 5a; the 'H_adv' / 'HA_T' output of compute_advection_fluxes()).

type H_adv:

array-like

param V_adv:

Vertical heat advection VAT [W/m^2] (Wang Eq. 6; the 'V_adv' / 'VAT' output of compute_advection_fluxes()). None (e.g. a horizontal-only call) is treated as all-zero (no vertical advection).

type V_adv:

array-like or None

param HA_Q:

Horizontal moisture advection HA_Q [W/m^2] (Wang Eq. 5b). None (default) is treated as all-zero, preserving the legacy three-argument call signature.

type HA_Q:

array-like or None, optional

param rn_min:

Net-radiation threshold [W/m^2] for Wang’s gate condition (1). Default 75 W/m^2.

type rn_min:

float, optional

returns:

Keys:

• 'Rn', 'G', 'H', 'LE' : the input series (float arrays).

• 'HA_T' / 'H_adv' : horizontal heat advection [W/m^2].

• 'HA_Q' : horizontal moisture advection [W/m^2].

• 'VAT' / 'V_adv' : vertical heat advection [W/m^2].

• 'H_plus_LE_orig' : uncorrected turbulent sum H + LE [W/m^2].

• 'H_plus_LE_corrected' : turbulent sum with the gated advective terms folded in [W/m^2].

• 'available_energy' : Rn - G [W/m^2] (unchanged by the correction; the target the turbulent sum is compared against).

• 'residual_orig' : Rn - G - H - LE [W/m^2] before correction.

• 'residual_corrected' : Rn - G - (H + LE)_corrected [W/m^2] after correction.

• 'included' : boolean mask, True where Wang’s gate passed and the advective terms were applied.

rtype:

dict

References

Moderow et al. (2021), OUT-positive sign convention and advection terms. Wang et al. (2024), conditional-inclusion rule and Eqs. 5a/5b/6.

advection.advect_detect.compute_advection_fluxes(main_data, upwind_data=None, detect_horizontal=None, detect_vertical=None, tower_distance=None)

Compute horizontal and vertical advection flux series for closure.

The horizontal terms are the physical gradient-based advection of Wang et al. (2024), Eqs. 5a/5b (Moderow et al. 2021, Term IV) – not a difference of two sensible-heat fluxes:

HA_T = rho * Cp     * u_bar * (dT/dx) * (zm - h)     # W/m^2  (Eq. 5a)
HA_Q = rho * lambda * u_bar * (dq/dx) * (zm - h)     # W/m^2  (Eq. 5b)

where rho is moist-air density (air_density()), Cp the moist-air specific heat (specific_heat_moist_air()), lambda the latent heat of vaporization (latent_heat_vaporization()) – all evaluated from the main-tower state (P, T, q) – u_bar the mean horizontal wind speed [m/s], and (zm - h) the depth [m] of the air layer between the canopy/lower height h and the measurement height zm.

Vertical mean-advection term (replaces the old residual fallback)

When the caller supplies the planar-fit mean vertical velocity w_bar (or a detect_vertical mask), the vertical heat advection is computed as the measured Lee (1998) / Wang et al. (2024) Eq. 6 term:

VAT = rho * Cp * w_bar * (T_zm - <T>)                # W/m^2  (Eq. 6)

where T_zm is the air temperature at the measurement height (the main tower 'T') and <T> = (1/zm) integral T dz is the column-mean temperature ('T_col' or a 'T_profile'/'z_profile' pair). This is not a closure residual: w_bar must be the planar-fit / tilt-corrected mean vertical velocity, never the raw sonic w (which is tilt-biased). If vertical advection is requested but w_bar or the column-mean temperature is missing, the function raises rather than falling back to the energy-balance residual (CLAUDE.md hard rule: never compute an advection term as the energy-balance residual).

Optionally, when two-level kinematic heat fluxes 'wT_zm' and 'wT_h' [K m/s] are provided, the vertical heat-flux divergence (Wang Eq. 12) is also returned:

VFD_T = -rho * Cp * (wT|zm - wT|h)                   # W/m^2  (Eq. 12)

Diagnostic residual (NOT an advective flux)

The energy-balance closure imbalance residual = (H + LE) - (Rn - G) is still computed and returned, but only as a closure diagnostic – it is not an advective flux and must never be added to the budget as one.

Sign derivation (read carefully)

Let the streamwise coordinate x increase downwind, i.e. point from the upwind tower toward the main tower. The along-wind temperature gradient is approximated by the two-tower finite difference:

dT/dx ~= (T_main - T_upwind) / tower_distance        # K/m

In the oasis case warm dry air sits upwind of a cool transpiring field, so T_upwind > T_main and therefore dT/dx < 0. With u_bar > 0 and (zm - h) > 0 this makes HA_T < 0. Under the Moderow et al. (2021) convention (positive = energy OUT of the control volume), a negative HA_T means sensible-heat energy is advected INTO the field – exactly the expected oasis behaviour. Defining dT/dx as (T_main - T_upwind)/tower_distance therefore yields the correct into-field (negative) sign with no extra negation; warmer upwind air gives HA_T < 0. The moisture gradient uses the identical downwind-positive convention, dq/dx = (q_main - q_upwind)/tower_distance, so drier upwind air (q_upwind < q_main) gives dq/dx > 0 and HA_Q > 0 – the drying signal of dry-air advection.

For the vertical term the canonical oasis has a stable internal boundary layer (cool wet surface, warm air aloft), so the measurement-height temperature exceeds the column mean (T_zm > <T>, i.e. T_zm - <T> > 0) while the planar-fit mean motion is subsidence (w_bar < 0). Then VAT = rho*Cp*w_bar*(T_zm - <T>) < 0: warm air is carried DOWN into the control volume, i.e. heat INTO the field – negative under the OUT-positive convention, the expected oasis vertical-advection sign.

Finally the diagnostic residual is defined here as (H + LE) - (Rn - G) – the negative of the Rn - G - H - LE imbalance written in CLAUDE.md, chosen so a positive value flags net energy advected IN (the legacy adv_in reading). In the oasis the latent flux exceeds the available energy (EF = LE/(Rn-G) > 1, the advective-input fingerprint), so H + LE > Rn - G and residual > 0. It is a closure diagnostic, never an advective flux.

Required inputs (a missing one RAISES, never silently returns zero)

Both main_data and each upwind tower must carry the fields the gradient terms need; tower_distance (or a per-tower 'distance') is required. If any required field is missing a ValueError is raised rather than producing a meaningless zero or a flux difference.

param main_data:

Main-tower series. Required keys:

• 'H' : sensible heat flux [W/m^2]

• 'LE' : latent heat flux [W/m^2]. Open-path LE is assumed already WPL (Webb-Pearman-Leuning 1980) density-corrected — a mandatory separate pre-step, not an advection fix (see advection.advection.wpl_latent_heat_flux()).

• 'Rn' : net radiation [W/m^2]

• 'G' : ground heat flux [W/m^2]

• 'T' : air temperature [°C or K]

• 'q' (or 'RH' [%]) : specific humidity [kg/kg]

• 'u' (or 'wind_speed') : mean horizontal wind speed u_bar [m/s]

• 'zm' : measurement height [m]

• 'h' : canopy / lower height [m] (must satisfy zm > h)

Optional keys: 'P' ambient pressure [Pa] (default 101325), 'wind_dir' [deg from north] for multi-tower fetch selection.

Optional vertical-advection keys (needed only when w_bar or detect_vertical is supplied; missing required ones then raise):

• 'w_bar' (aliases 'w_planar', 'w_mean') : planar-fit / tilt-corrected mean vertical velocity [m/s]. Never the raw sonic w. Required to compute VAT.

• 'T_col' (aliases 'T_column_mean', 'T_mean_col') : column-mean temperature <T> [°C or K]; or 'T_profile' + 'z_profile' (temperatures at heights [m]) from which <T> is approximated by trapezoidal integration. Required to compute VAT.

• 'wT_zm' and 'wT_h' : two-level kinematic heat fluxes [K m/s] at the measurement and canopy heights. When both are present the optional vertical heat-flux divergence VFD_T (Eq. 12) is returned.

'T', 'q'/'RH', 'u', 'zm', 'h', 'P', 'w_bar', 'T_col', 'wT_zm', 'wT_h' may each be a scalar (broadcast) or a length-n series.

type main_data:

dict

param upwind_data:

Upwind tower(s). Each must carry 'T' and 'q' (or 'RH') for the gradients; may carry 'bearing' [deg from north] for fetch selection and a per-tower 'distance' [m] overriding tower_distance. Required – horizontal advection cannot be computed without an upwind reference.

type upwind_data:

dict or list of dicts

param detect_horizontal:

Boolean mask from detect_horizontal_advection(). Where False, HA_T and HA_Q are set to 0 for that step (no advection event). Where None (default), the gradient terms are computed at every step.

type detect_horizontal:

np.ndarray, optional

param detect_vertical:

Boolean mask from detect_vertical_advection(). Supplying it engages the vertical term: VAT is then required (so w_bar and the column-mean temperature must be present, else a ValueError is raised), and where the mask is False VAT is set to 0 for that step (no vertical-advection event).

type detect_vertical:

np.ndarray, optional

param tower_distance:

Separation [m] between the main and upwind tower, used as dx in the gradients. Required unless every upwind dict supplies its own 'distance'.

type tower_distance:

float, optional

returns:

Keys:

• 'HA_T' : horizontal heat advection [W/m^2] (Eq. 5a). Negative = INTO the field (oasis, warm upwind air). length-n array.

• 'HA_Q' : horizontal moisture (latent-energy) advection [W/m^2] (Eq. 5b). length-n array.

• 'H_adv' : alias of 'HA_T' [W/m^2], kept for backward compatibility.

• 'VAT' : measured vertical heat advection rho*Cp*w_bar*(T_zm - <T>) [W/m^2] (Lee 1998; Wang Eq. 6) as a length-n array, or None when no vertical inputs were supplied (horizontal-only call). This is the proper mean-advection term, not a closure residual.

• 'V_adv' : backward-compatible alias of 'VAT'. (Previously a closure-residual estimate; it is now the measured Eq. 6 term.)

• 'VFD_T' : vertical heat-flux divergence -rho*Cp*(wT|zm - wT|h) [W/m^2] (Wang Eq. 12) as a length-n array when both 'wT_zm' and 'wT_h' were supplied; otherwise None.

• 'residual' : diagnostic energy-balance closure imbalance (H + LE) - (Rn - G) [W/m^2], length-n array. This is a closure diagnostic, NOT an advective flux – never add it to the budget as advection.

• 'adv_in' : deprecated alias of 'residual' (the old, misleadingly named key). Identical values.

rtype:

dict

raises ValueError:

If upwind_data is omitted, if tower_distance (and a per-tower 'distance') is missing, if any required gradient field ('T', 'q'/'RH', 'u', 'zm', 'h') is absent, or if the layer depth zm - h is not strictly positive. Additionally, when the vertical term is engaged (w_bar or detect_vertical supplied), if the planar-fit w_bar or the column-mean temperature (T_col / T_profile) is missing – the vertical term is never back-filled from the energy-balance residual.

References

Wang et al. (2024), Eqs. 5a, 5b, 6 and 12. Moderow et al. (2021), Term IV and the OUT-positive sign convention. Lee (1998), vertical mean-advection term.

advection.advect_detect.detect_horizontal_advection(main_flux, *, rn=None, g=None, le_main=None, temp_main=None, temp_upwind=None, humidity_main=None, humidity_upwind=None, wind_dir=None, upwind_dir=None, upwind_flux=None, rn_high=50.0, h_neg_threshold=0.0, ef_tol=1.05, temp_diff_threshold=1.0, humidity_diff_threshold=0.001, wind_sector_deg=45.0, upwind_h_excess=20.0, return_components=False)

Detect periods of horizontal advection influencing the main tower.

Fully vectorized (no Python per-timestep loop) and built from the documented oasis-regime detection signals. Each signal is an independent, parameterized boolean flag; the returned mask is their logical OR. Missing data is carried as np.nan and excluded from every comparison via numpy.isnan() (an x is None test, used by the previous version, never fires on a float ndarray and so silently failed to mask gaps).

Detection signals

1. Negative midday sensible heat / negative Bowen ratio (oasis fingerprint; Wang 2024 §2.2). Fires where H < h_neg_threshold while Rn > rn_high – warm dry air advected onto a cool transpiring surface drives a downward (negative) H in full daytime. Needs main_flux and rn.

2. Evaporative fraction ``EF = LE / (Rn - G) > 1`` (advective input; synthesis §2.3). Fires where LE > (Rn - G) * ef_tol and the available energy Rn - G is positive (EF is only meaningful by day). Needs le_main, rn and g.

3. Horizontal temperature / humidity gradient between the main and an upwind tower (warm and/or dry air upwind). The warm-upwind flag fires where temp_upwind > temp_main + temp_diff_threshold; the dry-upwind flag where humidity_main - humidity_upwind > humidity_diff_threshold (main moister than upwind). Both are evaluated only when the wind is blowing from the upwind tower – i.e. wind_dir is within ±wind_sector_deg of the bearing upwind_dir (a ±45° sector by default; ±20° is a common stricter choice). If either wind_dir or upwind_dir is omitted the sector gate is not applied.

In addition, when an upwind sensible-heat series upwind_flux is supplied, an upwind sensible-heat excess flag fires where upwind_flux > main_flux + upwind_h_excess (also subject to the wind-sector gate). upwind_h_excess is an absolute W/m^2 difference – this is the single, named, documented threshold that replaces the previous version’s contradiction between a docstring that said “20% higher” and code that added an absolute 20 W/m^2. The absolute reading is the one kept.

param main_flux:

Sensible heat flux H at the main tower [W/m^2].

type main_flux:

array-like

param rn:

Net radiation Rn at the main site [W/m^2].

type rn:

array-like, optional

param g:

Soil/ground heat flux G at the main site [W/m^2].

type g:

array-like, optional

param le_main:

Latent heat flux LE at the main tower [W/m^2].

type le_main:

array-like, optional

param temp_main:

Air temperature at the main and upwind towers [°C or K] (same units).

type temp_main:

array-like, optional

param temp_upwind:

Air temperature at the main and upwind towers [°C or K] (same units).

type temp_upwind:

array-like, optional

param humidity_main:

Specific humidity at the main and upwind towers [kg/kg].

type humidity_main:

array-like, optional

param humidity_upwind:

Specific humidity at the main and upwind towers [kg/kg].

type humidity_upwind:

array-like, optional

param wind_dir:

Wind direction at the main tower [deg from north].

type wind_dir:

array-like, optional

param upwind_dir:

Bearing from the main tower toward the upwind tower [deg from north]. Together with wind_dir this gates the gradient and upwind-flux signals to a ±wind_sector_deg sector about the fetch direction.

type upwind_dir:

float, optional

param upwind_flux:

Sensible heat flux H at an upwind reference tower [W/m^2].

type upwind_flux:

array-like, optional

param rn_high:

Net-radiation threshold [W/m^2] above which a period is treated as midday for the negative-H signal.

type rn_high:

float, default 50.0

param h_neg_threshold:

Sensible-heat threshold [W/m^2]; H below this (while Rn > rn_high) flags the oasis negative-H signal.

type h_neg_threshold:

float, default 0.0

param ef_tol:

Evaporative-fraction tolerance. LE > (Rn - G) * ef_tol flags advective input; the 5 % margin guards against measurement noise.

type ef_tol:

float, default 1.05

param temp_diff_threshold:

Minimum upwind-minus-main air-temperature excess [°C or K] for the warm-upwind gradient flag.

type temp_diff_threshold:

float, default 1.0

param humidity_diff_threshold:

Minimum main-minus-upwind specific-humidity excess [kg/kg] for the dry-upwind gradient flag.

type humidity_diff_threshold:

float, default 0.001

param wind_sector_deg:

Half-width [deg] of the wind-direction sector about upwind_dir within which the gradient and upwind-flux signals are evaluated.

type wind_sector_deg:

float, default 45.0

param upwind_h_excess:

Absolute amount [W/m^2] by which upwind_flux must exceed main_flux to flag an upwind sensible-heat excess. Absolute W/m^2, not a percentage.

type upwind_h_excess:

float, default 20.0

param return_components:

If True, also return a dict mapping each criterion name to its boolean mask (plus 'wind_aligned'), suitable for pandas.DataFrame(...).

type return_components:

bool, default False

returns:

• numpy.ndarray – Boolean mask (length n); True where any signal detects advection.

• dict, optional – Returned only when return_components is True, as (mask, components). Keys: 'negative_H', 'ef_gt_1', 'warm_upwind', 'dry_upwind', 'upwind_H_excess' and 'wind_aligned' (the sector gate), each a length-n boolean array.

References

Wang et al. (2024) §2.2 (negative Bowen / oasis fingerprint) and the EF > 1 advective-input synthesis §2.3. Moderow et al. (2021) OUT-positive sign convention.

advection.advect_detect.detect_vertical_advection(vertical_w=None, *, temp_profile_lower=None, temp_profile_upper=None, main_H=None, rn=None, g=None, w_threshold=0.05, rn_g_threshold=50.0, temp_grad_threshold=0.5, h_anomaly_threshold=20.0, use_h_anomaly=True, return_components=False)

Detect periods of energetically significant vertical heat advection.

Centred on the planar-fit mean vertical velocity w_bar (Lee 1998), not on the ad hoc “inverted profile” heuristic the previous version used as a gate. Lee (1998) shows that as |w_bar| approaches ~0.05 m/s the vertical advection term

\[VAT \approx \rho\, C_p\, \bar{w}\, (T_{zm} - \langle T\rangle)\]

becomes energetically significant (of order -100 W/m^2 at midday in the oasis regime; CLAUDE.md “Vertical heat advection”, Wang 2024 Eq. 6). A small mean vertical velocity is therefore, on its own, sufficient grounds to flag a period – it no longer has to co-occur with a temperature inversion.

Fully vectorized (no per-timestep Python loop). Missing data is carried as np.nan and excluded from every comparison with numpy.isnan() (an x is None test, used by the previous version, never fires on a float ndarray and so silently failed to mask gaps). The returned mask is the logical OR of the independent signals below.

Warning

vertical_w must be the planar-fit (or otherwise tilt-corrected) mean vertical velocity w_bar. The raw sonic w must not be used: instrument tilt relative to the mean streamline biases w by far more than the ~0.05 m/s threshold here, so an uncorrected w makes this detector fire on tilt, not on real subsidence/uplift.

Detection signals

1. Primary – significant mean vertical velocity (Lee 1998). Fires where |w_bar| > w_threshold (default 0.05 m/s). This is the dominant signal and the basis of the rewrite. Needs vertical_w.

2. Supporting – daytime vertical temperature gradient of the advective sign. Fires where the period is daytime (Rn - G > rn_g_threshold) and temperature increases with height (temp_profile_upper > temp_profile_lower + temp_grad_threshold). Warm air aloft over a cooler surface gives (T_zm - <T>) > 0; paired with the typical oasis subsidence (w_bar < 0) this yields a negative (downward, energy-IN) VAT – the oasis advective sign under the Moderow 2021 OUT-positive convention. Needs the two temperatures plus rn and g.

3. Optional, weak/secondary – anomalous daytime sensible heat. Fires where it is daytime and main_H < h_anomaly_threshold (very low or downward H when it would normally be positive). This is a screening heuristic only: anomalously low daytime H has many possible causes (cloud, advection, instrument issues), so it is corroborating evidence at best, never a quantitative VAT. Disable it with use_h_anomaly=False.

param vertical_w:

Planar-fit mean vertical velocity w_bar [m/s] – see the warning above. Scalars are broadcast. The primary detection signal.

type vertical_w:

array-like, optional

param temp_profile_lower:

Temperature near the surface/canopy [°C or K].

type temp_profile_lower:

array-like, optional

param temp_profile_upper:

Temperature at (or above) the measurement height [°C or K], same units as temp_profile_lower.

type temp_profile_upper:

array-like, optional

param main_H:

Sensible heat flux H at the main tower [W/m^2] (H-anomaly signal).

type main_H:

array-like, optional

param rn:

Net radiation Rn [W/m^2] (daytime gate).

type rn:

array-like, optional

param g:

Soil/ground heat flux G [W/m^2] (daytime gate).

type g:

array-like, optional

param w_threshold:

Magnitude of w_bar [m/s] above which vertical advection is treated as energetically significant (Lee 1998).

type w_threshold:

float, default 0.05

param rn_g_threshold:

Available-energy threshold Rn - G [W/m^2] separating daytime from nighttime/low-energy periods for the supporting and H-anomaly signals.

type rn_g_threshold:

float, default 50.0

param temp_grad_threshold:

Minimum upper-minus-lower temperature difference [°C or K] counted as a gradient of the advective sign (warm air aloft).

type temp_grad_threshold:

float, default 0.5

param h_anomaly_threshold:

Daytime H [W/m^2] below which the weak H-anomaly signal fires.

type h_anomaly_threshold:

float, default 20.0

param use_h_anomaly:

Whether to include the weak/secondary H-anomaly signal in the mask.

type use_h_anomaly:

bool, default True

param return_components:

If True also return the per-signal boolean arrays (see Returns).

type return_components:

bool, default False

returns:

• numpy.ndarray – Boolean mask of detected vertical-advection periods (length matches the input series).

• tuple of (numpy.ndarray, dict) – If return_components is True, (mask, components) with keys 'w_significant', 'temp_gradient', 'supporting', 'h_anomaly' and 'daytime' – each a length-n boolean array.

References

Lee, X. (1998), on the significance of planar-fit mean vertical velocity for vertical advection. Wang et al. (2024) Eq. 6 (VAT). Moderow et al. (2021) OUT-positive sign convention.