dased.helpers.srcloc¶

Source-localization helper utilities for DAS applications.

This subpackage provides forward operators, likelihood helpers and posterior tools used in Bayesian source localization. Notable exported symbols include ForwardHomogeneous, DataLikelihood, and calculate_posterior.

External types commonly referenced are torch.Tensor, xarray.DataArray, and geopandas.GeoDataFrame.

Classes¶

`ForwardBase`	Base class for seismic forward operators.
`ForwardHomogeneous`	Calculates travel times and SNRs assuming a constant velocity medium.
`ForwardLayeredLookup`	Uses precomputed travel times and incidence angles from lookup tables generated by ray tracing through layered earth models.
`MagnitudeRelation`	Magnitude-distance relationship for seismic signal amplitude modeling.
`DataLikelihood`	Computes likelihood of observed data given source model parameters and experimental design.
`SurfaceField_Distribution`	3D spatial distribution that follows surface topography with depth variation.
`PolygonUniform`	Uniform distribution within a Shapely polygon with custom vertical distribution.
`IndependentNormal`	Independent Normal distribution with element-wise independence.

Functions¶

calculate_posterior

Calculate Bayesian posterior distribution over a 2D or 3D grid.

Descriptions¶

class dased.helpers.srcloc.ForwardBase(wave_type, distance_relation, sensor_type='strain', direction_smoothing=0.0, incidence_max=None)¶

Base class for seismic forward operators.

This private base class provides core functionality for calculating travel times and signal-to-noise ratios for seismic waves between sources and sensors. It handles directional sensitivity calculations and coordinate transformations.

Parameters:

wave_type (str) – Type of wave (‘P’, ‘S’, ‘SV’, ‘SH’, ‘R’, ‘Rayleigh’, ‘Love’).
distance_relation (Callable) – Function to calculate distance-based attenuation.
sensor_type (str) – Type of sensor (‘strain’, ‘displacement’, ‘strain1C’, etc.).
direction_smoothing (float | Tuple[float, float]) – Standard deviation(s) (in radians) for directional sensitivity smoothing. If float, the same value is used for theta and phi smoothing. If tuple (sigma_theta, sigma_phi), specifies smoothing independently. Defaults to 0.0 (no smoothing).
incidence_max (float | None) – Maximum allowed phi angle (in radians) for directionality calculation. If None, no clipping is applied.

Note

This is a base class. Users should use the concrete implementations ForwardHomogeneous or ForwardLayeredLookup which inherit from this class.

__call__(m, design)¶

Calculate travel times and signal-to-noise ratios.

Parameters:

m (Tensor) – torch.Tensor Source coordinate tensor with shape (..., 3) containing XYZ coordinates in meters. Last dimension must be [x, y, z].
design (GeoDataFrame) – geopandas.GeoDataFrame GeoDataFrame with sensor locations containing a geometry column with Point geometries and an elevation column with elevations in meters. Optionally includes orientation columns u_x, u_y, u_z.

Returns:

Tuple[torch.Tensor, torch.Tensor]: Tuple of (travel_times, snr) with shapes matching the input source dimensions plus the number of sensors.

Return type:

Tuple[Tensor, Tensor]

snr(distance, wave_direction, sensor_direction)¶

Calculate SNR based on distance and directional sensitivity.

Parameters:

distance (Tensor) – torch.Tensor Distance between source and sensor.
wave_direction (Tensor) – torch.Tensor Direction vector of the wave.
sensor_direction (Tensor) – torch.Tensor Direction vector of the sensor.

Returns:

torch.Tensor: SNR values clamped to maximum of 100.0.

Return type:

Tensor

directional_sensitivity(m, design)¶

Calculate directional sensitivity between sources and sensors.

This is a convenience function that extracts the necessary data from the design GeoDataFrame, calculates wave directions, and returns the directional sensitivity factor for each source-sensor pair.

Parameters:

