commands¶

compare-weather¶

eap¶

This module implements the process described in VDI 3783 Part 16 [VDI3783p16] to find a substitute anemometer position (EAP) when wind measurements provided as input to the dispersion model AUSTAL [AST31] are not taken by an anemometer inside the AUSTAL model domain (i.e., on a nearby weather station or taken from a weather model).

The position is referred to as “EAP” since in German it is called “Ersatz-Anemometer-Position” (substitute anemometer position).

The module provides three approaches for calculating reference wind profiles:

General approach (calc_ref_geostrophic()): Uses geostrophic wind speeds from VDI 3783 Part 16 Table 1 prescribed at inversion height.
Adapted approach (calc_ref_adapted()): Uses frequency-weighted mean wind speeds from the meteorological time series at the effective anemometer height. This matches the approach used by AUSTAL/TALdia for wind library generation.
Austal approach (austal_ref()): Runs the AUSTAL model that creates a wind library for the current configuration but with terrain removed adn retrieves the reference profiles from this library

Approaches 1 & 2 use the two-layer wind profile model from VDI 3783 Part 8 [VDI3783p8] with Monin-Obukhov similarity theory in the surface layer and an Ekman spiral solution in the upper layer.

austaltools.eap.add_options(subparsers)¶

austaltools.eap.austal_ref(workdir, levels, dirs, tmproot=None, overwrite=False)¶

Generate reference wind profiles using AUSTAL/TALdia.

Creates reference profiles by running AUSTAL on flat terrain and extracting the wind profile at the model origin.

Parameters:

workdir (str or path-like) – Path to the working directory containing austal.txt.
levels (list of float) – Heights (m) to interpolate the wind profile to.
dirs (list of float) – Wind directions (degrees) for which to generate profiles.
tmproot (str or path-like or None, optional) – Directory for temporary files. Default is None.
overwrite (bool, optional) – Whether to overwrite existing reference file. Default is False.

Returns:

A tuple containing:

u_ref (numpy.ndarray) – Eastward reference wind components, shape (len(levels), N_CLASS, len(dirs)).
v_ref (numpy.ndarray) – Northward reference wind components, shape (len(levels), N_CLASS, len(dirs)).

Return type:

tuple(numpy.ndarray, numpy.ndarray)

austaltools.eap.calc_all_eap(g, mx_lvl=None)¶

Find the substitute anemometer position (EAP) for all vertical levels.

Parameters:

g (numpy.ndarray) – 3D array of quality measure values, shape (nx, ny, nz).
mx_lvl (int or None, optional) – Maximum level index to process. If None, all levels are processed. Default is None.

Returns:

A tuple containing:

eap_levels (list of list) – List containing, for each level, a list of EAP candidate coordinates (i, j).
g_upper_levels (list of list) – List containing, for each level, a list of summed quality measures G for each contiguous region.

Return type:

tuple(list, list)

Note

For levels beyond mx_lvl, empty lists are returned.

See also

find_eap()

Example:

>>> g = np.array([[[0.5, 0.8, 0.3],
...                [0.2, 0.7, 0.9],
...                [0.4, 0.6, 0.1]],
...               [[0.3, 0.6, 0.4],
...                [0.1, 0.5, 0.7],
...                [0.2, 0.8, 0.3]]])
>>> eap_levels, g_upper_levels = calc_all_eap(g, mx_lvl=1)
>>> print(eap_levels)
[[[(1, 2), (1, 1), (0, 1)], [(1, 2), (1, 1), (0, 1)]]]
>>> print(g_upper_levels)
[[2.4, 1.6, 1.3], [1.6, 1.3, 1.1]]

austaltools.eap.calc_quality_measure(u_grid, v_grid, u_ref, v_ref, nedge=3, minff=0.5, maxlev=-1)¶

Calculate the quality measure g according to VDI 3783 Part 16.

Compares an AUSTAL wind library to a reference profile and calculates quality criteria for wind direction (gd) and wind speed (gf), which are combined into an overall quality measure g = gd * gf.

Parameters:

u_grid (numpy.ndarray) – Eastward wind component from the wind library. Shape: (nx, ny, nz, nstab, ndir).
v_grid (numpy.ndarray) – Northward wind component from the wind library. Shape: (nx, ny, nz, nstab, ndir).
u_ref (numpy.ndarray) – Eastward reference wind component. Shape: (nz, nstab, ndir).
v_ref (numpy.ndarray) – Northward reference wind component. Shape: (nz, nstab, ndir).
nedge (int, optional) – Number of edge nodes to exclude along each boundary. Default is N_EGDE_NODES.
minff (float, optional) – Minimum wind speed threshold (m/s). Grid points with wind speed below this value are excluded. Default is MIN_FF.
maxlev (int, optional) – Maximum level index to evaluate. Negative values mean all levels are evaluated. Default is -1.

Returns:

A tuple containing:

g (numpy.ndarray) – Overall quality measure, shape (nx, ny, nz). Values in [0, 1], where 1 indicates perfect agreement.
gd (numpy.ndarray) – Quality measure for wind direction, shape (nx, ny, nz).
gf (numpy.ndarray) – Quality measure for wind speed, shape (nx, ny, nz).

Return type:

tuple(numpy.ndarray, numpy.ndarray, numpy.ndarray)

Raises:

ValueError – If grid shapes do not match or reference profile dimensions are incompatible with the wind library.

Note

The algorithm follows VDI 3783 Part 16, Section 6.1 [VDI3783p16] :

