Fields

Unified field API entrypoint.

class gwexpy.fields.ScalarField(data, unit=None, axis0=None, axis1=None, axis2=None, axis3=None, axis_names=None, axis0_domain: str = 'time', space_domain: str | dict[str, str] = 'real', **kwargs)[source]

Bases: FieldBase

4D Field with domain states and FFT operations.

This class extends Array4D to represent physical fields that can exist in different domains (time/frequency for axis 0, real/k-space for spatial axes 1-3).

Key feature: All indexing operations return a ScalarField, maintaining 4D structure. Integer indices result in axes with length 1 rather than being dropped.

Parameters:

  • data (array-like) – 4-dimensional input data.

  • unit (~astropy.units.Unit, optional) – Physical unit of the data.

  • axis0 (~astropy.units.Quantity or array-like, optional) – Index values for axis 0 (time or frequency).

  • axis1 (~astropy.units.Quantity or array-like, optional) – Index values for axis 1 (x or kx).

  • axis2 (~astropy.units.Quantity or array-like, optional) – Index values for axis 2 (y or ky).

  • axis3 (~astropy.units.Quantity or array-like, optional) – Index values for axis 3 (z or kz).

  • axis_names (iterable of str, optional) – Names for each axis (length 4). Defaults based on domain.

  • axis0_domain ({'time', 'frequency'}, optional) – Domain of axis 0. Default is ‘time’.

  • space_domain ({'real', 'k'} or dict, optional) – Domain of spatial axes. If str, applies to all spatial axes. If dict, maps axis names to domains. Default is ‘real’.

  • **kwargs – Additional keyword arguments passed to Array4D.

axis0_domain

Current domain of axis 0 (‘time’ or ‘frequency’).

Type:

str

space_domains

Mapping of spatial axis names to their domains (‘real’ or ‘k’).

Type:

dict

Examples

>>> import numpy as np
>>> from astropy import units as u
>>> from gwexpy.fields import ScalarField
>>> data = np.random.randn(100, 32, 32, 32)
>>> times = np.arange(100) * 0.01 * u.s
>>> x = np.arange(32) * 1.0 * u.m
>>> field = ScalarField(data, axis0=times, axis1=x, axis2=x, axis3=x,
...                 axis_names=['t', 'x', 'y', 'z'])
coherence_map(ref_point, plane='xy', at=None, **kwargs)[source]

Compute coherence map from a reference point.

This is a convenience wrapper around coherence_map().

Parameters:

  • ref_point (tuple of Quantity) – Reference point coordinates (x, y, z).

  • plane (str) – 2D plane to map: ‘xy’, ‘xz’, or ‘yz’. Default ‘xy’.

  • at (dict, optional) – Fixed value for the axis not in the plane.

  • **kwargs – Additional keyword arguments passed to coherence_map.

Returns:

Coherence map.

Return type:

ScalarField or FieldDict

See also

gwexpy.fields.signal.coherence_map

Full documentation.

compute_psd(point_or_region, **kwargs)[source]

Compute power spectral density using Welch’s method.

This is a convenience wrapper around compute_psd().

Parameters:

  • point_or_region (tuple, list of tuples, or dict) – Spatial location(s) to extract: - Single point: (x, y, z) tuple of Quantities - Multiple points: list of (x, y, z) tuples - Region dict: {'x': slice or value, 'y': ..., 'z': ...}

  • **kwargs – Additional keyword arguments passed to compute_psd: nperseg, noverlap, window, detrend, scaling, average.

Returns:

PSD estimate(s).

Return type:

FrequencySeries or FrequencySeriesList

See also

gwexpy.fields.signal.compute_psd

Full documentation.

compute_xcorr(point_a, point_b, **kwargs)[source]

Compute cross-correlation between two spatial points.

This is a convenience wrapper around compute_xcorr().

Parameters:

  • point_a (tuple of Quantity) – Spatial coordinates (x, y, z) for the two points.

  • point_b (tuple of Quantity) – Spatial coordinates (x, y, z) for the two points.

  • **kwargs – Additional keyword arguments passed to compute_xcorr.

Returns:

Cross-correlation function with lag axis.

Return type:

TimeSeries

See also

gwexpy.fields.signal.compute_xcorr

Full documentation.

diff(other, mode='diff')[source]

Compute difference or ratio between two ScalarField objects.

Parameters:

  • other (ScalarField) – The field to compare against.

  • mode (str, optional) – Comparison mode: - ‘diff’: Difference (self - other) - ‘ratio’: Ratio (self / other) - ‘percent’: Percentage difference ((self - other) / other * 100) Default is ‘diff’.

Returns:

Result field. For ‘diff’, unit is same as input. For ‘ratio’ and ‘percent’, unit is dimensionless.

Return type:

ScalarField

Raises:

  • ValueError – If mode is not recognized.

  • ValueError – If shapes are incompatible.

Examples

>>> diff_field = field1.diff(field2)
>>> ratio_field = field1.diff(field2, mode='ratio')
extract_points(points: Sequence[Sequence[Quantity]] | Sequence[tuple[Quantity, ...]], interp: str = 'linear') → Any[source]

Extract time series at specified spatial points.

Parameters:

  • points (list of tuple) – List of (x, y, z) coordinates.

  • interp (str, optional) – Interpolation method (‘linear’ or ‘nearest’). Default is ‘linear’.

Returns:

List of time series, one per point.

Return type:

TimeSeriesList

extract_profile(axis: str, at: dict[str, u.Quantity], reduce: str | None = None, interp: str = 'linear') → tuple[IndexLike, Any][source]

Extract a 1D profile along a specified axis.