m (Tensor) – torch.Tensor Source coordinate tensor with shape (..., 3) containing XYZ coordinates in meters. Last dimension must be [x, y, z].
design (GeoDataFrame) – geopandas.GeoDataFrame GeoDataFrame with sensor locations containing a geometry column with Point geometries and an elevation column with elevations in meters. Should include orientation columns u_x, u_y, u_z for directional sensors.

Returns:

torch.Tensor: Directional sensitivity values with shape matching the input source dimensions plus the number of sensors. Values range from 0 (no sensitivity) to 1 (maximum sensitivity).

Return type:

Tensor

calculate_distances(m_in, design_coords)¶

Calculate distances between sources and receivers.

Parameters:

m_in (Tensor) – torch.Tensor Source coordinates tensor with shape (B, 3).
design_coords (Tensor) – torch.Tensor Receiver coordinates tensor with shape (N, 3).

Returns:

Tuple[torch.Tensor, torch.Tensor]: Tuple of (total_distance, horizontal_distance) with shapes (B, N).

Return type:

Tuple[Tensor, Tensor]

abstract calculate_travel_times_and_directions(m_in, design_coords, distance=None)¶

Calculate travel times and wave directions. Must be implemented by subclasses.

Parameters:

m_in (Tensor) – torch.Tensor Source coordinates tensor.
design_coords (Tensor) – torch.Tensor Receiver coordinates tensor.
distance (Tensor | None) – Optional[torch.Tensor] Precalculated distances (optional).

Returns:

Tuple[torch.Tensor, torch.Tensor]: Tuple of (travel_times, wave_directions).

Return type:

Tuple[Tensor, Tensor]

class dased.helpers.srcloc.ForwardHomogeneous(velocity, wave_type, distance_relation, sensor_type='strain', direction_smoothing=0.0, incidence_max=None)¶

Calculates travel times and SNRs assuming a constant velocity medium. This is the simplest forward model, suitable for preliminary analysis or when detailed earth structure is unknown.

Parameters:

velocity (float) – Wave velocity in the medium (m/s). Should be positive.
wave_type (str) – Type of wave (‘P’, ‘S’, ‘SV’, ‘SH’, etc.). Affects directional sensitivity calculations for strain sensors.
distance_relation (Callable) – Function to calculate distance-based attenuation. Should accept distance in meters and return amplitude/SNR values.
sensor_type (str) – Type of sensor (‘strain’, ‘displacement’, ‘strain1C’, etc.). Defaults to ‘strain’ which maps to ‘strain1C’.
direction_smoothing (float | Tuple[float, float]) – Standard deviation(s) for directional sensitivity smoothing in radians. Can be float or tuple (sigma_theta, sigma_phi). Defaults to 0.0 (no smoothing).
incidence_max (float | None) – Maximum allowed phi angle (in radians) for directionality calculation. If None, no clipping is applied.

Examples

Basic P-wave forward model:

>>> from dased.helpers.srcloc import ForwardHomogeneous, MagnitudeRelation
>>> # Setup magnitude relation for SNR calculation
>>> mag_rel = MagnitudeRelation(magnitude_factor=1.0, reference_distance=1000)  # :class:`~dased.helpers.srcloc.magnitude_relation.MagnitudeRelation`
>>> # Create forward model with 5 km/s P-wave velocity
>>> forward = ForwardHomogeneous(
...     velocity=5000.0,
...     wave_type='P',
...     distance_relation=mag_rel
... )
>>> # Calculate travel times and SNR
>>> source = torch.tensor([[1000.0, 2000.0, 500.0]])  # x, y, z in meters, :class:`torch.Tensor`
>>> times, snr = forward(source, design_gdf)  # returns (:class:`torch.Tensor`, :class:`torch.Tensor`)

calculate_travel_times_and_directions(m_in, design_coords, distance=None)¶

Calculate travel times and wave directions in a homogeneous medium.

Parameters:

m_in (Tensor) – Source coordinates tensor (B, 3) or (1, 3).
design_coords (Tensor) – Receiver coordinates tensor (N, 3).
distance (Tensor | None) – Precalculated distances (optional) (B, N) or (1, N).

Returns:

Tuple of (travel_times, wave_directions) with shapes (B, N) and (B, N, 3).

Return type:

Tuple[Tensor, Tensor]

class dased.helpers.srcloc.ForwardLayeredLookup(lookup_path, wave_type, distance_relation, sensor_type='strain', direction_smoothing=0.0, incidence_max=None)¶

Uses precomputed travel times and incidence angles from lookup tables generated by ray tracing through layered earth models. This provides more accurate modeling for realistic earth structure.

Parameters:

lookup_path (str) – Path to the lookup table file in NetCDF format. The file must contain coordinates ‘receiver_depth’, ‘distance’, ‘source_depth’ and data variables ‘arrival_time’, ‘incidence_angle’.
wave_type (str) – Type of wave (‘P’, ‘S’, ‘SV’, ‘SH’, etc.).
distance_relation (Callable) – Function to calculate distance-based attenuation.
sensor_type (str) – Type of sensor. Defaults to ‘strain’.
direction_smoothing (float | Tuple[float, float]) – Standard deviation(s) for directional sensitivity smoothing. Defaults to 0.0.
incidence_max (float | None) – Maximum allowed phi angle (in radians) for directionality calculation. If None, no clipping is applied.

Raises:

FileNotFoundError – If lookup table file is not found.
ValueError – If lookup table is missing required coordinates or variables.
IOError – If lookup table cannot be loaded or initialized.

Examples

Using a precomputed lookup table:

>>> forward = ForwardLayeredLookup(
...     lookup_path='travel_times_P.nc',
...     wave_type='P',
...     distance_relation=mag_rel  # :class:`~dased.helpers.srcloc.magnitude_relation.MagnitudeRelation`
... )
>>> times, snr = forward(source_locations, design_gdf)  # returns (:class:`torch.Tensor`, :class:`torch.Tensor`)

Note

The lookup table assumes depths are positive downward from the surface. Coordinates in the lookup table should use the same units as the source and receiver coordinates (typically meters).

calculate_travel_times_and_directions(m_in, design_coords, distance=None)¶

Calculate travel times and wave directions using lookup tables.

Parameters:

m_in (Tensor) – Source coordinates tensor (B, 3) or (1, 3).
design_coords (Tensor) – Receiver coordinates tensor (N, 3).
distance (Tensor | None) – Precalculated total distances (optional) (B, N) or (1, N).

Returns:

Tuple of (travel_times, wave_directions) with shapes (B, N) and (B, N, 3).

Return type:

Tuple[Tensor, Tensor]

class dased.helpers.srcloc.MagnitudeRelation(log_coeff=0, magnitude_factor=1.0, reference_distance=None, reference_magnitude=None, reference_relation=None)¶

Magnitude-distance relationship for seismic signal amplitude modeling.

Implements the relationship:

log10(Amplitude) = magnitude_factor * M_l + log_coeff * log10(R)

Where R is the distance in kilometers and M_l is the local magnitude. This relationship is commonly used in seismology to relate observed signal amplitudes to earthquake magnitude and distance.

Parameters:

log_coeff (float) – Coefficient for the log10(distance) term. Controls how amplitude decreases with distance. Typical values range from -1.0 to -3.0. Defaults to 0 (no distance dependence).
magnitude_factor (float) – Coefficient for the magnitude term. Controls how amplitude scales with magnitude. Defaults to 1.0.
reference_distance (float | None) – Reference distance in meters for relative amplitude calculations. If provided, :pyattr:`reference_magnitude` is calculated assuming amplitude=1 at this distance. Cannot be used with reference_magnitude or reference_relation.
reference_magnitude (float | None) – Explicit reference magnitude value. Overrides calculation from reference_distance. Cannot be used with reference_distance or reference_relation.
reference_relation (MagnitudeRelation | None) – Another MagnitudeRelation instance to inherit the reference magnitude from. Overrides both reference_magnitude and reference_distance.

Raises:

ValueError – If magnitude_factor is zero, or if none of the reference parameters are provided, or if multiple reference parameters are given.

Examples

Standard magnitude-distance relationship:

>>> # Richter magnitude relation with -2.0 distance decay
>>> relation = MagnitudeRelation(
...     log_coeff=-2.0,
...     magnitude_factor=1.0,
...     reference_distance=1000.0  # 1 km reference
... )
>>> amplitude = relation(distance=5000.0)  # 5 km distance, returns :class:`torch.Tensor`

Using magnitude and distance to calculate amplitude:

>>> magnitude = 4.5
>>> distance = 2000.0  # 2 km
>>> amplitude = relation.get_amplitude(magnitude, distance)  # :class:`torch.Tensor`

Converting observed amplitude to magnitude:

>>> observed_amplitude = 0.01
>>> distance = 3000.0  # 3 km
>>> estimated_magnitude = relation.get_magnitude(observed_amplitude, distance)  # :class:`torch.Tensor`

get_magnitude(amplitude, distance)¶

Calculate the magnitude (M_l) based on observed amplitude and distance.

Uses the inverse relationship:

M_l = (log10(Amplitude) - log_coeff * log10(R)) / magnitude_factor

Parameters:

amplitude (Tensor | float) – Observed amplitude(s) of the signal. Must be positive.
distance (Tensor | float) – Distance(s) from the source in meters.

Returns:

Calculated magnitude(s) as a torch.Tensor.

Return type:

Tensor

Note

Distance is converted to kilometers internally for the calculation.

get_amplitude(magnitude, distance)¶

Calculate the amplitude based on magnitude (M_l) and distance.

Uses the relationship:

Amplitude = 10^(magnitude_factor * M_l + log_coeff * log10(R))

Parameters:

magnitude (Tensor | float) – Local magnitude (M_l) value(s).
distance (Tensor | float) – Distance(s) from the source in meters.

Returns:

Calculated amplitude(s) as a torch.Tensor.

Return type:

Tensor

Note

Distance is converted to kilometers internally for the calculation.

__call__(distance)¶

Calculate amplitude for given distance using the stored reference magnitude.

This represents the amplitude relative to a baseline defined by the reference magnitude, useful for signal-to-noise ratio calculations.

Parameters:: distance (Tensor | float) – Distance(s) from the source in meters.
Returns:: Calculated amplitude(s) for the reference magnitude as a torch.Tensor.
Return type:: Tensor

class dased.helpers.srcloc.DataLikelihood(forward_function, std_corr, std_uncorr, cor_length=0.0, std_cutoff=10.0, f_max=10.0, K_sh=10.0, data_type_design_map=None, remove_mean=True)¶

Computes likelihood of observed data given source model parameters and experimental design.

This class handles computation of data likelihood for Bayesian source localization. It supports correlated and uncorrelated noise models and multiple data types. Docstrings reference types using Sphinx roles such as torch.Tensor, pandas.DataFrame and geopandas.GeoDataFrame.

Parameters:

forward_function (Dict[str, Callable] | Callable) –
Forward model(s) for each data type. Can be:
- Single callable: Used for all data (single data type)
- Dict of callables: Keys are data type names, values are forward functions
Each forward function should accept (model_samples, design) and return either predictions or (predictions, snr_values) tuple.
std_corr (Dict[str, Callable] | Callable | float) –
Correlated standard deviation specification. Can be:
- Float: Constant correlated uncertainty for all data
- Callable: Function of predicted data returning std values
- Dict: Mapping data types to float/callable specifications
std_uncorr (Dict[str, Callable] | Callable | float) – Uncorrelated standard deviation specification. Same format options as std_corr.
cor_length (float) – Correlation length for spatial correlation in meters. If 0.0, assumes uncorrelated measurements. Defaults to 0.0.
std_cutoff (float) – Maximum allowed standard deviation for SNR conversion. Used to prevent numerical issues. Defaults to 10.0.
f_max (float) – Maximum frequency for SNR to standard deviation conversion using Shannon-Hartley theorem. Defaults to 10.0 Hz.
K_sh (float) – Constant for Shannon-Hartley SNR conversion. Defaults to 10.0.
data_type_design_map (Dict[str, List[int]] | None) – Optional mapping from data type names to lists of design row indices. If None, all design rows are used for all data types. Useful when different data types use different subsets of sensors.

