Photon Point-Kernel Example

This example evaluates Cs-137 PKI doses at arbitrary 3-D coordinates without being limited to the G-14 hotspot selector (for short term scenarios).

For

paths at a tall stack, see the runnable compare_submersion_models.py script. For a visual vertical plume slice with dose curves, see plot_pki_plume_slices.py. The complete list of reproducible diagnostics is in the plots README.

Arbitrary Short-Term Points

The short-term runner pki_nuclides option is intentionally a hotspot comparison. For arbitrary points, use ScenarioEvaluator.evaluate_short_term_pki instead. It reuses the source and orography already configured for the evaluator, does not run hotspot selection, and evaluates the requested weather category directly:

import numpy as np
from astropy.units import Bq, m, s

from ensi_g14_calculator.averaging.weather_averager import SourceParams
from ensi_g14_calculator.dose import PointKernelConfig, PointKernelIntegrator
from ensi_g14_calculator.orography import Orography
from ensi_g14_calculator.short_term_scenarios.runner import ScenarioEvaluator


source = SourceParams(
    H_a=20 * m,
    H_G=15 * m,
    D_i=1.2 * m,
    D_e=1.2 * m,
    W_0=5.0 * m / s,
)
orography = Orography(flat=True)
evaluator = ScenarioEvaluator(source, orography)

points_x = np.array([300.0, 700.0, 1200.0]) * m
points_y = np.array([0.0, 100.0, -100.0]) * m
Q = 1.0e6 * Bq

point_results = evaluator.evaluate_short_term_pki(
    nuclides_bq={"Cs-137": Q},
    x=points_x,
    y=points_y,
    weather_category="F",
    wind_direction=225.0,
    pki_config=PointKernelConfig(  # not strictly necessary, defaults work well. 
        radial_order=12,
        polar_order=12,
        azimuthal_order=16,
    ),
)

point_results is a DoseCollection containing only non-reference pki_photon point records. Query it by coordinate index, interval, nuclide, and age group:

adult_0_8h = point_results.get_point_doses(
    weather="F",
    nuclide="Cs-137",
    interval="0-8h",
    age_group="adult",
    model="pki_photon",
)

for record in adult_0_8h:
    print(record.grid_index, record.x, record.y, record.z, record.dose)

The result has one record for each requested coordinate, nuclide, interval, and age group. Coordinates may be scalar or broadcastable arrays. z is optional; when omitted, observers are placed at local terrain height. Explicit z values are global coordinates relative to the chimney base. No min_radius_direct or min_radius_residence values are involved.

The connector uses the standard short-term release split: Q / 3 for 0-8h, 2 * Q / 3 for 8-24h, and the G-14 half-dispersion factor for the latter interval. Set wind_speed to override the evaluator's default 1 m/s. By default the parent decays inside the plume, every cloud point at its own travel time \(x/U\) (ScenarioEvaluator.decay_in_plume = True). Pass T_fz explicitly to use a fixed flight time instead, or set decay_in_plume = False for the G-14 flight time (and then optionally calculate_time_of_flight = True for \(r/U\)).

For a custom weather distribution or a custom dose calculation, the lower- level PointKernelIntegrator and WeatherAverager interfaces remain available:

import pandas as pd

from ensi_g14_calculator.averaging.weather_averager import WeatherAverager, WeatherStats
from ensi_g14_calculator.dose import PointKernelIntegrator
from ensi_g14_calculator.dose.short_term_dose import ShortTermDoseCalculator

weather = WeatherStats.from_dataframe(
    pd.DataFrame(
        {
            "weather_category": ["F"],
            "wind_speed": [2.0],
            "wind_direction": [225.0],
            "probability": [1.0],
        }
    ),
    randomize_wind_direction_within_bins=False,
)
averager = WeatherAverager(orography=orography, source_params=source)
averager.correct_wind_speed_depending_on_height = False
integrator = PointKernelIntegrator()
dose_calculator = ShortTermDoseCalculator()

def dispersion_at_points(x, y, z):
    return averager.average_dispersion_at_points(x, y, z, weather)

coefficient = integrator.dose_coefficient(
    "Cs-137", (300.0 * m, 0.0 * m, 0.0 * m), dispersion_at_points
)
point_dose = dose_calculator.submersion_pki(
    "Cs-137", Q / 3.0, coefficient, k_s=1.0
)

The point-kernel coefficient uses the packaged air-attenuation (\(\exp{-\mu s}\)), buildup (\(B(E)\)), and ISO fluence-to-effective-dose (\(K(E)\)) data and is truncated at the configured mean-free-path limit (default 5).

Long-Term Points

The long-term runner automatically evaluates PKI at the arbitrary exposed-group coordinates supplied to LongTermEvaluator.run. To use a different set of points without treating them as G-14 exposed groups, use the same pattern with long-term weather and LongTermDoseCalculator.submersion_pki:

from ensi_g14_calculator.dose.long_term_dose import LongTermDoseCalculator

long_term_averager = WeatherAverager(orography=orography, source_params=source)


def long_term_dispersion_at_points(x, y, z):
    return long_term_averager.average_dispersion_at_points(x, y, z, weather)


long_term_coefficient = integrator.dose_coefficient(
    "Cs-137", (points_x[0], points_y[0], 0.0 * m), long_term_dispersion_at_points
)

long_term_dose = LongTermDoseCalculator().submersion_pki(
    "Cs-137", Q, long_term_coefficient
)

The long-term PKI result is also integrated Sv, and the coefficient is in Sv/Bq. For a weather distribution, pass all 3-parameter weather bins to WeatherStats and the same average_dispersion_at_points method will apply their probability-weighted physical dispersion.