Parameters:

  • axis (str) – Axis name to extract along (‘x’, ‘y’, ‘z’, or their k-variants).

  • at (dict) – Dictionary specifying fixed values for the other axes. Must include all axes except the extraction axis. Example: {'t': 0.5 * u.s, 'y': 2.0 * u.m, 'z': 0.0 * u.m}

  • reduce (None) – Reserved for future averaging support. Currently ignored.

  • interp (str, optional) – Interpolation method (‘linear’ or ‘nearest’). Default is ‘linear’.

Returns:

(axis_index, values): Both are ~astropy.units.Quantity arrays. axis_index is the coordinate along the extraction axis, values is the data values along that axis.

Return type:

tuple

fft_space(axes: Iterable[str] | None = None, n: Sequence[int] | None = None, overwrite: bool = False) → ScalarField[source]

Compute FFT along spatial axes.

This method uses two-sided FFT (scipy.fft.fftn) and produces angular wavenumber (k = 2π·fftfreq).

Parameters:

  • axes (iterable of str, optional) – Axis names to transform (e.g., [‘x’, ‘y’]). If None, transforms all spatial axes in ‘real’ domain.

  • n (tuple of int, optional) – FFT lengths for each axis.

  • overwrite (bool, optional) – If True, perform FFT in-place on a temporary copy of the data to reduce peak memory usage. Default is False.

Returns:

Transformed field with specified axes in ‘k’ domain.

Return type:

ScalarField

Raises:

  • ValueError – If any specified axis is not in ‘real’ domain.

  • ValueError – If any specified axis is not uniformly spaced.

Notes

Angular Wavenumber Convention

The wavenumber axis is computed as k = 2π * fftfreq(n, d=dx), satisfying k = 2π / λ. This is the standard angular wavenumber definition in physics, with units of [rad/length].

Note: This is NOT the cycle wavenumber (1/λ) commonly used in some fields. To convert: k_cycle = k_angular / (2π).

Sign Convention for Descending Axes (dx < 0)

If the spatial axis is descending (dx < 0), the k-axis is sign-flipped to preserve physical consistency with the phase factor convention e^{+ikx}. This ensures that positive k corresponds to waves propagating in the positive x direction, regardless of the data storage order.

This convention differs from the standard FFT behavior (which ignores axis direction) but maintains physical consistency for interferometer simulations and wave propagation analysis.

This formula has been validated through unit tests and independent technical review (2026-02-01). The 2π factor is correctly applied, and units are properly set as 1/dx_unit.

References

fft_time(nfft: int | None = None) → ScalarField[source]

Compute FFT along time axis (axis 0).

This method applies the same normalization as GWpy’s TimeSeries.fft(): rfft / nfft, with DC-excluded bins multiplied by 2 (except Nyquist bin for even nfft).

Parameters:

nfft (int, optional) – Length of the FFT. If None, uses the length of axis 0.

Returns:

Transformed field with axis0_domain='frequency'.

Return type:

ScalarField

Raises:

  • ValueError – If axis0_domain is not ‘time’.

  • ValueError – If time axis length < 2 or is irregularly spaced.

  • TypeError – If input data is complex-valued.

See also

gwpy.timeseries.TimeSeries.fft

The reference implementation.

filter(*args, **kwargs) → ScalarField[source]

Apply a filter along the time axis (axis 0).

Parameters:

  • *args – Filter specification. Supports same arguments as gwpy.timeseries.TimeSeries.filter().

  • **kwargs – Filter specification. Supports same arguments as gwpy.timeseries.TimeSeries.filter().

Returns:

Filtered field.

Return type:

ScalarField

freq_space_map(axis, at=None, **kwargs)[source]

Compute frequency-space map along a spatial axis.

This is a convenience wrapper around freq_space_map().

Parameters:

  • axis (str) – Spatial axis to scan along (‘x’, ‘y’, or ‘z’).

  • at (dict, optional) – Fixed values for the other two spatial axes.

  • **kwargs – Additional keyword arguments passed to freq_space_map.

Returns:

2D frequency-space map.

Return type:

ScalarField

See also

gwexpy.fields.signal.freq_space_map

Full documentation.

classmethod from_pyroomacoustics_field(room: Any, *, grid_shape: tuple[int, ...], source: int = 0, mode: str = 'rir', unit: Any | None = None) → ScalarField[source]

Create from pyroomacoustics room with grid-placed microphones.

Parameters:

  • room (pyroomacoustics.Room) – Room with microphones on a regular spatial grid.

  • grid_shape (tuple of int) – Spatial grid shape (nx, ny, nz) or (nx, ny).

  • source (int, default 0) – Source index (for mode='rir').

  • mode ({'rir', 'signals'}) – 'rir' for impulse responses, 'signals' for mic signals.

  • unit (str or astropy.units.Unit, optional) – Unit to assign to the data.

Return type:

ScalarField

ifft_space(axes: Iterable[str] | None = None, n: Sequence[int] | None = None, overwrite: bool = False) → ScalarField[source]

Compute inverse FFT along k-space axes.

Parameters:

  • axes (iterable of str, optional) – Axis names to transform (e.g., [‘kx’, ‘ky’]). If None, transforms all spatial axes in ‘k’ domain.

  • n (tuple of int, optional) – Output lengths for each axis.

  • overwrite (bool, optional) – If True, perform IFFT in-place on a temporary copy.

Returns:

Transformed field with specified axes in ‘real’ domain.

Return type:

ScalarField

Raises:

ValueError – If any specified axis is not in ‘k’ domain.

ifft_time(nout: int | None = None) → ScalarField[source]

Compute inverse FFT along frequency axis (axis 0).

This method applies the inverse normalization of fft_time() / GWpy’s FrequencySeries.ifft().

Parameters:

nout (int, optional) – Length of the output time series. If None, computed as (n_freq - 1) * 2.

Returns:

Transformed field with axis0_domain='time'.

Return type:

ScalarField

Raises:

  • ValueError – If axis0_domain is not ‘frequency’.

  • ValueError – If frequency axis length < 2 or is irregularly spaced.

See also