Examples

Basic setup with single forward model:

>>> from dased.helpers.srcloc import DataLikelihood, ForwardHomogeneous
>>> from dased.helpers.srcloc import MagnitudeRelation
>>>
>>> # Setup forward model
>>> mag_rel = MagnitudeRelation(magnitude_factor=1.0, reference_distance=1000)
>>> forward = ForwardHomogeneous(velocity=5000, wave_type='P', distance_relation=mag_rel)  # :class:`~dased.helpers.srcloc.forward.ForwardHomogeneous`
>>>
>>> # Create likelihood with uncorrelated noise
>>> likelihood = DataLikelihood(  # :class:`~dased.helpers.srcloc.data_likelihood.DataLikelihood`
...     forward_function=forward,
...     std_corr=0.1,     # 0.1 second correlated uncertainty
...     std_uncorr=0.05,  # 0.05 second uncorrelated uncertainty
...     cor_length=0.0    # No spatial correlation
... )

Multiple data types with different forward models:

>>> forward_models = {
...     'P_wave': forward_p,
...     'S_wave': forward_s
... }
>>> std_specs = {
...     'P_wave': 0.1,
...     'S_wave': 0.2
... }
>>> likelihood = DataLikelihood(
...     forward_function=forward_models,
...     std_corr=std_specs,
...     std_uncorr=0.05
... )

__call__(model_samples, design, posterior=None)¶

Compute the data likelihood distribution for given model samples and design.

Parameters:

model_samples (Tensor) – Model parameter samples with shape (…, n_params). Last dimension typically contains [x, y, z] coordinates.
design – Design DataFrame containing sensor information. Must include ‘geometry’ column with Point geometries and ‘elevation’ column. For correlated case, geometries are used for distance calculations.

Returns:

Likelihood distribution. Returns IndependentNormal for uncorrelated case or MultivariateNormal for correlated case.

Return type:

WeightedIndependentNormal | WeightedMultivariateNormal

snr2std(snr)¶

Convert SNR to standard deviation using the Shannon-Hartley theorem.

Parameters:: snr (Tensor | None) – SNR tensor or None.
Returns:: Standard deviation tensor.
Return type:: Tensor

class dased.helpers.srcloc.SurfaceField_Distribution(distribution, topo_data, depth=200.0)¶

3D spatial distribution that follows surface topography with depth variation.

This distribution generates 3D points where the horizontal coordinates follow a specified 2D distribution and the vertical coordinate is distributed between the surface topography and a maximum depth below the surface.

Parameters:

distribution (Distribution) – 2D probability distribution for horizontal coordinates. Should be a PyTorch distribution that samples (x, y) pairs.
topo_data (DataArray) – xarray.DataArray containing topography data with 'x' and 'y' coordinates and elevation values. Used for elevation at sampled horizontal coordinates.
depth (float) – Maximum depth below surface in meters. Depth is sampled uniformly between the surface and this maximum depth. Must be positive. Defaults to 200 meters.

Examples

Create distribution over mountainous terrain:

>>> import torch.distributions as dist
>>> # Define horizontal distribution
>>> horizontal_dist = dist.Uniform(
...     torch.tensor([0.0, 0.0]),
...     torch.tensor([10000.0, 10000.0])
... )
>>> # Assume topo_data is loaded as :class:`~xarray.DataArray`
>>> spatial_dist = SurfaceField_Distribution(
...     distribution=horizontal_dist,
...     topo_data=topo_data,
...     depth=500.0  # Up to 500m below surface
... )
>>> samples = spatial_dist.sample((1000,))  # Sample 1000 3D points

sample(sample_shape=1)¶

Sample 3D points following surface topography.

Parameters:: sample_shape (int | tuple | Size) – Shape of samples to generate. Can be int, tuple, or Size.
Returns:: Tensor of 3D points with shape sample_shape + (3,).
Return type:: Tensor