Exclude edge nodes (nedge from each boundary).
Reject points where wind doesn’t rotate consistently or wind speed is below minff.
Calculate correlation-based direction criterion gd.
Calculate speed ratio criterion gf.
Combine: g = gd * gf.

See also

find_eap()

austaltools.eap.calc_ref_adapted(levels: list[float], dirs: list[float], working_dir: str | None = None, z0: float | None = None, overwrite: bool = False) → tuple[ndarray, ndarray]¶

Calculate reference wind profiles adapted to the meteorological data.

Uses the two-layer wind profile model from VDI 3783 Part 8 with frequency-weighted mean wind speeds from the meteorological time series prescribed at the effective anemometer height for each stability class.

This approach matches what AUSTAL/TALdia uses internally when generating wind libraries from dispersion class statistics (AKS).

Parameters:

levels (list of float) – Heights above ground (m) for the output profile.
dirs (list of float) – Wind directions (degrees, meteorological convention) for which to generate profiles.
working_dir (str or path-like or None, optional) – Working directory containing austal.txt and the time series file. Default is current directory.
z0 (float or None, optional) – Roughness length (m). Default is VDI_DEFAULT_ROUGHNESS.
overwrite (bool, optional) – Whether to overwrite existing reference file. Default is False.

Returns:

A tuple containing:

u_ref (numpy.ndarray) – Eastward reference wind components, shape (len(levels), N_CLASS, len(dirs)).
v_ref (numpy.ndarray) – Northward reference wind components, shape (len(levels), N_CLASS, len(dirs)).

Return type:

tuple(numpy.ndarray, numpy.ndarray)

Raises:

FileNotFoundError – If the time series file specified in austal.txt is not found.
ValueError – If no time series file is defined in austal.txt or all wind data are invalid.

Note

This is the “adapted” approach that uses:

The effective anemometer height \(h_a\) from the time series file header (dependent on roughness length class).
The frequency-weighted mean wind speed for each stability class from the time series data.

For stability classes with no data, the wind speed is estimated by scaling the VDI geostrophic wind values proportionally.

austaltools.eap.calc_ref_geostrophic(levels: list[float], dirs: list[float], z0: float | None = None, overwrite: bool = False) → tuple[ndarray, ndarray]¶

Calculate reference wind profiles using geostrophic wind values.

Uses the two-layer wind profile model from VDI 3783 Part 8 with geostrophic wind speeds from VDI 3783 Part 16 Table 1 prescribed at the inversion height for each stability class.

Parameters:

levels (list of float) – Heights above ground (m) for the output profile.
dirs (list of float) – Wind directions (degrees, meteorological convention) for which to generate profiles.
z0 (float or None, optional) – Roughness length (m). Default is VDI_DEFAULT_ROUGHNESS.
overwrite (bool, optional) – Whether to overwrite existing reference file. Default is False.

Returns:

A tuple containing:

u_ref (numpy.ndarray) – Eastward reference wind components, shape (len(levels), N_CLASS, len(dirs)).
v_ref (numpy.ndarray) – Northward reference wind components, shape (len(levels), N_CLASS, len(dirs)).

Return type:

tuple(numpy.ndarray, numpy.ndarray)

Note

This is the “general” approach that uses standardized geostrophic wind values independent of the actual meteorological data.

The friction velocity \(u_*\) is calculated iteratively to match the prescribed geostrophic wind at the top of the boundary layer.

austaltools.eap.calc_vdi3783_8(levels: list, dirs: list, z0: float = None, u_a_classes: list[float] | None = None, h_a_classes: list[float] | None = None, overwrite: bool = False)¶

Calculate wind profiles using the VDI 3783 Part 8 two-layer model.

Implements the two-layer boundary layer wind profile model consisting of a Monin-Obukhov surface layer with linear direction turning and an Ekman spiral solution in the upper layer.

Parameters:

levels (array-like) – Heights above ground (m) for the output profile, must be positive and increasing.
dirs (array-like) – Wind directions at reference height (degrees, meteorological convention, 0° = North, 90° = East).
z0 (float or None, optional) – Roughness length (m). Default is VDI_DEFAULT_ROUGHNESS.
u_a_classes (list of float or None, optional) – Reference wind speed (m/s) for each stability class. If None, uses VDI_GEOSTROPIC_WIND (geostrophic wind at inversion height). Default is None.
h_a_classes (list of float or None, optional) – Reference height (m) for each stability class where u_a_classes is prescribed. If None, uses VDI_INVERSION_HEIGHT (inversion height). Default is None.
overwrite (bool, optional) – Whether to overwrite existing output file. Default is False.

Returns:

A tuple containing:

u_ref (numpy.ndarray) – Eastward wind components, shape (len(levels), N_CLASS, len(dirs)).
v_ref (numpy.ndarray) – Northward wind components, shape (len(levels), N_CLASS, len(dirs)).

Return type:

tuple(numpy.ndarray, numpy.ndarray)

Note

The two-layer model from VDI 3783 Part 8 [VDI3783p8] consists of:

Lower layer (\(z \leq h_1\)):

Surface layer following Monin-Obukhov similarity with wind speed:

\[u_1(z) = \frac{u_*}{\kappa} \left[ \ln\frac{z}{z_0} - \psi_m\left(\frac{z}{L}\right) \right]\]

and linear direction turning with gradient \(a = -0.2 A\).

Upper layer (\(z > h_1\)):

Ekman spiral solution with exponentially decaying oscillations:

\[ \begin{align}\begin{aligned}\tilde{u}(z) = u_1(h_1) c_1 + \frac{1}{2A}[(1-c_z)p + s_z q]\\\tilde{v}(z) = u_1(h_1) s_1 + \frac{1}{2A}[(c_z-1)q + s_z p]\end{aligned}\end{align} \]

where \(A = \sqrt{|f_c|/(2K)}\) is the Ekman parameter.

The layer interface height \(h_1\) is calculated from Eq. (A19):

Stable: \(h_1 = \frac{L}{20}\left(\sqrt{1 + \frac{10 h_m}{3\alpha L}} - 1\right)\)
Unstable/neutral: \(h_1 = \frac{h_m}{12\alpha}\)

See also

_calc_h1(), _calc_Km(), _calc_ekman_layer(), _calc_u_star_from_vg()

austaltools.eap.contiguous_areas(array: ndarray) → tuple[ndarray, int]¶

Identify and label contiguous areas in a 2D binary array.

Assigns a unique label to each contiguous region of adjacent True values using 4-connectivity (top, bottom, left, right neighbors).

Parameters:

array (numpy.ndarray) – A 2D boolean array where True represents cells belonging to a contiguous region and False represents background.

Returns:

A tuple containing:

labels (numpy.ndarray) – A 2D integer array of the same shape as array where each contiguous region is labeled with a unique non-negative integer. Background cells are labeled with -1.
num_areas (int) – The number of unique contiguous areas found.

Return type:

tuple(numpy.ndarray, int)

Note

The function uses the union-find algorithm with path compression for efficient region labeling in a two-pass approach.

Example:

>>> arr = np.array([[1, 0, 0], [1, 1, 0], [0, 1, 1]]).astype(bool)
>>> labels, num = contiguous_areas(arr)
>>> print(labels)
[[ 0 -1 -1]
 [ 0  0 -1]
 [-1  0  0]]
>>> print(num)
1

austaltools.eap.find_eap(g_lower: ndarray)¶

Find the substitute anemometer position (EAP) from quality measure.

Identifies the optimal grid point for the substitute anemometer position based on the quality measure g at a single vertical level.

Parameters:

g_lower (numpy.ndarray) – 2D array of quality measure values for each (x, y) grid point.

Returns:

A tuple containing:

eap (list of tuple) – List of EAP candidate coordinates (i, j) in the grid, sorted by decreasing quality (best candidate first).
g_upper (list of float) – List of corresponding summed quality measures G for each contiguous region, sorted in decreasing order.

Return type:

tuple(list, list)

Note

The algorithm follows VDI 3783 Part 16, Section 6.1:

Within each contiguous region of valid points, sum the quality measures to get G.
In the region with the largest G, find the point with the largest individual g.
This point is defined as the EAP.

Example:

>>> g_lower = np.array([[0.5, 0.8, 0.3],
...                     [0.2, 0.7, 0.9],
...                     [0.4, 0.6, 0.1]]).astype(float)
>>> eap, g_upper = find_eap(g_lower)
>>> print(eap[0])  # Best EAP location
(1, 2)

austaltools.eap.interpolate_wind(u_in: list, v_in: list, z_in: list, levels: list)¶

Interpolate wind components to specified heights.

Uses logarithmic interpolation for wind speed and linear interpolation for wind direction.

Parameters:

u_in (list of float) – Eastward wind component values at input heights.
v_in (list of float) – Northward wind component values at input heights.
z_in (list of float) – Heights (m) corresponding to input wind values.
levels (list of float) – Target heights (m) to interpolate to.

Returns:

A tuple containing:

u_out (list of float) – Interpolated eastward wind components.
v_out (list of float) – Interpolated northward wind components.

Return type:

tuple(list, list)

Raises:

ValueError – If u_in, v_in, and z_in do not have the same length.

Example:

>>> u_in = [1.0, 2.0, 3.0]
>>> v_in = [0.5, 1.0, 1.5]
>>> z_in = [10.0, 50.0, 100.0]
>>> levels = [25.0, 75.0]
>>> u_out, v_out = interpolate_wind(u_in, v_in, z_in, levels)

austaltools.eap.main(args)¶

Main entry point for the EAP analysis.

Reads a wind library, calculates reference profiles, computes quality measures, finds optimal EAP locations, and optionally creates plots and updates the AUSTAL configuration.

Parameters:

args (dict) –

Command line arguments dictionary with keys:

working_dir: Path to working directory.
grid: Grid ID to evaluate.
reference: Reference profile method (‘general’, ‘simple’, ‘file’, or ‘austal’).
overwrite: Whether to overwrite existing files.
max_height: Maximum evaluation height.
edge_nodes: Number of edge nodes to exclude.
min_ff: Minimum wind speed threshold.
height: Target height for EAP selection.
report: Whether to print detailed report.
austal: Whether to update austal.txt.
plot: Plot output specification.

See also

add_options()

austaltools.eap.print_report(args: dict, g: ndarray, gd: ndarray, gf: ndarray, eaps: list[list[tuple]], g_upper: list[list[float]], axes: dict[str, list])¶

Print a detailed report of EAP analysis results.

Outputs a formatted report mimicking the style of the VDI 3783 Part 16 reference implementation (TAL-Anemo.zip).

Parameters:

args (dict) – Command line arguments dictionary.
g (numpy.ndarray) – Overall quality measure, shape (nx, ny, nz).
gd (numpy.ndarray) – Direction quality measure, shape (nx, ny, nz).
gf (numpy.ndarray) – Speed quality measure, shape (nx, ny, nz).
eaps (list of list of tuple) – EAP coordinates for each level, as returned by calc_all_eap().
g_upper (list of list of float) – Summed quality measures for each level.
axes (dict) – Dictionary with keys ‘x’, ‘y’, ‘z’ containing grid coordinates.