gwpy.frequencyseries.FrequencySeries.ifft

Reference implementation.

plot_map2d(plane='xy', at=None, mode='real', method='pcolormesh', ax=None, add_colorbar=True, vmin=None, vmax=None, title=None, cmap=None, **kwargs)[source]

Plot a 2D map (heatmap) of the field.

Parameters:

  • plane (str, optional) – The plane to plot: ‘xy’, ‘xz’, ‘yz’, ‘tx’, ‘ty’, ‘tz’. Default is ‘xy’.

  • at (dict, optional) – Dictionary specifying fixed values for axes not in the plane. If None, axes with length=1 are used automatically.

  • mode (str, optional) – Component to extract from complex data: ‘real’, ‘imag’, ‘abs’, ‘angle’, ‘power’. Default is ‘real’.

  • method (str, optional) – Plot method: ‘pcolormesh’ or ‘imshow’. Default is ‘pcolormesh’.

  • ax (matplotlib.axes.Axes, optional) – Axes to plot on. If None, creates new figure.

  • add_colorbar (bool, optional) – Whether to add a colorbar. Default is True.

  • vmin (float, optional) – Color scale limits.

  • vmax (float, optional) – Color scale limits.

  • title (str, optional) – Plot title.

  • cmap (str or Colormap, optional) – Colormap to use.

  • **kwargs – Additional arguments passed to the plot method.

Returns:

(fig, ax): The matplotlib figure and axes objects.

Return type:

tuple

Examples

>>> fig, ax = field.plot_map2d('xy', at={'t': 0.5 * u.s, 'z': 0.0 * u.m})
plot_profile(axis, at, mode='real', ax=None, label=None, **kwargs)[source]

Plot a 1D profile along a specified axis.

Parameters:

  • axis (str) – Axis name to plot along.

  • at (dict) – Dictionary specifying fixed values for other axes.

  • mode (str, optional) – Component to extract: ‘real’, ‘imag’, ‘abs’, ‘angle’, ‘power’. Default is ‘real’.

  • ax (matplotlib.axes.Axes, optional) – Axes to plot on. If None, creates new figure.

  • label (str, optional) – Line label for legend.

  • **kwargs – Additional arguments passed to plot.

Returns:

(fig, ax): The matplotlib figure and axes objects.

Return type:

tuple

Examples

>>> fig, ax = field.plot_profile(
...     'x', at={'t': 0.5 * u.s, 'y': 0.0 * u.m, 'z': 0.0 * u.m}
... )
plot_time_space_map(axis='x', at=None, mode='real', method='pcolormesh', ax=None, add_colorbar=True, vmin=None, vmax=None, title=None, cmap=None, **kwargs)[source]

Plot a 2D time-space map (t vs one spatial axis).

Parameters:

  • axis (str, optional) – Spatial axis name. Default is ‘x’.

  • at (dict, optional) – Fixed values for other spatial axes.

  • mode (str, optional) – Component to extract. Default is ‘real’.

  • method (str, optional) – Plot method. Default is ‘pcolormesh’.

  • ax (matplotlib.axes.Axes, optional) – Axes to plot on.

  • add_colorbar (bool, optional) – Whether to add colorbar. Default is True.

  • vmin (float, optional) – Color scale limits.

  • vmax (float, optional) – Color scale limits.

  • title (str, optional) – Plot title.

  • cmap (str or Colormap, optional) – Colormap to use.

  • **kwargs – Additional plot arguments.

Returns:

(fig, ax): The matplotlib figure and axes objects.

Return type:

tuple

Examples

>>> fig, ax = field.plot_time_space_map('x', at={'y': 0*u.m, 'z': 0*u.m})
plot_timeseries_points(points, labels=None, interp='nearest', ax=None, legend=True, **kwargs)[source]

Plot time series extracted at specified spatial points.

Parameters:

  • points (list of tuple) – List of (x, y, z) coordinates.

  • labels (list of str, optional) – Labels for each time series. If None, auto-generated.

  • interp (str, optional) – Interpolation method. Default is ‘nearest’.

  • ax (matplotlib.axes.Axes, optional) – Axes to plot on. If None, creates new figure.

  • legend (bool, optional) – Whether to show legend. Default is True.

  • **kwargs – Additional arguments passed to plot.

Returns:

(fig, ax): The matplotlib figure and axes objects.

Return type:

tuple

Examples

>>> points = [(1.0 * u.m, 2.0 * u.m, 3.0 * u.m)]
>>> fig, ax = field.plot_timeseries_points(points)
psd(**kwargs)[source]

Compute power spectral density along time axis.

Convenience method equivalent to spectral_density(axis=0). Uses Welch’s method by default for robust PSD estimation.

Parameters:

**kwargs –

Keyword arguments passed to spectral_density() or compute_psd() (if point_or_region is used). Common options:

  • point_or_regiontuple or list, optional

    If provided, computes PSD at specific spatial point(s) or region average instead of full field. Returns FrequencySeries(List).

  • method : {‘welch’, ‘fft’}, default ‘welch’

  • fftlengthfloat, optional

    Segment length in seconds (time-based specification)

  • nfftint, optional

    Number of samples per segment (sample-based specification)

  • overlapfloat, optional

    Overlap in seconds

  • noverlapint, optional

    Number of overlapping samples

  • window : str, window function

  • scaling : {‘density’, ‘spectrum’}

Returns:

PSD field with axis0_domain=’frequency’.

Return type:

ScalarField

Examples

>>> from astropy import units as u
>>> # Time-based specification
>>> psd_field = field.psd(fftlength=1.0, overlap=0.5)
>>> # Or sample-based specification
>>> psd_field = field.psd(nfft=512, noverlap=256)
>>> psd_field.shape  # (n_freq, nx, ny, nz)

See also

spectral_density

Generalized spectral density for any axis.

resample(rate, **kwargs) → ScalarField[source]