log_prob(value, fast_eval=True)¶

Calculate log probability density for given points.

Parameters:

value (Tensor) – Points to evaluate with shape (…, 3).
fast_eval (bool) – If True, uses approximate elevation for speed. If False, interpolates actual elevation for each point (slower but more accurate).

Returns:

Log probability densities with shape (…,).

Return type:

Tensor

class dased.helpers.srcloc.PolygonUniform(polygon, vertical_dist)¶

Uniform distribution within a Shapely polygon with custom vertical distribution.

This distribution generates 3D points where the horizontal coordinates are uniformly distributed within a specified polygon and the vertical coordinate follows a given distribution.

Parameters:

polygon (Polygon) – Shapely Polygon object defining the horizontal bounds.
vertical_dist (Distribution) – PyTorch distribution for the vertical (z) component. Should be a 1D distribution.

Examples

Create uniform distribution within rectangular area:

>>> from shapely.geometry import Polygon
>>> import torch.distributions as dist
>>> # Define rectangular polygon
>>> coords = [(0, 0), (1000, 0), (1000, 500), (0, 500), (0, 0)]
>>> polygon = Polygon(coords)
>>> # Define vertical distribution (uniform from 0 to 1000m depth)
>>> vertical_dist = dist.Uniform(0, 1000)
>>> # Create 3D distribution
>>> spatial_dist = PolygonUniform(polygon, vertical_dist)
>>> samples = spatial_dist.sample((500,))  # :class:`torch.Tensor` of shape (500, 3)

Note

This distribution uses rejection sampling for horizontal coordinates, which may be slow for polygons with low area-to-bounding-box ratios. The implementation includes adaptive oversampling to improve efficiency.

sample(sample_shape=())¶

Sample points uniformly from within the polygon.

Parameters:: sample_shape (Size) – Shape of samples to generate.
Returns:: Tensor of 3D points with shape sample_shape + batch_shape + (3,).
Return type:: Tensor

log_prob(value)¶

Calculate log probability density for points.

Parameters:: value (Tensor) – Points to evaluate with shape (…, 3).
Returns:: Log probability densities.
Return type:: Tensor

class dased.helpers.srcloc.IndependentNormal(loc, scale)¶

Independent Normal distribution with element-wise independence.

This distribution represents a multivariate normal distribution where each component is independent (diagonal covariance matrix). It’s more efficient than a full multivariate normal when correlations are not needed.

Parameters:

loc (Tensor) – Mean of the normal distribution with shape (…, n).
scale (Tensor) – Standard deviation of the normal distribution with shape (…, n).

Examples

Create independent normal distribution:

>>> loc = torch.tensor([0.0, 1.0, -0.5])  # :class:`torch.Tensor`
>>> scale = torch.tensor([1.0, 2.0, 0.5])  # :class:`torch.Tensor`
>>> dist = IndependentNormal(loc, scale)
>>> samples = dist.sample((100,))  # Sample 100 points

log_prob(value)¶

Compute log probability density.

sample(sample_shape=())¶

Sample from the distribution.

rsample(sample_shape=())¶

Reparameterized sampling for gradient flow.

expand(batch_shape, _instance=None)¶

Returns a new distribution with expanded batch dimensions.

property mean: Tensor¶: Compute the mean of the distribution.

property variance: Tensor¶: Compute the variance of the distribution.

property stddev: Tensor¶: Compute the standard deviation of the distribution.

dased.helpers.srcloc.calculate_posterior(true_source, design, data_likelihood, x_range=(-5000.0, 5000.0), y_range=(-5000.0, 5000.0), z_range=None, grid_size=100, prior_dist=None, clean_data=True, downsample=None, seed=0, posterior_std_obs=0.02)¶

Calculate Bayesian posterior distribution over a 2D or 3D grid.

This function computes the posterior probability distribution for source locations by evaluating the likelihood function over a regular grid and combining it with a prior distribution (uniform if not specified).

Types in this module use Sphinx roles, e.g. torch.Tensor for the grid and DataLikelihood for likelihood objects. The design argument typically expects a geopandas.GeoDataFrame.