austaltools.eap.read_ref(file: str, levels: list[float], dirs: list[float], linear_interpolation: bool = False)¶

Read reference wind profiles from file.

Reads wind profiles in the format of Ref1d.dat from the VDI 3783 Part 16 reference implementation and interpolates to the requested heights and directions.

Parameters:

file (str) – Path to the reference profile file.
levels (list of float) – Target heights (m) to interpolate to.
dirs (list of float) – Target wind directions (degrees) to extract.
linear_interpolation (bool, optional) – If True, use linear interpolation for wind speed (for comparison with VDI reference implementation). If False, use logarithmic interpolation. Default is False.

Returns:

A tuple containing:

u_ref (numpy.ndarray) – Eastward reference wind components, shape (len(levels), N_CLASS, len(dirs)).
v_ref (numpy.ndarray) – Northward reference wind components, shape (len(levels), N_CLASS, len(dirs)).

Return type:

tuple(numpy.ndarray, numpy.ndarray)

Raises:

ValueError – If no matching profile is found for a stability class.

See also

write_ref()

austaltools.eap.run_austal(workdir, tmproot=None)¶

Create a reference wind library using AUSTAL/TALdia.

Invokes AUSTAL with the -l parameter to generate a wind library for flat terrain with the anemometer at the model origin.

Parameters:

workdir (str or path-like) – Path to the working directory containing austal.txt.
tmproot (str or path-like or None, optional) – Directory for temporary files. If None, uses workdir. Default is None.

Returns:

A tuple containing:

u_tmp (numpy.ndarray) – Eastward wind component grid.
v_tmp (numpy.ndarray) – Northward wind component grid.
ax_tmp (dict) – Dictionary containing grid axes and metadata.

Return type:

tuple(numpy.ndarray, numpy.ndarray, dict)

Raises:

ValueError – If austal.txt is not found or AUSTAL fails.
OSError – If the AUSTAL executable is not found.

Note

This function creates a temporary directory, modifies the AUSTAL configuration for flat terrain, runs AUSTAL, extracts the results, and cleans up the temporary files.

austaltools.eap.same_sense_rotation(val, ref)¶

Check if wind directions rotate in the same sense.

Determines whether the wind directions in val and ref both rotate in the same direction (both clockwise or both counter-clockwise) as the input wind direction varies.

Parameters:

val (array-like) – Tested wind directions (degrees, meteorological convention).
ref (array-like) – Reference wind directions (degrees, meteorological convention).

Returns:

True if both arrays rotate in the same sense, False otherwise.

Return type:

bool

Note

This function is used in the EAP algorithm to reject grid points where the wind does not rotate consistently with the reference profile as required by VDI 3783 Part 16, Section 6.1, criterion 2 [VDI3783p16] .

Write reference wind profiles to file.

Writes wind profiles in the format of Ref1d.dat from the VDI 3783 Part 16 reference implementation (TAL-Anemo.zip).

Parameters:

file (str) – Output file path.
out_levels (array-like) – Heights (m) to include in output.
out_dirs (array-like) – Wind directions (degrees) to include in output.
u_ref (numpy.ndarray) – Eastward wind components, shape (nz, N_CLASS, ndir).
v_ref (numpy.ndarray) – Northward wind components, shape (nz, N_CLASS, ndir).
axes_ref (tuple) – Tuple of (levels, stability_classes, directions) arrays corresponding to the dimensions of u_ref and v_ref.
overwrite (bool or None, optional) – If True, overwrite existing file. If False, raise FileExistsError. If None, prompt user interactively. Default is None.

Raises:

FileExistsError – If file exists and overwrite is False.

See also

read_ref()

austaltools.eap.AUSTAL_ROUGHNESS = 0.01¶

float: Roughness length \(z_0\) (m) used for reference wind profile calculation.

Corresponds to CORINE class 231 “short grass”, according to VDI 3783 Part 8 [VDI3783p8] .

austaltools.eap.MAX_HEIGHT = 100.0¶: float: Maximum height (m) above ground to which wind data are included in the EAP search algorithm.

austaltools.eap.MIN_FF = 0.5¶: float: Minimum wind speed (m/s) for which wind data are included in the EAP search algorithm.

austaltools.eap.N_CLASS = 6¶: int: Number of Klug-Manier stability classes (I through V, with III split into III1 and III2).

austaltools.eap.N_EGDE_NODES = 3¶: int: Number of model nodes along each side of the model domain that should be excluded to avoid edge effects in the EAP search algorithm.

austaltools.eap.VDI_DEFAULT_ROUGHNESS = 0.1¶

float: Default roughness length \(z_0\) (m) for wind profile calculation.

Value of 0.1 m is used instead of the original VDI value (0.02 m for LBM-DE class 231 “Wiesen und Weiden”) since 2023, according to UBA TEXTE 144/2023 “Weiterentwicklung ausgewählter methodischer Grundlagen der Schornsteinhöhenbestimmung und der Ausbreitungsrechnung nach TA Luft”.

austaltools.eap.VDI_GEOSTROPIC_WIND = [1.6, 2.5, 7.8, 5.6, 4.2, 3.8]¶

list of float: Geostrophic wind speed \(v_g\) (m/s) for each stability class.

Values from VDI 3783 Part 16, Table 1.