Resample the field along the time axis (axis 0).

Parameters:

  • rate (float, Quantity) – The new sampling rate (e.g., in Hz).

  • **kwargs – Additional arguments passed to gwpy.timeseries.TimeSeries.resample().

Returns:

Resampled field.

Return type:

ScalarField

classmethod simulate(method: str, *args: Any, **kwargs: Any) → ScalarField[source]

Generate a simulated ScalarField.

Parameters:

  • method (str) – Name of the generator from gwexpy.noise.field. (e.g., ‘gaussian’, ‘plane_wave’).

  • *args – Arguments passed to the generator.

  • **kwargs – Arguments passed to the generator.

Returns:

Generated field.

Return type:

ScalarField

Examples

>>> from gwexpy.fields import ScalarField
>>> field = ScalarField.simulate('gaussian', shape=(100, 10, 10, 10))
slice_map2d(plane='xy', at=None)[source]

Extract a 2D slice (map) from the 4D field.

Parameters:

  • plane (str, optional) – The plane to extract: ‘xy’, ‘xz’, ‘yz’, ‘tx’, ‘ty’, ‘tz’. Default is ‘xy’.

  • at (dict, optional) – Dictionary specifying fixed values for axes not in the plane. If None, axes with length=1 are used automatically.

Returns:

A ScalarField with the non-plane axes having length=1.

Return type:

ScalarField

Raises:

  • ValueError – If plane specification is invalid.

  • ValueError – If at is None and there is ambiguity about which axes to fix.

Examples

>>> # Extract xy plane at a specific time and z
>>> field_2d = field.slice_map2d('xy', at={'t': 0.5 * u.s, 'z': 0.0 * u.m})
>>> field_2d.plot_map2d()
spectral_density(axis=0, **kwargs)[source]

Compute spectral density along any axis.

Generalized spectral density function that works on time axis (0) or spatial axes (1-3). Returns a new ScalarField with the transformed axis in spectral domain.

Parameters:

  • axis (int or str) – Axis to transform. Default 0 (time axis).

  • **kwargs – Additional arguments passed to spectral_density(). See that function for full parameter list.

Returns:

Spectral density field with transformed axis.

Return type:

ScalarField

Examples

>>> # Time PSD
>>> psd_field = field.spectral_density(axis=0)
>>> psd_field.axis0_domain  # 'frequency'

>>> # Spatial wavenumber spectrum
>>> kx_spec = field.spectral_density(axis='x')

See also

gwexpy.fields.signal.spectral_density

Full documentation.

psd

Convenience alias for time-axis PSD.

time_delay_map(ref_point, plane='xy', at=None, **kwargs)[source]

Compute time delay map from a reference point.

This is a convenience wrapper around time_delay_map().

Parameters:

  • ref_point (tuple of Quantity) – Reference point coordinates (x, y, z).

  • plane (str) – 2D plane to map: ‘xy’, ‘xz’, or ‘yz’. Default ‘xy’.

  • at (dict, optional) – Fixed value for the axis not in the plane.

  • **kwargs – Additional keyword arguments passed to time_delay_map.

Returns:

Time delay map.

Return type:

ScalarField

See also

gwexpy.fields.signal.time_delay_map

Full documentation.

time_space_map(axis='x', at=None, mode='real', reduce=None)[source]

Extract a 2D time-space map (t vs one spatial axis).

Parameters:

  • axis (str, optional) – Spatial axis name (‘x’, ‘y’, ‘z’ or k-variants). Default is ‘x’.

  • at (dict, optional) – Fixed values for the other two spatial axes.

  • mode (str, optional) – Component to extract: ‘real’, ‘imag’, ‘abs’, ‘angle’, ‘power’. Default is ‘real’.

  • reduce (None) – Reserved for future averaging support. Currently ignored.

Returns:

(t_axis, space_axis, data_2d): Quantity arrays for axes and 2D numpy array for the data.

Return type:

tuple

Raises:

  • ValueError – If axis is not valid.

  • ValueError – If at dictionary is missing required axes.

Examples

>>> from astropy import units as u
>>> t, x, data = field.time_space_map('x', at={'y': 0 * u.m, 'z': 0 * u.m})
time_stat_map(stat='mean', t_range=None, plane='xy', at=None)[source]

Compute a time-aggregated 2D map.

Parameters:

  • stat (str, optional) – Statistic to compute: ‘mean’, ‘std’, ‘rms’, ‘max’, ‘min’. Default is ‘mean’.

  • t_range (tuple of Quantity, optional) – Time range (t_start, t_end) to aggregate over. If None, uses the entire time axis.

  • plane (str, optional) – Spatial plane to visualize: ‘xy’, ‘xz’, ‘yz’. Default is ‘xy’.

  • at (dict, optional) – Fixed values for axes not in the plane.

Returns:

Result field with time axis reduced to length=1.

Return type:

ScalarField

Raises:

ValueError – If stat is not recognized.

Examples

>>> from astropy import units as u
>>> mean_map = field.time_stat_map('mean', t_range=(0 * u.s, 1 * u.s))
>>> mean_map.plot_map2d('xy')
wavelength(axis: str | int) → Quantity[source]

Compute wavelength from wavenumber axis.

Parameters:

axis (str or int) – The k-domain axis name or index.

Returns:

Wavelength values (\(\lambda = 2\pi / |k|\)). k=0 returns inf.

Return type:

~astropy.units.Quantity

Raises:

ValueError – If the axis is not in ‘k’ domain.

zscore(baseline_t=None)[source]

Compute z-score normalized field using a baseline period.

The z-score is computed as (data - mean) / std, where mean and std are computed from the baseline period along axis 0 (time).

Parameters:

baseline_t (tuple of Quantity, optional) – Time range (t_start, t_end) for computing baseline statistics. If None, uses the entire time axis.

Returns:

Z-score normalized field (dimensionless).

Return type:

ScalarField

Raises:

  • ValueError – If axis0_domain is not ‘time’.

  • ValueError – If baseline time range is outside the available data.

Examples

>>> from astropy import units as u
>>> zscore_field = field.zscore(baseline_t=(0 * u.s, 1 * u.s))
gwexpy.fields.make_demo_scalar_field(pattern: ~typing.Literal['gaussian', 'sine', 'standing', 'noise'] = 'gaussian', nt: int = 100, nx: int = 16, ny: int = 16, nz: int = 16, dt: ~astropy.units.quantity.Quantity = <Quantity 0.01 s>, dx: ~astropy.units.quantity.Quantity = <Quantity 0.1 m>, dy: ~astropy.units.quantity.Quantity | None = None, dz: ~astropy.units.quantity.Quantity | None = None, unit: ~astropy.units.core.Unit = Unit("m / s"), seed: int = 42, **kwargs) → ScalarField[source]

Generate a deterministic demo ScalarField with known waveforms.

This function provides reproducible ScalarField data for testing, tutorials, and examples. All generated fields have proper axis metadata and units.

Parameters:

  • pattern ({'gaussian', 'sine', 'standing', 'noise'}) – Type of demo pattern: - ‘gaussian’: Propagating Gaussian pulse (default) - ‘sine’: Plane sinusoidal wave - ‘standing’: Standing wave pattern - ‘noise’: White noise (for PSD validation)

  • nt (int) – Number of time points. Default is 100.

  • nx (int) – Number of spatial points per axis. Default is 16.

  • ny (int) – Number of spatial points per axis. Default is 16.

  • nz (int) – Number of spatial points per axis. Default is 16.

  • dt (Quantity) – Time step. Default is 0.01 s.

  • dx (Quantity) – Spatial step for x-axis. Default is 0.1 m.

  • dy (Quantity, optional) – Spatial steps for y and z. If None, uses dx.

  • dz (Quantity, optional) – Spatial steps for y and z. If None, uses dx.

  • unit (Unit) – Physical unit of the data. Default is m/s.

  • seed (int) – Random seed for reproducibility. Default is 42.

  • **kwargs – Additional keyword arguments passed to the specific pattern generator: - For ‘gaussian’: velocity, sigma_t, sigma_x - For ‘sine’: frequency, k_direction - For ‘standing’: mode_x, mode_y, mode_z, omega

Returns:

Deterministic demo field with proper metadata.

Return type:

ScalarField

Examples

>>> from gwexpy.fields.demo import make_demo_scalar_field
>>> field = make_demo_scalar_field('gaussian', nt=50, nx=32)
>>> field.shape
(50, 32, 16, 16)
>>> field.axis0_domain
'time'

>>> # White noise for PSD testing
>>> noise_field = make_demo_scalar_field('noise', seed=123)
gwexpy.fields.make_propagating_gaussian(**kw) → ScalarField[source]

Generate a propagating Gaussian demo field. See make_demo_scalar_field.

gwexpy.fields.make_sinusoidal_wave(**kw) → ScalarField[source]

Generate a sinusoidal wave demo field. See make_demo_scalar_field.

gwexpy.fields.make_standing_wave(**kw) → ScalarField[source]

Generate a standing wave demo field. See make_demo_scalar_field.

gwexpy.fields.compute_psd(field: ScalarField, point_or_region: tuple[Quantity, Quantity, Quantity] | list[tuple[Quantity, Quantity, Quantity]] | dict[str, Any], *, fftlength: float | None = None, overlap: float | None = None, nfft: int | None = None, noverlap: int | None = None, window: str = 'hann', detrend: str | bool = 'constant', scaling: Literal['density', 'spectrum'] = 'density', average: Literal['mean', 'median'] = 'mean', return_onesided: bool = True, **kwargs) → FrequencySeries | FrequencySeriesList[source]

Compute power spectral density using Welch’s method.

Extracts time series from one or more spatial points/regions and estimates their PSD using scipy.signal.welch.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • point_or_region (tuple, list of tuples, or dict) –

    Spatial location(s) to extract:

    • Single point: (x, y, z) tuple of Quantities

    • Multiple points: list of (x, y, z) tuples

    • Region dict: {'x': slice or value, 'y': ..., 'z': ...} If region, averages over all points in the region.

  • fftlength (float, optional) – Segment length in seconds. Default is min(256, len(time_axis)) samples.

  • overlap (float, optional) – Overlap in seconds. Default is fftlength / 2.

  • window (str) – Window function name. Default is ‘hann’.

  • detrend (str or bool) – Detrending method: ‘constant’, ‘linear’, or False. Default ‘constant’.

  • scaling ({'density', 'spectrum'}) – ‘density’ for PSD (V²/Hz), ‘spectrum’ for power spectrum (V²). Default is ‘density’.

  • average ({'mean', 'median'}) – Averaging method for segments. Default is ‘mean’.

  • return_onesided (bool) – If True, return one-sided spectrum for real data. Default True.

  • **kwargs – Additional keyword arguments. Passing the removed nperseg or noverlap parameters will raise TypeError.

Returns:

  • Single point: FrequencySeries with PSD values

  • Multiple points: FrequencySeriesList

Return type:

FrequencySeries or FrequencySeriesList

Raises:

ValueError – If time axis is not regularly spaced or axis0_domain != ‘time’.

Examples

>>> from astropy import units as u
>>> psd = compute_psd(field, (1.0*u.m, 2.0*u.m, 0.0*u.m))
>>> psd.frequencies  # Frequency axis with units

>>> # Multiple points
>>> points = [(0*u.m, 0*u.m, 0*u.m), (1*u.m, 0*u.m, 0*u.m)]
>>> psd_list = compute_psd(field, points)
gwexpy.fields.freq_space_map(field: ScalarField, axis: str, at: dict[str, Quantity] | None = None, *, method: Literal['welch', 'fft'] = 'welch', fftlength: float | None = None, overlap: float | None = None, nfft: int | None = None, noverlap: int | None = None, window: str = 'hann', detrend: str | bool = 'constant', scaling: Literal['density', 'spectrum'] = 'density', **kwargs) → ScalarField[source]