Parameters:

true_source (Tensor) – torch.Tensor True source location tensor with shape (3,) containing [x, y, z] coordinates in meters. Used to generate synthetic data and as reference for evaluation.
design – geopandas.GeoDataFrame Design GeoDataFrame containing sensor information. Must include a geometry column with Point geometries and an elevation column. Used by the data likelihood function.
data_likelihood – DataLikelihood DataLikelihood instance for computing likelihood values. Should be properly configured with forward models and uncertainty specifications.
x_range (Tuple[float, float]) – Tuple of (min, max) values for x-axis in meters. Defines the spatial extent of the grid in the x-direction. Defaults to (-5000.0, 5000.0).
y_range (Tuple[float, float]) – Tuple of (min, max) values for y-axis in meters. Defines the spatial extent of the grid in the y-direction. Defaults to (-5000.0, 5000.0).
z_range (Tuple[float, float] | None) – Optional tuple of (min, max) values for z-axis in meters. If None, computes 2D posterior at the true source depth. If provided, computes 3D posterior over the specified depth range.
grid_size (int | Tuple[int, ...]) – Number of grid points in each dimension. Can be:
Int (-) – Same size for all dimensions (nx=ny=nz=grid_size)
Tuple (-) –
Specify size for each dimension (nx, ny) or (nx, ny, nz)

Larger values give higher resolution but increase computation time. Defaults to 100.
prior_dist (Callable | None) – Optional prior distribution. Should be a callable that accepts grid points and returns log prior probabilities. If None, uses uniform prior. Defaults to None.
clean_data (bool) – Whether to use clean synthetic data (True) or add noise (False). Clean data is useful for testing and validation. Defaults to True.
downsample (int | None) – Optional downsampling factor for sensors. If provided, uses every downsample-th sensor to reduce computation. Useful for large sensor arrays. Defaults to None (use all sensors).
seed (int) – Random seed for reproducibility of noise generation and sampling. Defaults to 0.
posterior_std_obs (float) – Fill-in value for observation uncertainty in posterior computation. Defaults to 0.02.
Returns –
For 2D case: Tuple of (torch.Tensor, torch.Tensor, torch.Tensor, dict) For 3D case: Tuple of (torch.Tensor, torch.Tensor, torch.Tensor, torch.Tensor, dict)

Where: - x_grid, y_grid, z_grid: Grid coordinate tensors as torch.Tensor - log_posterior: Log posterior probabilities as torch.Tensor with same shape as coordinate grids - post_info: Dictionary with additional information including design (geopandas.GeoDataFrame),

observed data (torch.Tensor), data likelihood distribution (DataLikelihood), etc.

Examples

Basic 2D posterior calculation:

>>> import torch
>>> # Setup true source location
>>> true_source = torch.tensor([1000.0, 2000.0, 500.0])
>>> # Calculate 2D posterior
>>> x_grid, y_grid, log_post, info = calculate_posterior(
...     true_source=true_source,
...     design=design_gdf,
...     data_likelihood=likelihood,
...     x_range=(-2000, 4000),
...     y_range=(-1000, 5000),
...     grid_size=50
... )
>>> # Find maximum likelihood location
>>> max_idx = torch.argmax(log_post)
>>> max_x = x_grid[max_idx // log_post.shape[1]]
>>> max_y = y_grid[max_idx % log_post.shape[1]]

3D posterior with custom prior:

>>> import torch
>>> # Define prior favoring shallow sources
>>> class ShallowPrior:
...     def log_prob(self, points):
...         z_values = points[..., 2]
...         return -0.001 * z_values
>>>
>>> shallow_prior = ShallowPrior()
>>> x_grid, y_grid, z_grid, log_post, info = calculate_posterior(
...     true_source=true_source,
...     design=design_gdf,
...     data_likelihood=likelihood,
...     z_range=(0, 2000),  # 0 to 2 km depth
...     grid_size=(40, 40, 20),  # Higher resolution in x,y
...     prior_dist=shallow_prior
... )