Index 0 corresponds to Class I (very stable), Index 5 corresponds to Class V (very unstable).

austaltools.eap.VDI_INVERSION_HEIGHT = [250, 250, 800, 800, 1100, 1100]¶

list of int: Mixing layer / inversion height \(h_m\) (m) for each stability class.

Values from VDI 3783 Part 8 (2002), Table 4.

Index 0 corresponds to Class I (very stable), Index 5 corresponds to Class V (very unstable).

austaltools.eap.VDI_THETA_GRADIENT = [0.008, 0.0057, 0.0032, 0.0012, 0.0003, 0.0]¶

list of float: Potential temperature vertical gradient (K/m) for each stability class.

Values from VDI 3783 Part 16, Table 1.

Index 0 corresponds to Class I (very stable), Index 5 corresponds to Class V (very unstable).

fill_timeseries¶

This module allows to create time-dependent source strenght timeseries as input for simulations with the German regulatory dispersion model AUSTAL [AST31]

austaltools.fill_timeseries.add_options(subparsers)¶

austaltools.fill_timeseries.expand_cycles(yinfo)¶

Processes a dictionary of cycle information, applying templates to cycles as needed.

This function validates the provided yinfo dictionary to ensure it contains the correct data structure, extracts templates, and applies them to the cycles. If a cycle specifies a template, the template is applied, including any emission factor calculations based on specified substances.

Parameters:

yinfo (dict) – A dictionary containing cycle information. The keys represent cycle IDs and the values are dictionaries with specific cycle information, which can include ‘column’, ‘source’, ‘template’, and ‘factors’.

Raises:

ValueError –

If yinfo is not a dictionary, or
if it contains invalid structure such as null at the top level,
missing template definitions for requested cycles,
missing emission factors for specified substances, or
if emission factors are present without a selected substance.

Returns:

A dictionary of processed cycles where each cycle has necessary attributes set such as ‘multiplier’, ‘emissionfactor’, and ‘substance’. The keys are cycle IDs and the values are dictionaries containing the expanded cycle information.

Return type:

dict

Example:

Consider a set of cycle information with one defined template and two cycles:

>>> yinfo = {
...     'template1': {'column': None, 'source': None, 'factors': {'NOX': 1.0}},
...     'cycle1': {'column': '01.nox', 'template': {'name': 'template1', 'substance': 'NOX'}},
...     'cycle2': {'column': '01.xx'},
...     'cycle3': {'column': '02.nox', 'multiplier': 2.5}
... }
>>> expand_cycles(yinfo)
{
    'cycle1': {'column': '01.nox', 'source': None, 'substance': 'NOX', 'emissionfactor': 1.0, 'multiplier': 1.0},
    'cycle2': {'column': '01.xx', 'multiplier': 1.0, 'emissionfactor': 1.0, 'substance': None},
    'cycle3': {'column': '02.nox', 'multiplier': 2.5, 'emissionfactor': 1.0, 'substance': None}}

This example demonstrates how the specified template is applied to cycle1 and cycle2 is processed without a template.

austaltools.fill_timeseries.get_timeseries(file: str, time: DatetimeIndex)¶

Parse yaml file containing cycle(s) information and generate an emission time series.

This funtion is essentially a wrapper that applies for parse_cycle() to a yaml file.

Parameters:

file (str) – filename (optionally containing a path)
time (pandas.Series) – Time series

Returns:

time series of emssions of all emissions descrcibed in file

Return type:

pandas.Dataframe with time as index and column-ids as colums

Example:

>>> yaml_text = '''
... meinname:
...   column: 01.so2
...   start:
...     at:
...       time: 1-11/2
...       unit: month
...     offset:
...       time: 1,3
...       unit: week
...   sequence:
...   - ramp:
...       time: 1
...       unit: day
...       value: 9.0
...   - const:
...       time: 36
...       unit: hour
...       value: 1.1
...'''
>>> with open("cycle.yaml, "w") as f:
>>>     f.write(yaml_text)
>>> time = pandas.date_range("2000-01-01 00:00",
...                          "2000-01-02 00:00", freq="1h")
>>> get_cycle(file, time)
    ('01.so2',
     2000-01-01 00:00:00    0.0
     2000-01-01 01:00:00    0.0
     2000-01-01 02:00:00    0.0
     2000-01-01 03:00:00    0.0
     2000-01-01 04:00:00    0.0
                           ...
     2000-12-24 07:00:00    1.1
     2000-12-24 08:00:00    1.1
     2000-12-24 09:00:00    1.1
     2000-12-24 10:00:00    1.1
     2000-12-24 11:00:00    1.1
     Name: foo, Length: 745, dtype: float64)

Note:

The format of the yaml file is described under variable values

austaltools.fill_timeseries.main(args)¶

Process the data file based on the provided arguments.

Parameters:

args (dict) – Dictionary containing the following keys:
args["action"] – (str) – The action to perform. Possible values are ‘list’, ‘week-5’, ‘week-6’, or ‘cycle’.
args["cycle_file"] – (str) – The name of the cycle file (required for ‘cycle’ action).
args["holiday_month"] – (*list, optional) – List of months (1-12) considered as holidays.
args["holiday_week"] – (*list, optional) – List of weeks (1-52) considered as holidays.
args["hour_begin"] – (int, optional) – The daily start of the working time, i.e. the first hour of each working day the source emits pollutants (evaluated for ‘week-5’ and ‘week-6’ actions). Defaults to DEFAULT__BEGIN.
args["hour_end"] – (int, optional) – The daily end of the working time, i.e. the last hour of each working day the source emits pollutants (evaluated for ‘week-5’ and ‘week-6’ actions). Defaults to DEFAULT_END .
args["column_id"] – (str) – The column ID to process (required for ‘week-5’ and ‘week-6’ actions).
args["output"] – (list) – The source strength (in g/s) when the source is emitting (required for ‘week-5’ and ‘week-6’ actions).
args["working_dir"] – (str) – The path to the directory containing the data file. The datafile is named zeitreihe.dmna or timeseries.dmna, depending on the language setting of the AUSTAL model.

Raises:

ValueError – If the data file is not in DMNA timeseries format.
ValueError – If the action is unknown.
ValueError – If required arguments are missing or invalid.

Note:

the datafile zeitreihe.dmna/timeseries.dmna must be created by invoking AUSTAL with paramter -z

austaltools.fill_timeseries.parse_cycle(c_id: str, c_info: dict, time: DatetimeIndex) → Series¶

Parse cycle information and generate an emission time series.

Parameters:

c_id (str) – Cycle identifier
c_info (dict) –
Cycle information dictionary. Must contain the keys:
- ”column”: str, column identifier (must not be equal to c_id)
- ”start”: dict, must contain: - “at”: str, start time information - “offset” (optional): str, offset time information
- ”sequence” or “list”: list, sequence or list of values
- ”unit” (optional): str, unit information in the format “<mass unit>/<time interval>”
time (pandas.Series) – Time series

Raises:

ValueError –

If required keys are missing or invalid values are found in c_info. Possible errors include:

if time is an invalid type or time series does not have a unique interval
if c_info does not contain the referred column name
if the cycle name c_id is equal to the column name
if c_info has neithert none or both of a cycle or list entry
if c_info has not start entry
if the start entry is not a dict or does not contain an at entry
’sequence’ item contains more or less than one entry or the entry cannot be parsed
c_info['list'] does not contain a list
the unit info in c_info['unit']s cannot be parsed
the mass unit in c_info['unit'] is not al valid weight unit
the time interval in c_info['unit'] is not a valid time unit

Returns:

Column identifier and generated cycle series

Return type:

tuple (str, pandas.Series)

Example:

>>> import pandas as pd
>>> c_id = "foo"
>>> c_info = {'column': '01.so2',
...   'start': {'at': {'time': '1-11/2', 'unit': 'month'},
...   'offset': {'time': '1,3', 'unit': 'week'}},
...   'sequence': [
...     {'ramp': {'time': 1, 'unit': 'day', 'value': 9.0}
...    },
...    {'const': {'time': 36, 'unit': 'hour', 'value': 1.1}}]}
>>> time = pd.date_range("2000-01-01 00:00",
...                          "2000-01-02 00:00", freq="1h")
>>> fill_timeseries.parse_cycle(c_id, c_info, time)
    ('01.so2',
     2000-01-01 00:00:00    0.0
     2000-01-01 01:00:00    0.0
     2000-01-01 02:00:00    0.0
     2000-01-01 03:00:00    0.0
     2000-01-01 04:00:00    0.0
                           ...
     2000-12-24 07:00:00    1.1
     2000-12-24 08:00:00    1.1
     2000-12-24 09:00:00    1.1
     2000-12-24 10:00:00    1.1
     2000-12-24 11:00:00    1.1
     Name: foo, Length: 745, dtype: float64)

austaltools.fill_timeseries.parse_time(info, name='', multi=True)¶

Parse time information from a given dictionary.

The dictionary info must contain the following keys: - ‘time’: A string representing the time information. - ‘unit’: A string representing the unit of time.

Parameters:

info (dict) – Dictionary containing time information.
name (str) – Optional name for the time info, used in error messages.
multi (bool) – Flag indicating whether multiple times are allowed.

Raises:

ValueError – If ‘time’ or ‘unit’ keys are missing in the info dictionary.
ValueError – If multiple times are defined when multi is False.

Returns:

A tuple containing the parsed time count and unit.

Return type:

tuple

austaltools.fill_timeseries.parse_time_unit(string)¶

Parse a string and determine which time unit it describes: - ‘month’, ‘months’, ‘mon’ for months - ‘day’, ‘days’, ‘d’ for days - ‘hour’, ‘hours’, ‘hr’, ‘hrs’, ‘h’ for hours

Parameters:: string (str) – the string to parse
Returns:: the parsed time unit
Return type:: str

austaltools.fill_timeseries.DEFAULT_BEGIN = 8¶: Default staring hour for a workday (first hour during which emsssions are created)

austaltools.fill_timeseries.DEFAULT_END = 17¶: Default end hour for a workday (last hour during which emsssions are created)

heating¶

import_buildings¶

This module provides funtions to processes a GeoJSON file containing building data, extracts the corner points, fit rectangles to corner points, plot the buildings and to write building inforamtion the ‘austal.txt’ configuration file.

austaltools.import_buildings.add_options(subparsers)¶

austaltools.import_buildings.building_corners(build: Building) → list[tuple[float, float]]¶

Return the four corner positions of a rectangle with the properties:

Parameters:: build – Building object defining lower-left corner, rectangle extensions and rotation in degrees counterclockwise from the x-axis.
Returns:: list of corner positions
Return type:: list[tuple[float, float]]

austaltools.import_buildings.building_new()¶: return a new Building object :return: empty building object :rtype: _tools.Building()

austaltools.import_buildings.check_tolerances(tolerance: float, build: Building, points: list[tuple[float, float]]) → bool¶

Check if the given points are within the specified tolerance from the building corners.

This function calculates the minimum distance from each point to the building corners and checks if all distances are within the specified tolerance. It also ensures that all four corners of the building are represented by the closest points.

Parameters:

tolerance (float) – The maximum allowable distance from the points to the building corners.
build (_tools.Building) – The building object containing the corner coordinates.
points (list[tuple[float, float]]) – A list of tuples representing the coordinates of the points to be checked.

Returns:

True if all points are within the tolerance and all corners are represented, False otherwise.

Return type:

bool

Example:

>>> building = _tools.Building(corners=[(0, 0), (0, 10), (10, 0), (10, 10)])
>>> points = [(1, 1), (1, 9), (9, 1), (9, 9)]
>>> check_tolerances(2.0, building, points)
True

austaltools.import_buildings.deduplicate(points, tolerance=None)¶

Removes duplicate points from a list of points.

Parameters:

points (list[tuple[float, float]]:) – A list of (x, y) coordinate tuples. tolerance (float, optional): Maximum distance for considering points as duplicates. If provided, points within this distance are considered duplicates.

Returns:

A list of unique points after removing duplicates.

Return type:

list[tuple[float, float]]:

Note:

If tolerance is specified, points within the tolerance distance are considered duplicates.
The dist_points function (not defined here) calculates the distance between two points.

Example:

>>> points = [(1.0, 2.0), (3.0, 4.0), (1.0, 2.0), (5.0, 6.0)]
>>> deduplicate(points)
[(1.0, 2.0), (3.0, 4.0), (5.0, 6.0)]

austaltools.import_buildings.dist_points(p1: tuple[float, float], p2: tuple[float, float]) → float¶

Calulate distance between two points in a 2D a cartesian coordinate system

Parameters:

p1 (tuple[float, float]) – point 1
p2 (tuple[float, float]) – point 2

Returns:

distance

Return type:

float

austaltools.import_buildings.dist_to_line(a: float, b: float, s: float, p: tuple[float, float]) → float¶

returns distance of point p to line with slope b ant offset a

Parameters:

a (float) – offset
b (float) – slope
s – (rotation) sense (see austaltools.austal_buildings_geojson.line_through())
p (tuple[float, float]) – point

Returns:

distance

Return type:

float

austaltools.import_buildings.extract_polygons(features, origin)¶

Extracts polygons from a list of features.

This function processes a list of features and extracts polygons from them. It checks the type of each feature, validates the geometry, and converts coordinates to a model coordinate system based on the specified origin.

Parameters:

features (list[dict]) – A list of feature dictionaries.
origin (tuple[float, float]) – The origin point for coordinate conversion.

Returns:

A list of polygons, where each polygon is represented by a tuple containing: - Feature index - Polygon index within the feature - List of points (x, y) in the model coordinate system

Return type:

list[tuple[int, int, list[tuple[float, float]]]]

Note:

If a feature has unsupported geometry type, it will be skipped.
For MultiPolygon features, only the exterior ring (first set of coordinates) is considered.
Holes in polygons are ignored.
The logger is used to report errors and warnings.

Example:

>>> features = [
...     {'type': 'Feature', 'geometry': {
...         'type': 'Polygon',
...         'coordinates': [[(0, 0), (1, 0), (1, 1), (0, 1)]]}
...     },
...     {'type': 'Feature', 'geometry': {
...         'type': 'MultiPolygon',
...         'coordinates': [[[(2, 2), (3, 2), (3, 3), (2, 3)]]]}
...     },
... ]
>>> origin = (0, 0)
>>> extract_polygons(features, origin)
[(0, 0, [(0.0, 0.0), (1.0, 0.0), (1.0, 1.0), (0.0, 1.0)]),
 (1, 0, [(2.0, 2.0), (3.0, 2.0), (3.0, 3.0), (2.0, 3.0)])]

austaltools.import_buildings.find_building_around(points: list[tuple[float, float]], tolerance: float) → Building | None¶

Find the minimal rectagle encircling the points. Returns lower left corner as x and y coordinate, the exetensions of the rectangle, width in x-direction and depth in y-direction and its rotation angle in degrees counterclockwise from the x-axis.

Parameters:

points (list[tuple[float, float]]) – list of the points positions
tolerance (float) – minimum distance between points to consider them as different positions

Returns:

Building object defining x, y, width, depth and angle or none if finding fails

Return type:

_tools.Building or None

austaltools.import_buildings.is_rectangle(points, tolerance=0.1)¶

Check if the given points form a rectangle within the specified tolerance.

This function calculates the diagonals of the quadrilateral formed by the points and checks if the difference between the diagonals is within the specified tolerance value.

Parameters:

points (list[tuple[float, float]]) – A list of tuples representing the points.
tolerance (float, optional) – The maximum allowed difference between the diagonals.

Returns:

True if the points form a rectangle within the tolerance, False otherwise.

Return type:

bool

Example:

>>> points = [(0, 0), (0, 2), (2, 0), (2, 2)]
>>> is_rectangle(points)
True

austaltools.import_buildings.line_through(p1: tuple[float, float], p2: tuple[float, float]) -> (<class 'float'>, <class 'float'>)¶

Returns parameters of the line through two points: slope and offset of the linear equation and (rotation) sense: - +1 if p2 is to the right (positive x axis) of p1 - -1 if p2 is to the left (negative x axis) of p1

Parameters:

p1 (tuple[float, float]) – first point
p2 (tuple[float, float]) – second point

Returns:

intercept and slope of the line

Return type:

float, float

Example:

>>> p1 = (1, 2)
>>> p2 = (3, 4)
>>> line_through(p1, p2)
(1.0, 1.0, 1)

austaltools.import_buildings.main(args)¶

Main entry point: extract buildings from a GeoJSON file and write them to the config file ‘austal.txt’.

This function processes a GeoJSON file containing building data, extracts the relevant information, and writes it to a configuration file for further use. The function also supports optional plotting of building shapes.

Parameters:

args (dict) – A dictionary containing the following keys: - ‘zvalue’: (optional) The name of the JSON variable denoting building height. - ‘height’: (optional) A fixed height value for all buildings. - ‘tolerance’: The tolerance value for checking if the points form a rectangle. - ‘wdir’: The working directory where the ‘austal.txt’ file is located. - ‘file’: The name of the GeoJSON file containing building data. - ‘dry_run’: A boolean flag indicating whether to perform a dry run (no file output). - ‘plot’: A boolean flag indicating whether to plot the building shapes.

Raises:

ValueError – If the GeoJSON file is not of type ‘FeatureCollection’ or if the CRS is not ‘EPSG:31463’.
ValueError – If neither GaussKrueger nor UTM coordinates are found in the configuration.
ValueError – If no height information is available for a building.

Example:

>>> args = {
>>>     'zvalue': 'height',
>>>     'height': None,
>>>     'tolerance': 0.1,
>>>     'wdir': '/path/to/working/directory',
>>>     'file': 'buildings.geojson',
>>>     'dry_run': True,
>>>     'plot': False
>>> }
>>> main(args)

austaltools.import_buildings.nearest_point_on_line(a: float, b: float, p: tuple[float, float]) → tuple[float, float]¶

Returns position of the point on the line of slope b ant offset a that is closest to point p.

Parameters:

a (float) – offset
b (float) – slope
p (tuple[float, float]) – point

Returns:

distance

Return type:

tuple[float, float]

austaltools.import_buildings.plot_building_shapes(args: dict, polygons: list[tuple], buildings: list[Building], topo: str = None)¶

Plot buildings and polygon shapes from geojson file

Parameters:

args (dict) – command line arguments
polygons (list[tuple]) – list of tuple (#feature, # polygon, list of points)
buildings (list[_tools.Building]) – Building objects
topo (str (optional)) – Name of topography file (*.grid)

austaltools.import_buildings.rotating_caliper(points: list[tuple[float, float]]) -> (<class 'float'>, <class 'float'>, <class 'float'>, tuple[float, float])¶

Return the equation of the one of all lines through two adjacent points, for which all other points are closest to the line, as well as dististance to and postion of the most distant point

Parameters:: points (tuple[float, float]) – point positions
Returns:: offset and slope of the line, most distant point distance and position of the first base point
Return type:: float, float, float, tuple[float, float]

austaltools.import_buildings.sort_anticlock(points: list[tuple[float, float]]) → list[tuple[float, float]]¶

Sort points anticlockwise around the center point

Parameters:: points (list[tuple[float,float]]) – point positions to sort
Returns:: sorted point positions
Return type:: list[tuple[float,float]]

austaltools.import_buildings.DEFAULT_FILE = 'haeuser.geojson'¶: default name of the geojson value that indicates building height

austaltools.import_buildings.DEFAULT_ZVALUE = 'height'¶: allowed difference between geojson polygon corners and the rectangle fitted to them in m

austaltools.import_buildings.logger = <RootLogger root (WARNING)>¶: default name of the geojson file that contains building data

input_terrain¶

input_weather¶

plot¶

Module containing functions to create a basic plot from austal result data

Plots can be shon interactively if the user operates on a terminal that has an X-server running. For example Linux with a running desktop environment, an Anaconda environment or a ssh connection with active X-forwarding and a local X-server running.

austaltools.plot.add_options(subparsers)¶

austaltools.plot.main(args)¶

This is the main working function

Parameters:

args (dict) – The command line arguments as a dictionary.
args['working_dir'] (str) – The working directory where files are located (i.e. where austal.txt is stored).
args['file'] (str) – The input file name. If it doesn’t have a ‘.dmna’ extension, it will be added.
args['buildings'] (bool) – A flag indicating whether to plot buildings from the configuration.
args['stdvs'] (float) – The standard deviation value to mark (additional) concentrations as significant by overlaying dots..
args['plot'] (str or None) – The plot file name. If None or ‘-’, the plot will be shown interactively. If ‘__default__’, the name of the displayed data file with extension .png will beused.

Raises:

OSError – If the configuration file cannot be found or read.
ValueError – If the data shape is not understood or if the standard deviation shape does not match the data shape.

austaltools.plot.parse_austal_outputname(filename: str)¶

analyze name of austal output file

Parameters:

filename – str

Returns:

information about file contents:

substance: name of pollutant (xx for unknown/not specified)
averaging: duration of averaging interval (accumulation, year, day or hour)
rank: rank of output value in list of all averages of the same length
kind: type of output (load, stdev or index)
grid: number of grid. 0 if not given / no staggered grids.

Return type:

dict

select-year¶

steepness¶

create basic plot for austal result files

austaltools.steepness.add_options(subparsers)¶

austaltools.steepness.main(args)¶

commands¶

compare-weather¶

eap¶

fill_timeseries¶

heating¶

import_buildings¶

input_terrain¶

input_weather¶

plot¶

select-year¶

steepness¶

transform¶

volout¶

windfield¶

windrose¶

Table of Contents

Previous topic

Next topic