Compute frequency-space map along a spatial axis.

Scans along the specified spatial axis, computing PSD at each position while fixing other spatial coordinates. Produces a 2D map with frequency on one axis and spatial coordinate on the other.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • axis (str) – Spatial axis to scan along (‘x’, ‘y’, or ‘z’).

  • at (dict, optional) – Fixed values for the other two spatial axes. Example: {'y': 0*u.m, 'z': 0*u.m} when axis=’x’.

  • method ({'welch', 'fft'}) – PSD estimation method. Default is ‘welch’.

  • fftlength (float, optional) – Segment length in seconds. Default is min(256, nt) samples.

  • overlap (float, optional) – Overlap in seconds. Default is fftlength / 2.

  • window (str) – Window function. Default is ‘hann’.

  • detrend (str or bool) – Detrending method. Default is ‘constant’.

  • scaling ({'density', 'spectrum'}) – PSD scaling. Default is ‘density’.

  • **kwargs – Additional keyword arguments. Passing the removed nperseg or noverlap parameters will raise TypeError.

Returns:

2D frequency-space map stored as ScalarField with: - axis0: frequency (with axis0_domain=’frequency’) - axis1: spatial coordinate - axis2, axis3: length-1 dummy axes

Return type:

ScalarField

Examples

>>> from astropy import units as u
>>> fsmap = freq_space_map(field, 'x', at={'y': 0*u.m, 'z': 0*u.m})
>>> fsmap.shape  # (n_freq, n_x, 1, 1)

See also

compute_psd

Single-point PSD computation.

compute_freq_space

Alias for this function.

gwexpy.fields.compute_freq_space(field: ScalarField, axis: str, at: dict[str, Quantity] | None = None, *, method: Literal['welch', 'fft'] = 'welch', fftlength: float | None = None, overlap: float | None = None, nfft: int | None = None, noverlap: int | None = None, window: str = 'hann', detrend: str | bool = 'constant', scaling: Literal['density', 'spectrum'] = 'density', **kwargs) → ScalarField

Compute frequency-space map along a spatial axis.

Scans along the specified spatial axis, computing PSD at each position while fixing other spatial coordinates. Produces a 2D map with frequency on one axis and spatial coordinate on the other.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • axis (str) – Spatial axis to scan along (‘x’, ‘y’, or ‘z’).

  • at (dict, optional) – Fixed values for the other two spatial axes. Example: {'y': 0*u.m, 'z': 0*u.m} when axis=’x’.

  • method ({'welch', 'fft'}) – PSD estimation method. Default is ‘welch’.

  • fftlength (float, optional) – Segment length in seconds. Default is min(256, nt) samples.

  • overlap (float, optional) – Overlap in seconds. Default is fftlength / 2.

  • window (str) – Window function. Default is ‘hann’.

  • detrend (str or bool) – Detrending method. Default is ‘constant’.

  • scaling ({'density', 'spectrum'}) – PSD scaling. Default is ‘density’.

  • **kwargs – Additional keyword arguments. Passing the removed nperseg or noverlap parameters will raise TypeError.

Returns:

2D frequency-space map stored as ScalarField with: - axis0: frequency (with axis0_domain=’frequency’) - axis1: spatial coordinate - axis2, axis3: length-1 dummy axes

Return type:

ScalarField

Examples

>>> from astropy import units as u
>>> fsmap = freq_space_map(field, 'x', at={'y': 0*u.m, 'z': 0*u.m})
>>> fsmap.shape  # (n_freq, n_x, 1, 1)

See also

compute_psd

Single-point PSD computation.

compute_freq_space

Alias for this function.

gwexpy.fields.compute_xcorr(field: ScalarField, point_a: tuple[Quantity, Quantity, Quantity], point_b: tuple[Quantity, Quantity, Quantity], *, max_lag: int | Quantity | None = None, mode: Literal['full', 'same', 'valid'] = 'full', normalize: bool = True, detrend: bool = True, window: str | None = None) → TimeSeries[source]

Compute cross-correlation between two spatial points.

Calculates the normalized (or unnormalized) cross-correlation function between time series extracted at two spatial locations.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • point_a (tuple of Quantity) – Spatial coordinates (x, y, z) for the two points.

  • point_b (tuple of Quantity) – Spatial coordinates (x, y, z) for the two points.

  • max_lag (int or Quantity, optional) – Maximum lag to compute. If Quantity, converted to samples. If None, uses all available lags.

  • mode ({'full', 'same', 'valid'}) – Correlation mode (see numpy.correlate). Default is ‘full’.

  • normalize (bool) – If True, normalize to [-1, 1] range. Default True.

  • detrend (bool) – If True, remove mean before correlation. Default True.

  • window (str, optional) – Window function to apply before correlation. Default None.

Returns:

Cross-correlation function with lag axis in time units. Positive lag means point_b leads point_a.

Return type:

TimeSeries

Examples

>>> from astropy import units as u
>>> xcorr = compute_xcorr(field, (0*u.m, 0*u.m, 0*u.m), (1*u.m, 0*u.m, 0*u.m))
>>> xcorr.times  # Lag axis
>>> xcorr.value  # Correlation values

Notes

The correlation is computed as correlate(a, b, mode), where positive values in the output correspond to b leading a.

gwexpy.fields.time_delay_map(field: ScalarField, ref_point: tuple[Quantity, Quantity, Quantity], plane: str = 'xy', at: dict[str, Quantity] | None = None, *, max_lag: int | Quantity | None = None, stride: int | Quantity = 1, roi: dict[str, slice] | None = None, normalize: bool = True, detrend: bool = True) → ScalarField[source]

Compute time delay map from a reference point to a 2D slice.

For each point in the specified plane, computes cross-correlation with the reference point and extracts the lag at maximum correlation.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • ref_point (tuple of Quantity) – Reference point coordinates (x, y, z).

  • plane (str) – 2D plane to map: ‘xy’, ‘xz’, or ‘yz’. Default ‘xy’.

  • at (dict, optional) – Fixed value for the axis not in the plane. Example: {'z': 0*u.m} when plane=’xy’.

  • max_lag (int or Quantity, optional) – Maximum lag to search for peak. Default uses all lags.

  • stride (int or Quantity, optional) – Spatial subsample stride. If int, treated as grid point stride. If Quantity, treated as physical distance stride and converted to grid points based on spatial resolution. Default 1.

  • roi (dict, optional) – Region of interest as dict of slices. Example: {'x': slice(10, 50), 'y': slice(20, 80)}

  • normalize (bool) – Normalize correlation. Default True.

  • detrend (bool) – Detrend before correlation. Default True.

Returns:

Time delay map as ScalarField slice with delay values in time units. Shape matches the input plane dimensions (with stride applied).

Return type:

ScalarField

Examples

>>> from astropy import units as u
>>> delay_map = time_delay_map(
...     field,
...     ref_point=(0*u.m, 0*u.m, 0*u.m),
...     plane='xy',
...     at={'z': 0*u.m},
...     stride=0.05*u.m  # spatial stride in meters
... )

Notes

Positive delay means the point leads the reference.

gwexpy.fields.coherence_map(field: ScalarField, ref_point: tuple[Quantity, Quantity, Quantity], *, plane: str = 'xy', at: dict[str, Quantity] | None = None, band: tuple[Quantity, Quantity] | None = None, fftlength: float | None = None, overlap: float | None = None, nfft: int | None = None, noverlap: int | None = None, window: str = 'hann', stride: int | Quantity = 1, **kwargs) → ScalarField | FieldDict[source]

Compute magnitude-squared coherence map from a reference point.

Calculates coherence between the reference time series and all points in the specified 2D slice.

Parameters:

  • field (ScalarField) – Input 4D field with axis0_domain=’time’.

  • ref_point (tuple of Quantity) – Reference point coordinates (x, y, z).

  • plane (str) – 2D plane to map: ‘xy’, ‘xz’, or ‘yz’. Default ‘xy’.

  • at (dict, optional) – Fixed value for the axis not in the plane.

  • band (tuple of Quantity, optional) – Frequency band (fmin, fmax) to average coherence over. If None, returns coherence at all frequencies.

  • fftlength (float, optional) – Segment length in seconds. Default min(256, nt) samples. Cannot be used with nfft.

  • overlap (float, optional) – Overlap in seconds. Default fftlength / 2. Cannot be used with noverlap.

  • nfft (int, optional) – Number of samples per segment (alternative to fftlength). Cannot be used with fftlength.

  • noverlap (int, optional) – Number of overlapping samples (alternative to overlap). Cannot be used with overlap.

  • window (str) – Window function. Default ‘hann’.

  • stride (int or Quantity, optional) – Spatial subsample stride. If int, treated as grid point stride. If Quantity, treated as physical distance stride and converted to grid points based on spatial resolution. Default 1.

  • **kwargs – Additional keyword arguments.

Returns:

  • If band is specified: ScalarField with scalar coherence values (0-1).

  • If band is None: FieldDict with keys as frequency indices, or a single ScalarField with frequency as axis0.

Return type:

ScalarField or FieldDict

Examples

>>> from astropy import units as u
>>> # Compute coherence map with time-based FFT parameters
>>> coh_map = coherence_map(
...     field,
...     ref_point=(0*u.m, 0*u.m, 0*u.m),
...     plane='xy',
...     at={'z': 0*u.m},
...     band=(10*u.Hz, 100*u.Hz),
...     fftlength=1.0,    # 1.0 seconds
...     overlap=0.5,      # 0.5 seconds
...     window='hann'
... )

>>> # Or with sample-based FFT parameters
>>> coh_map = coherence_map(
...     field,
...     ref_point=(0*u.m, 0*u.m, 0*u.m),
...     nfft=1024,        # samples
...     noverlap=512,     # samples
...     stride=0.1*u.m    # spatial stride in meters
... )
class gwexpy.fields.FieldList(items=None, validate=False)[source]

Bases: list

List-like collection for ScalarField objects with batch operations.

fft_time_all(**kwargs)[source]

Apply fft_time to all fields, returning FieldList.

ifft_time_all(**kwargs)[source]

Apply ifft_time to all fields, returning FieldList.

fft_space_all(axes=None, **kwargs)[source]

Apply fft_space to all fields, returning FieldList.

ifft_space_all(axes=None, **kwargs)[source]

Apply ifft_space to all fields, returning FieldList.

resample_all(rate, **kwargs)[source]

Apply resample to all fields, returning FieldList.

filter_all(*args, **kwargs)[source]

Apply filter to all fields, returning FieldList.

sel_all(**kwargs)[source]

Apply sel to all fields, returning FieldList.

isel_all(**kwargs)[source]

Apply isel to all fields, returning FieldList.

class gwexpy.fields.FieldDict(items=None, validate=False)[source]

Bases: dict

Dict-like collection for ScalarField objects with batch operations.

copy() → FieldDict[source]

Return a copy of this FieldDict.

fft_time_all(**kwargs)[source]

Apply fft_time to all fields, returning FieldDict.

ifft_time_all(**kwargs)[source]

Apply ifft_time to all fields, returning FieldDict.

fft_space_all(axes=None, **kwargs)[source]

Apply fft_space to all fields, returning FieldDict.

ifft_space_all(axes=None, **kwargs)[source]

Apply ifft_space to all fields, returning FieldDict.

resample_all(rate, **kwargs)[source]

Apply resample to all fields, returning FieldDict.

filter_all(*args, **kwargs)[source]

Apply filter to all fields, returning FieldDict.

sel_all(**kwargs)[source]

Apply sel to all fields, returning FieldDict.

isel_all(**kwargs)[source]

Apply isel to all fields, returning FieldDict.

class gwexpy.fields.VectorField(components=None, basis: Literal['cartesian', 'custom'] = 'cartesian', validate: bool = True)[source]

Bases: FieldDict

Vector-valued field as a collection of ScalarField components.

This class maintains a collection of ScalarField objects, constrained to have identical axis metadata. This structure allows performing physical operations (FFT, filtering) on each component independently while ensuring geometrical consistency.

Parameters:

  • components (dict[str, ScalarField], optional) – Dictionary mapping component names (e.g., ‘x’, ‘y’) to ScalarField objects.

  • basis ({'cartesian', 'custom'}, optional) – The coordinate basis of the vector. Default is ‘cartesian’.

  • validate (bool, optional) – Whether to validate that all components have consistent axes and units. Default is True.

Examples

>>> from gwexpy.fields import ScalarField, VectorField
>>> fx = ScalarField(np.random.randn(100, 4, 4, 4))
>>> fy = ScalarField(np.random.randn(100, 4, 4, 4))
>>> vf = VectorField({'x': fx, 'y': fy})
>>> vnorm = vf.norm()
copy() → VectorField[source]

Return a copy of this VectorField.

to_array() → ndarray[source]

Convert components to a single 5D ndarray.

The last axis (axis 4) corresponds to the vector components.

Returns:

A 5D array of shape (axis0, axis1, axis2, axis3, n_components).

Return type:

numpy.ndarray

norm() → ScalarField[source]

Compute the L2 norm of the vector field.

Returns:

A scalar field representing the magnitude \(\sqrt{\sum |E_i|^2}\).

Return type:

ScalarField

dot(other: VectorField) → ScalarField[source]

Compute the scalar (dot) product with another VectorField.

Parameters:

other (VectorField) – The other vector field to dot with. Must have identical components.

Returns:

The resulting dot product field.

Return type:

ScalarField

cross(other: VectorField) → VectorField[source]

Compute the cross product with another VectorField (3-components).

Parameters:

other (VectorField) – The other vector field. Both must have exactly ‘x’, ‘y’, ‘z’ components.

Returns:

The resulting vector field.

Return type:

VectorField

project(direction: VectorField) → ScalarField[source]

Project the vector field onto a given direction.

Result is (self . direction) / norm(direction).

Parameters:

direction (VectorField) – The direction to project onto.

Returns:

The projected scalar field.

Return type:

ScalarField

plot_magnitude(x=None, y=None, **kwargs)[source]

Plot the magnitude (norm) of the vector field.

Parameters:

  • x (str, optional) – X-axis name.

  • y (str, optional) – Y-axis name.

  • **kwargs – Fixed coordinates and arguments passed to ScalarField.plot.

quiver(x=None, y=None, **kwargs)[source]

Plot the vector field using arrows (quiver).

Parameters:

  • x (str, optional) – X-axis name.

  • y (str, optional) – Y-axis name.

  • **kwargs – Fixed coordinates and arguments passed to FieldPlot.add_vector.

streamline(x=None, y=None, **kwargs)[source]

Plot the vector field using streamlines.

Parameters:

  • x (str, optional) – X-axis name.

  • y (str, optional) – Y-axis name.

  • **kwargs – Fixed coordinates and arguments passed to FieldPlot.add_vector.

plot(x=None, y=None, **kwargs)[source]

Plot magnitude and overlay quiver.

This is a convenience method combining plot_magnitude and quiver.

class gwexpy.fields.TensorField(components: dict[tuple[int, ...], ScalarField] | None = None, rank: int | None = None, validate: bool = True)[source]

Bases: FieldDict

Tensor-valued field as a collection of ScalarField components.

Keys are typically tuples of indices representing the tensor components, e.g., (0, 0) for the \(T_{00}\) component.

Parameters:

  • components (dict[tuple, ScalarField], optional) – Dictionary mapping component index tuples to ScalarField objects.

  • rank (int, optional) – The rank (order) of the tensor. If not provided, it is inferred from the keys.

  • validate (bool, optional) – Whether to validate component consistency. Default is True.

copy() → TensorField[source]

Return a copy of this TensorField.

trace() → ScalarField[source]

Compute the trace of a rank-2 tensor.

Returns:

The sum of diagonal components \(\sum T_{ii}\).

Return type:

ScalarField

Raises:

ValueError – If the tensor rank is not 2.

to_array(order: Literal['first', 'last'] = 'last') → ndarray[source]

Convert components into a single NumPy array.

Parameters:

order ({'first', 'last'}, optional) – If ‘first’ (default for VectorField), the tensor dimensions are at the beginning. If ‘last’ (preferred for linalg), they are at the end.

Returns:

Shape (dim1, dim2, N0, N1, N2, N3) if order=’first’. Shape (N0, N1, N2, N3, dim1, dim2) if order=’last’.

Return type:

ndarray

det() → ScalarField[source]

Compute the determinant of a rank-2 tensor.

Returns:

The determinant field.

Return type:

ScalarField

symmetrize() → TensorField[source]

Symmetrize a rank-2 tensor: \(S_{ij} = (T_{ij} + T_{ji}) / 2\).

Returns:

The symmetrized tensor.

Return type:

TensorField

Raises:

ValueError – If the tensor rank is not 2.

plot_components(x=None, y=None, **kwargs)[source]

Plot the tensor components in a grid layout.

Parameters:

  • x (str, optional) – X-axis name.

  • y (str, optional) – Y-axis name.

  • **kwargs – Fixed coordinates and arguments passed to FieldPlot.add_scalar.