from __future__ import annotations
import hashlib
import math
from dataclasses import asdict, dataclass, field, fields
from pathlib import Path
from typing import TYPE_CHECKING, Any
import numpy as np
import xarray as xr
from geopandas import GeoSeries
from numpy.random import Generator
from shapely.ops import transform
from tqdm import tqdm
from vtk import vtkPolyData
from cratermaker.constants import _VBIG, PairOfFloats
from cratermaker.core.base import ComponentBase, CratermakerBase, import_components
from cratermaker.utils.general_utils import format_large_units, validate_and_normalize_location
if TYPE_CHECKING:
from cratermaker.components.projectile import Projectile
from cratermaker.components.scaling import Scaling
from cratermaker.components.surface import Surface
from cratermaker.components.target import Target
_TALLY_LONG_NAME = "Unique crater identification number"
[docs]
@dataclass(frozen=True, slots=True)
class CraterFixed:
id: np.uint32 = None
"""A unique identifyer for the crater, which can be used to track it across time steps in a simulation. This is automatically generated when a crater is emplaced in a simulation as a hash of its fixed attributes."""
semimajor_axis: float | None = None
"""The semimajor axis of the crater in meters. For a circular crater, this is equal to the semiminor axis and the radius."""
semiminor_axis: float | None = None
"""The semiminor axis of the crater in meters. For a circular crater, this is equal to the semimajor axis and the radius."""
orientation: float | None = None
"""Orientation of the crater in degrees, measured clockwise from north. For a circular crater, this value is not meaningful, but will be set to the same value as `projectile_orientation`."""
transient_diameter: float | None = None
"""The transient diameter of the crater in meters."""
projectile_diameter: float | None = None
"""The diameter of the projectile in meters."""
projectile_velocity: float | None = None
"""The velocity of the projectile in meters per second."""
projectile_angle: float | None = None
"""The impact angle of the projectile relative to the surfacein degrees."""
projectile_mass: float | None = None
"""The mass of the projectile in kg."""
location: PairOfFloats | None = None
"""The location of the crater center in (longitude, latitude) in degrees."""
morphology_type: str | None = None
"""The morphology type of the crater, i.e. "simple", "complex", etc."""
time: float | None = None
"""The time in My before present at which the crater was emplaced or observed."""
radius: float | None = field(default=None, init=False)
"""The radius of the crater in meters. For a circular crater, this is equal to the semimajor and semiminor axes. For a non-circular crater, this is the geometric mean of the semimajor and semiminor axes."""
diameter: float | None = field(default=None, init=False)
"""The diameter of the crater in meters. For a circular crater, this is equal to twice the semimajor and semiminor axes. For a non-circular crater, this is twice the geometric mean of the semimajor and semiminor axes."""
transient_radius: float | None = field(default=None, init=False)
"""The transient radius of the crater in meters, which is half the transient diameter."""
projectile_radius: float | None = field(default=None, init=False)
"""The radius of the projectile in meters, which is half the projectile diameter."""
projectile_density: float | None = field(default=None, init=False)
"""The density of the projectile in kg/m³."""
projectile_vertical_velocity: float | None = field(default=None, init=False)
"""The vertical component of the projectile velocity in meters per second."""
projectile_direction: float | None = field(default=None, init=False)
"""The direction of the projectile in degrees, measured clockwise from north. This is the same as `orientation`."""
def __post_init__(self):
# Compute derived attributes that depend on other fields and cannot be set directly
if self.semimajor_axis is not None and self.semiminor_axis is not None:
object.__setattr__(self, "radius", math.sqrt(self.semimajor_axis * self.semiminor_axis))
object.__setattr__(self, "diameter", self.radius * 2)
if self.transient_diameter is not None:
object.__setattr__(self, "transient_radius", self.transient_diameter / 2.0)
if self.projectile_diameter is not None:
object.__setattr__(self, "projectile_radius", self.projectile_diameter / 2.0)
if self.projectile_mass is not None:
object.__setattr__(
self, "projectile_density", (3.0 * self.projectile_mass) / (4.0 * math.pi * self.projectile_radius**3)
)
if self.projectile_velocity is not None and self.projectile_angle is not None:
object.__setattr__(
self, "projectile_vertical_velocity", self.projectile_velocity * math.sin(math.radians(self.projectile_angle))
)
if self.orientation is not None:
object.__setattr__(self, "projectile_direction", self.orientation)
[docs]
class CraterVariable:
def __init__(
self,
measured_diameter: float | None = None,
measured_radius: float | None = None,
measured_semimajor_axis: float | None = None,
measured_semiminor_axis: float | None = None,
measured_orientation: float | None = None,
measured_location: PairOfFloats | None = None,
measured_rim_height: float | None = None,
measured_floor_depth: float | None = None,
degradation_state: float | None = None,
production_time: tuple[float, float] | float | None = None,
production_ND: tuple[float, float, float] | tuple[float, float] | None = None,
production_sequence: int | None = None,
name: str | None = None,
**kwargs: Any,
):
object.__setattr__(self, "_measured_semimajor_axis", None)
object.__setattr__(self, "_measured_semiminor_axis", None)
object.__setattr__(self, "_measured_radius", None)
object.__setattr__(self, "_measured_orientation", None)
object.__setattr__(self, "_measured_location", None)
object.__setattr__(self, "_measured_rim_height", None)
object.__setattr__(self, "_measured_floor_depth", None)
object.__setattr__(self, "_degradation_state", None)
object.__setattr__(self, "_production_time", None)
object.__setattr__(self, "_production_ND", None)
object.__setattr__(self, "_production_sequence", None)
object.__setattr__(self, "_name", None)
if measured_diameter is not None:
self.measured_diameter = measured_diameter
if measured_radius is not None:
self.measured_radius = measured_radius
if measured_semimajor_axis is not None:
self.measured_semimajor_axis = measured_semimajor_axis
if measured_semiminor_axis is not None:
self.measured_semiminor_axis = measured_semiminor_axis
if measured_orientation is not None:
self.measured_orientation = measured_orientation
if measured_location is not None:
self.measured_location = measured_location
if measured_rim_height is not None:
self.measured_rim_height = measured_rim_height
if measured_floor_depth is not None:
self.measured_floor_depth = measured_floor_depth
if degradation_state is not None:
self.degradation_state = degradation_state
if production_time is not None:
self.production_time = production_time
if production_ND is not None:
self.production_ND = production_ND
if production_sequence is not None:
self.production_sequence = production_sequence
if name is not None:
self.name = name
return
def __repr__(self):
return (
f"CraterVariable(measured_diameter={self.measured_diameter}, "
f"measured_radius={self.measured_radius}, "
f"measured_semimajor_axis={self.measured_semimajor_axis}, "
f"measured_semiminor_axis={self.measured_semiminor_axis}, "
f"measured_orientation={self.measured_orientation}, "
f"measured_location={self.measured_location}, "
f"measured_rim_height={self.measured_rim_height}, "
f"measured_floor_depth={self.measured_floor_depth}, "
f"degradation_state={self.degradation_state}, "
f"production_time={self.production_time}, "
f"production_ND={self.production_ND}, "
f"production_sequence={self.production_sequence}, "
f"name={self.name} "
)
[docs]
def as_dict(self):
"""
Convert the variable attributes to a dictionary for serialization, excluding any attributes that are None.
If both measured_diameter and measured_semimajor_axis are provided and consistent, only include measured_diameter. If both measured_radius and measured_semimajor_axis are provided and consistent, only include measured_radius. This is to avoid redundancy in the output.
"""
dict_repr = {
"measured_orientation": self.measured_orientation,
"measured_location": self.measured_location,
"measured_rim_height": self.measured_rim_height,
"measured_floor_depth": self.measured_floor_depth,
"degradation_state": self.degradation_state,
}
if self.measured_semimajor_axis is not None:
if self.measured_radius is not None and self.measured_radius == self.measured_semimajor_axis:
dict_repr["measured_diameter"] = self.measured_diameter
else:
dict_repr["measured_semimajor_axis"] = self.measured_semimajor_axis
dict_repr["measured_semiminor_axis"] = self.measured_semiminor_axis
if self.production_time is not None:
dict_repr["production_time"] = self.production_time
if self.production_ND is not None:
dict_repr["production_ND"] = self.production_ND
if self.name is not None:
dict_repr["name"] = self.name
if self.production_sequence is not None:
dict_repr["production_sequence"] = self.production_sequence
return dict_repr
@property
def name(self) -> str | None:
"""Optional name of the crater, which does not need to be unique."""
return self._name
@name.setter
def name(self, value):
if value is not None:
self._name = str(value)
@property
def measured_diameter(self) -> float | None:
"""Final diameter of the crater in meters."""
return self.measured_radius * 2 if self.measured_radius is not None else None
@measured_diameter.setter
def measured_diameter(self, value: float | None):
if value is None:
self.measured_radius = None
else:
self.measured_radius = value / 2.0
return
@property
def measured_radius(self) -> float | None:
"""
The measured radius of the crater in meters. For a circular crater, this is equal to the measured semimajor and semiminor axes.
For a non-circular crater, this is the geometric mean of the measured semimajor and measured semiminor axes. This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed semimajor axis attribute.
"""
if self.measured_semimajor_axis is not None and self.measured_semiminor_axis is not None:
return math.sqrt(self.measured_semimajor_axis * self.measured_semiminor_axis)
return None
@measured_radius.setter
def measured_radius(self, value: float | None):
if value is None:
self.measured_semimajor_axis = None
self.measured_semiminor_axis = None
else:
self.measured_semimajor_axis = value
self.measured_semiminor_axis = value
return
@property
def measured_semimajor_axis(self) -> float | None:
"""
The measured semimajor axis of the crater in meters.
For a circular crater, this is equal to the measured semiminor axis and the measured radius. For a non-circular crater, this is the longest axis of the crater. This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed semimajor axis attribute.
"""
return self._measured_semimajor_axis
@measured_semimajor_axis.setter
def measured_semimajor_axis(self, value: float | None):
if value is None:
self._measured_semimajor_axis = None
return
if value <= 0:
raise ValueError("measured_semimajor_axis must be positive.")
if self._measured_semiminor_axis is not None and value < self._measured_semiminor_axis:
self._measured_semimajor_axis = self._measured_semiminor_axis
self._measured_semiminor_axis = float(value)
else:
self._measured_semimajor_axis = float(value)
return
@property
def measured_semiminor_axis(self) -> float | None:
"""
The measured semiminor axis of the crater in meters.
For a circular crater, this is equal to the measured semimajor axis and the measured radius. This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed semimajor axis attribute.
"""
return self._measured_semiminor_axis
@measured_semiminor_axis.setter
def measured_semiminor_axis(self, value: float | None):
"""
The measured semiminor axis of the crater in meters.
For a circular crater, this is equal to the measured semimajor axis and the measured radius. This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed semimajor axis attribute.
"""
if value is None:
self._measured_semiminor_axis = None
return
if value <= 0:
raise ValueError("measured_semiminor_axis must be positive.")
if self._measured_semimajor_axis is not None and value > self._measured_semimajor_axis:
self._measured_semiminor_axis = self._measured_semimajor_axis
self._measured_semimajor_axis = float(value)
else:
self._measured_semiminor_axis = float(value)
return
@property
def measured_orientation(self) -> float | None:
"""
The measured orientation of the crater in degrees, measured clockwise from north.
This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed orientation attribute.
"""
return self._measured_orientation
@measured_orientation.setter
def measured_orientation(self, value: float | None):
if value is None:
self._measured_orientation = None
return
self._measured_orientation = float(value) % 360
return
@property
def measured_location(self) -> tuple[float, float] | None:
return self._measured_location
@measured_location.setter
def measured_location(self, value: tuple[float, float] | None):
"""
The measured location of the crater as a tuple of (longitude, latitude) in degrees.
This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of the fixed location attribute.
"""
if value is None:
self._measured_location = None
return
self._measured_location = validate_and_normalize_location(value)
return
@property
def measured_rim_height(self) -> float | None:
return self._measured_rim_height
@measured_rim_height.setter
def measured_rim_height(self, value: float | None):
"""
The measured rim height above the local reference plane of the crater in meters.
This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of any fixed attributes.
"""
if value is None:
self._measured_rim_height = None
return
self._measured_rim_height = float(value)
return
@property
def measured_floor_depth(self) -> float | None:
"""
The measured floor depth below the local reference plane of the crater in meters.
This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of any fixed attributes.
"""
return self._measured_floor_depth
@measured_floor_depth.setter
def measured_floor_depth(self, value: float | None):
if value is None:
self._measured_floor_depth = None
return
self._measured_floor_depth = float(value)
return
@property
def degradation_state(self) -> float | None:
"""The degradation state of the crater in m², which is a measure of how much the crater has degraded diffusively from its original state.
This is a variable attribute that can be updated as the crater degrades or as measurements are refined, and is independent of any fixed attributes.
"""
return self._degradation_state
@degradation_state.setter
def degradation_state(self, value: float | None):
if value is None:
self._degradation_state = None
return
self._degradation_state = float(value)
return
@property
def production_time(self) -> float | None:
"""The range of ages of the crater in Myr before present, used by the quasi-monte carlo sampling method to emplace a user-defined crater within a time period."""
return self._production_time
@production_time.setter
def production_time(self, value: tuple[float, float] | float | None):
if value is not None:
if isinstance(value, (int, float)):
value = (value, 0)
if len(value) != 2:
raise ValueError("production_time must be a single float or a tuple of (time, time_stdev).")
if value[0] is None:
return
self._production_time = float(value[0]), float(value[1])
return
@property
def production_ND(self) -> tuple[float, float, float] | tuple[float, float] | None:
"""A tuple of diameter and cumulative number values, in the form of a (D, N), or (D, N, N_stdev) used by the quasi-monte carlo sampling method to emplace a user-defined crater within a number-diameter range."""
return self._production_ND
@production_ND.setter
def production_ND(self, value: tuple[float, float, float] | tuple[float, float] | None):
if value is not None:
if len(value) != 3 and len(value) != 2:
raise ValueError("production_ND must be a tuple of floats in the form of either (D, N) or (D, N, N_stdev).")
if value[0] is None:
return
if len(value) == 2:
self._production_ND = (float(value[0]), float(value[1]), 0.0)
else:
self._production_ND = (float(value[0]), float(value[1]), float(value[2]))
return
@property
def production_sequence(self) -> int | None:
"""The production sequence number of the crater, used by the quasi-monte carlo sampling method to emplace a user-defined crater in a relative sequence order with other quasi-monte carlo sampled craters."""
return self._production_sequence
@production_sequence.setter
def production_sequence(self, value: int | None):
if value is not None:
self._production_sequence = int(value)
return
[docs]
class Crater(ComponentBase):
_registry: dict[str, Crater] = {}
def __init__(
self,
crater: Crater | None = None,
fixed_cls: type[CraterFixed] = CraterFixed,
variable_cls: type[CraterVariable] = CraterVariable,
**kwargs: Any,
):
if crater is not None:
if not isinstance(crater, Crater):
raise TypeError("crater must be an instance of Crater or None")
fixed_fields = asdict(crater._fixed)
for f in fields(crater._fixed):
if not f.init:
fixed_fields.pop(f.name)
var_fields = crater._var.as_dict()
# Be sure to scrub any fields that are redundant with any potential kwargs
kwargs = self._scrub_measured_input_size(**kwargs)
if "measured_semimajor_axis" in kwargs or "measured_semiminor_axis" in kwargs:
var_fields.pop("measured_diameter", None)
var_fields.pop("measured_radius", None)
else:
fixed_fields = {}
var_fields = {}
for f in fields(fixed_cls):
if f.name in kwargs and f.init and (f.name not in fixed_fields or fixed_fields[f.name] is None):
fixed_fields[f.name] = kwargs.pop(f.name)
if f.name not in fixed_fields and f.name in kwargs and f.init:
fixed_fields[f.name] = kwargs.pop(f.name)
var_fields.update({k: v for k, v in kwargs.items() if k not in fixed_fields})
self._fixed = fixed_cls(**fixed_fields)
self._var = variable_cls(**var_fields)
return
def __getattr__(self, name: str):
"""
Delegates attribute access to the fixed and variable fields.
If an attribute is not found in the Crater class, it will first check the fixed fields, then the variable fields.
"""
if hasattr(self._fixed, name):
return getattr(self._fixed, name)
if hasattr(self._var, name):
return getattr(self._var, name)
raise AttributeError(f"{type(self).__name__} has no attribute {name!r}")
def __setattr__(self, name: str, value: Any):
"""
Delegates attribute setting to the variable fields.
Only attributes in the variable fields can be set. Fixed fields are immutable.
"""
if name in ("_fixed", "_var", "_fixed_fields", "_var_fields"):
object.__setattr__(self, name, value)
return
if name in CraterFixed.__dataclass_fields__:
raise AttributeError(f"{name!r} is a fixed attribute and cannot be modified.")
# Call the setters in CraterVariable if they exist
if hasattr(self._var, name):
setattr(self._var, name, value)
return
[docs]
def as_dict(self, ignore_keys: list[str] = [], **kwargs) -> dict:
"""
Expose delegated attributes for serialization.
"""
base = {}
base.update(asdict(self._fixed))
base.update(self._var.as_dict())
if ignore_keys is not None:
for key in ignore_keys:
base.pop(key, None)
return base
def __dir__(self):
"""
Expose delegated attributes for autocomplete/introspection.
"""
base = set(super().__dir__())
base.update(dir(self.fixed))
base.update(dir(self.var))
return sorted(base)
def __str__(self):
str_repr = super().__str__()
str_repr += f"ID: {self.id}\n"
if self.name is not None and self.name != "":
str_repr += f"Name: {self.name}\n"
if self.semimajor_axis != self.semiminor_axis:
str_repr += f"Elliptical size {format_large_units(self.semimajor_axis, quantity='length')} x {format_large_units(self.semiminor_axis, quantity='length')}\nMean Diameter: {format_large_units(self.diameter, quantity='length')}\n"
else:
str_repr += f"Diameter: {format_large_units(self.diameter, quantity='length')}\n"
if self.measured_semimajor_axis is not None and self.measured_semimajor_axis != self.semimajor_axis:
if self.measured_semimajor_axis != self.measured_semiminor_axis:
str_repr += f"Measured elliptical size {format_large_units(self.measured_semimajor_axis, quantity='length')} x {format_large_units(self.measured_semiminor_axis, quantity='length')}\nMeasured Mean Diameter: {format_large_units(self.measured_diameter, quantity='length')}\n"
else:
str_repr += f"Measured diameter: {format_large_units(self.measured_diameter, quantity='length')}\n"
str_repr += f"morphology_type: {self.morphology_type}\n"
str_repr += f"transient_diameter: {format_large_units(self.transient_diameter, quantity='length')}\n"
str_repr += f"projectile_diameter: {format_large_units(self.projectile_diameter, quantity='length')}\n"
str_repr += f"projectile_mass: {self.projectile_mass:.4e} kg\n"
str_repr += f"projectile_density: {self.projectile_density:.0f} kg/m³\n"
str_repr += f"projectile_velocity: {format_large_units(self.projectile_velocity, quantity='velocity')}\n"
str_repr += f"projectile_angle: {self.projectile_angle:.1f}°\n"
str_repr += f"projectile_direction: {self.projectile_direction:.1f}°\n"
if self.measured_orientation is not None and self.measured_orientation != self.orientation:
str_repr += f"Measured orientation: {self.measured_orientation:.1f}°\n"
str_repr += f"location (lon,lat): ({self.location[0]:.4f}°, {self.location[1]:.4f}°)\n"
if self.measured_location is not None and self.measured_location != self.location:
str_repr += f"Measured location (lon, lat): ({self.measured_location[0]:.4f}°, {self.measured_location[1]:.4f}°)\n"
if self.measured_rim_height is not None and np.abs(self.measured_rim_height) < _VBIG:
str_repr += f"Measured rim height: {format_large_units(self.measured_rim_height, quantity='length')}\n"
if self.measured_floor_depth is not None and np.abs(self.measured_floor_depth) < _VBIG:
str_repr += f"Measured floor depth: {format_large_units(self.measured_floor_depth, quantity='length')}\n"
if self.production_sequence is not None:
str_repr += f"Production sequence: {self.production_sequence}\n"
if self.production_ND is not None:
str_repr += f"Production N({self.production_ND[0]:.0f}): {self.production_ND[1]:.4f} ± {self.production_ND[2]:.4f}\n"
if self.production_time is not None:
str_repr += f"Production time: {format_large_units(self.production_time[0], quantity='time')} ± {format_large_units(self.production_time[1], quantity='time')}\n"
if self.time is not None:
str_repr += f"Emplacement time: {format_large_units(self.time, quantity='time')}\n"
if self.degradation_state is not None:
str_repr += f"Degradation state: {format_large_units(self.degradation_state, quantity='area')}\n"
return str_repr
def _scrub_measured_input_size(self, **kwargs: Any):
if "measured_diameter" in kwargs:
measured_diameter = kwargs.pop("measured_diameter")
kwargs["measured_semimajor_axis"] = measured_diameter / 2.0
kwargs["measured_semiminor_axis"] = measured_diameter / 2.0
elif "measured_radius" in kwargs:
measured_radius = kwargs.pop("measured_radius")
kwargs["measured_semimajor_axis"] = measured_radius
kwargs["measured_semiminor_axis"] = measured_radius
return kwargs
[docs]
@classmethod
def maker(
cls: type[Crater],
crater: Crater | None = None,
scaling: str | Scaling | None = None,
target: str | Target | None = None,
projectile: str | Projectile | None = None,
diameter: float | None = None,
radius: float | None = None,
final_diameter: float | None = None,
final_radius: float | None = None,
transient_diameter: float | None = None,
transient_radius: float | None = None,
semimajor_axis: float | None = None,
semiminor_axis: float | None = None,
orientation: float | None = None,
projectile_diameter: float | None = None,
projectile_radius: float | None = None,
projectile_mass: float | None = None,
projectile_density=None,
projectile_velocity: float | None = None,
projectile_mean_velocity: float | None = None,
projectile_vertical_velocity=None,
projectile_angle: float | None = None,
projectile_direction: float | None = None,
projectile_location: PairOfFloats | None = None,
location: PairOfFloats | None = None,
longitude: float | None = None,
latitude: float | None = None,
time: float | None = None,
name: str | None = None,
measured_semimajor_axis: float | None = None,
measured_semiminor_axis: float | None = None,
measured_orientation: float | None = None,
measured_diameter: float | None = None,
measured_radius: float | None = None,
measured_location: PairOfFloats | None = None,
measured_rim_height: float | None = None,
measured_floor_depth: float | None = None,
degradation_state: float | None = None,
production_time: tuple[float, float] | float | None = None,
production_ND: tuple[float, float, float] | tuple[float, float] | None = None,
production_sequence: int | None = None,
simdir: str | Path | None = None,
rng: Generator = None,
rng_seed: str | int | None = None,
rng_state: dict | None = None,
check_redundant_inputs: bool = True,
**kwargs: Any,
):
"""
Create a Crater object with the given parameters.
Parameters
----------
crater : Crater, optional
A Crater object to copy parameters from.
scaling : str or Scaling, optional
A string key or instance of a scaling model. If none provided, a default will be used
target : str or Target, optional
A string key or instance of a target model. If none provided, a default will be used.
projectile : str or Projectile, optional
A string key or instance of a projectile model. If none provided, a default will be used.
diameter : float, optional
The final diameter of the crater in meters.
final_diameter : float, optional
The final diameter of the crater in meters. This is an alias of "diameter"
radius : float, optional
The final radius of the crater in meters.
final_radius : float, optional
The final radius of the crater in meters. This is an alias of "radius"
transient_diameter : float, optional
The transient diameter of the crater in meters.
transient_radius : float, optional
The transient radius of the crater in meters.
semimajor_axis : float, optional
The semimajor axis of the crater in meters. Same as radius if the crater is circular.
semiminor_axis : float, optional
The semiminor axis of the crater in meters. Same as radius if the crater is circular.
orientation : float, optional
The orientation of the crater in degrees if it is elliptical.
projectile_diameter : float, optional
The diameter of the projectile in meters.
projectile_radius : float, optional
The radius of the projectile in meters.
projectile_mass : float, optional
The mass of the projectile in kilograms.
projectile_density : float, optional
The density of the projectile in kg/m³. If not provided, it will be defined through the projectile population model provided.
projectile_velocity : float, optional
The total impact velocity of the projectile in m/s.
projectile_mean_velocity : float, optional
The mean velocity from which to sample a projectile velocity.
projectile_vertical_velocity : float, optional
The vertical component of the velocity in m/s.
projectile_angle : float, optional
The impact angle in degrees (0-90).
projectile_direction : float, optional
The direction of the impact in degrees (0-360).
projectile_location : pair of float, optional
The (longitude, latitude) location of the projectile impact. This is equivalent to `location`, which takes precedence
location : pair of floats, optional
The (longitude, latitude) location of the crater on the target surface in degrees.
longitude : float, optional
The longitude of the crater on the target surface in degrees. If used, you must also pass `latitude` and cannot pass `location`.
latitude : float, optional
The latitude of the crater on the target surface in degrees. If used, you must also pass `longitude` and cannot pass `location`.
time : float, optional
The time of the crater impact in Myr before present.
measured_semimajor_axis : float, optional
The measured semimajor axis of the crater in meters.
measured_semiminor_axis : float, optional
The measured semiminor axis of the crater in meters.
measured_orientation : float, optional
The measured orientation of the crater in degrees.
measured_diameter : float, optional
The measured diameter of the crater in meters.
measured_radius : float, optional
The measured radius of the crater in meters.
measured_location : pair of floats, optional
The measured (longitude, latitude) location of the crater.
measured_rim_height : float, optional
The measured rim height of the crater in meters.
measured_floor_depth : float, optional
The measured floor depth of the crater in meters.
degradation_state : float, optional
The current degradation state of the crater in meters squared.
production_ND: tuple[float, float, float] | tuple[float, float], optional
The production N(D) value for quasi-Monte Carlo emplacement. Can either be the form of a tuple, of [D, N, Nsigma] or [D, N], where D is the reference diameter in km and N is number of craters per 1 million sq. km, with an optional 1 sigma undertainty Nsigma
production_time: tuple[float, float] | float, optional
The production time for quasi-Monte Carlo emplacement. Can either be in the form of a tuple of [t,tsigma] or a float of t, where t is the emplacement time in My with an optional 1 sigma uncertainty tsigma.
production_sequence: int | None = None,
A positive integer value used to constrain the sequence of quasi-Monte Carlo craters. Craters with a sequence number can only form after craters with a lower sequence number.
simdir : str | Path
|simdir|
rng : numpy.random.Generator | None
|rng|
rng_seed : Any type allowed by the rng_seed argument of numpy.random.Generator, optional
|rng_seed|
rng_state : dict, optional
|rng_state|
check_redundant_inputs : bool, optional
If True, check for redundant inputs such as providing both diameter and radius. Default is True.
**kwargs : Any
|kwargs|
Returns
-------
Crater
A frozen Crater dataclass with derived attributes.
Notes
-----
- Unless `check_redundant_inputs` is False, the following rules apply:
- Exactly one of the following must be provided: `diameter`, `radius`, both of `semimajor_axis` and `semiminor_axis`, `transient_diameter`, `transient_radius`, `projectile_diameter`, `projectile_radius`, or `projectile_mass`.
- Velocity may be specified in one of these ways:
- `projectile_mean_velocity` alone (samples a velocity)
- Any two of (`projectile_velocity`, `projectile_vertical_velocity`, `projectile_angle`). the third is inferred.
- `projectile` is mutually exclusive with velocity-related inputs; if provided, it overrides velocity, angle, direction, and density unless explicitly set.
- The `scaling`, and `rng` models are required for scaling and density inference, but are not stored in the returned Crater object.
- If providing measured properties, either both `measured_semimajor_axis` and `measured_semiminor_axis`, or one of `measured_diameter` or `measured_radius` must be provided.
- if `check_redundant_inputs` is True and any of the rules are violated, then it is assumed that the user wants a copy. In that case, some parameters will not be computed and only the user-input arguments will be used
"""
from cratermaker.components.projectile import Projectile
from cratermaker.components.scaling import Scaling
from cratermaker.components.target import Target
make_copy = crater is not None
def _set_id(**kwargs: Any) -> str:
"""
Sets the hash id of the crater based on input parameters.
This is used as a unique identifier for the crater. To reduce storage requirements, the id is a 32-bit unsigned integer derived from the hash of the input parameters. There is a very small chance of hash collisions, but this is acceptable for the purposes of crater identification.
Parameters
----------
**kwargs : Any
|kwargs|
"""
id_args = [
"time",
"semimajor_axis",
"projectile_angle",
"location",
"orientation",
"transient_diameter",
"projectile_diameter",
"projectile_velocity",
"projectile_mass",
"semiminor_axis",
]
combined_args = [k for k in id_args if k in kwargs and kwargs[k] is not None]
combined = "".join(f"{kwargs[k]}" for k in combined_args)
hexid = hashlib.shake_256(combined.encode()).hexdigest(4)
return np.uint32(int(f"0x{hexid}", 16))
# Convert from the old API "final_diameter/final_radius" to "diameter/radius"
if final_diameter is not None:
diameter = final_diameter
if final_radius is not None:
radius = final_radius
if radius is not None and diameter is not None:
if check_redundant_inputs:
raise ValueError("Only one of diameter or radius may be set.")
else:
diameter = None # Give priority to radius if both are set and we aren't checking for redundant inputs
if latitude is not None or longitude is not None:
if latitude is None or longitude is None:
raise ValueError("Both longitude and latitude must be passed or location must be passed.")
if check_redundant_inputs and (location is not None or projectile_location is not None):
raise ValueError("location cannot be used with longitude and latitude as separate arguments")
location = [longitude, latitude]
if location is None:
location = projectile_location
elif projectile_location is not None:
if check_redundant_inputs:
raise ValueError("Only one of location or projectile_location may be set.")
else:
projectile_location = None # Give priority to location if both are set and we aren't checking for redundant inputs
# Validate ellipticity parameters
if semimajor_axis is not None or semiminor_axis is not None:
if semimajor_axis is None or semiminor_axis is None: # We will assume circular if only one is set
radius = semimajor_axis if semiminor_axis is None else semiminor_axis
semimajor_axis = radius
semiminor_axis = radius
if check_redundant_inputs and (diameter is not None or radius is not None):
raise ValueError(
"diameter or radius cannot be used for elliptical craters; use semimajor_axis and semiminor_axis instead."
)
# Turn any circular crater arguments into elliptial ones
if semimajor_axis is None and semiminor_axis is None:
if radius is not None:
semimajor_axis = radius
semiminor_axis = radius
elif diameter is not None:
semimajor_axis = diameter / 2.0
semiminor_axis = diameter / 2.0
# For now, crater orientation and projectile direction are linked
if check_redundant_inputs and orientation is not None and projectile_direction is not None:
raise ValueError("Only one of orientation or projectile_direction may be set for elliptical craters.")
if orientation is not None and projectile_direction is None:
projectile_direction = orientation
if projectile_direction is not None and orientation is None:
orientation = projectile_direction
velocity_inputs = {
"projectile_velocity": projectile_velocity,
"projectile_vertical_velocity": projectile_vertical_velocity,
"projectile_angle": projectile_angle,
}
n_velocity_inputs = sum(x is not None for x in velocity_inputs.values())
if n_velocity_inputs > 2:
if check_redundant_inputs:
raise ValueError(f"Only two of {', '.join(k for k, v in velocity_inputs.items() if v is not None)} may be set.")
else:
make_copy = True
if projectile_mean_velocity is not None and (projectile_velocity is not None or projectile_vertical_velocity is not None):
if check_redundant_inputs:
raise ValueError("projectile_mean_velocity cannot be used with projectile_velocity or projectile_vertical_velocity")
else:
make_copy = True
size_inputs = {
"semimajor_axis": semimajor_axis,
"transient_diameter": transient_diameter,
"transient_radius": transient_radius,
"projectile_diameter": projectile_diameter,
"projectile_radius": projectile_radius,
"projectile_mass": projectile_mass,
}
# Assemble initial arguments
args = {
"semimajor_axis": semimajor_axis,
"semiminor_axis": semiminor_axis,
"orientation": orientation,
"transient_diameter": transient_diameter,
"projectile_diameter": projectile_diameter,
"projectile_mass": projectile_mass,
"projectile_velocity": projectile_velocity,
"projectile_angle": projectile_angle,
"projectile_density": projectile_density,
"morphology_type": "Not Set",
"location": location,
"time": time,
}
if measured_diameter is not None:
args["measured_diameter"] = measured_diameter
if measured_radius is not None:
args["measured_radius"] = measured_radius
if measured_semimajor_axis is not None:
args["measured_semimajor_axis"] = measured_semimajor_axis
if measured_semiminor_axis is not None:
args["measured_semiminor_axis"] = measured_semiminor_axis
if measured_orientation is not None:
args["measured_orientation"] = measured_orientation
if measured_location is not None:
args["measured_location"] = measured_location
if measured_rim_height is not None:
args["measured_rim_height"] = measured_rim_height
if measured_floor_depth is not None:
args["measured_floor_depth"] = measured_floor_depth
if degradation_state is not None:
args["degradation_state"] = degradation_state
if production_time is not None:
args["production_time"] = production_time
if production_ND is not None:
args["production_ND"] = production_ND
if production_sequence is not None:
args["production_sequence"] = production_sequence
n_size_inputs = sum(v is not None for v in size_inputs.values())
# Process input crater object if provided
if crater is not None:
if not isinstance(crater, Crater):
raise TypeError("crater must be a Crater object.")
old_parameters = {}
for field in crater.as_dict():
if field in locals() and locals()[field] is None and getattr(crater, field) is not None:
old_parameters[field] = getattr(crater, field)
if (
n_size_inputs == 0
): # The user has not passed any size parameters, so we will use the size parameters from the crater object
n_size_inputs = 1 # Make sure we don't trigger the error below
make_copy = True
else: # The user is passing a size parameter, so we cannot use any of the values from the crater object
make_copy = False
for field in size_inputs:
old_parameters.pop(field, None)
if projectile_mean_velocity is None:
# Be sure to keep only two velocitity components, with a preference for angle over vertical velocity
if n_velocity_inputs == 0:
old_parameters.pop("projectile_vertical_velocity", None)
elif n_velocity_inputs == 1:
for k, v in velocity_inputs.items():
if v is not None:
old_parameters.pop(k, None)
if k == "projectile_velocity" or k == "projectile_angle":
old_parameters.pop("projectile_vertical_velocity", None)
elif k == "projectile_vertical_velocity":
old_parameters.pop("projectile_velocity", None)
elif n_velocity_inputs == 2:
for field in velocity_inputs:
old_parameters.pop(field, None)
# Make sure we don't override the projectile_velocity and projectile_angle
n_velocity_inputs = 2
# Now set the local arguments to be what's left from the old_parameters
for field in old_parameters:
if field in args and args[field] is None or field not in args:
args[field] = old_parameters[field]
if n_size_inputs != 1:
if check_redundant_inputs:
raise ValueError(f"Exactly one of {', '.join(k for k, v in size_inputs.items() if v is not None)} must be set.")
else:
make_copy = True
if make_copy:
# The following will change the id
a = args.pop("semimajor_axis")
b = args.pop("semiminor_axis")
td = args.pop("transient_diameter")
pdir = args.pop("orientation")
pd = args.pop("projectile_diameter")
pm = args.pop("projectile_mass")
pv = args.pop("projectile_velocity")
pang = args.pop("projectile_angle")
prho = args.pop("projectile_density")
location = args.pop("location")
if "measured_diameter" in args:
measured_diameter = args.pop("measured_diameter")
measured_semimajor_axis = measured_diameter / 2.0
measured_semiminor_axis = measured_diameter / 2.0
measured_diameter = None
measured_radius = None
elif "measured_radius" in args:
measured_radius = args.pop("measured_radius")
measured_semimajor_axis = measured_radius
measured_semiminor_axis = measured_radius
measured_diameter = None
measured_radius = None
elif "measured_semimajor_axis" in args or "measured_semiminor_axis" in args:
measured_semimajor_axis = args.pop("measured_semimajor_axis", None)
measured_semiminor_axis = args.pop("measured_semiminor_axis", None)
measured_orientation = args.pop("measured_orientation", None)
measured_location = args.pop("measured_location", None)
name = args.pop("name", None)
if crater is not None:
mt = crater.morphology_type
else:
mt = kwargs.pop("morphology_type", "UNKNOWN")
else:
# --- Merge the CommonArgs with any arguments passed here using CratermakerBase ---
argproc = CratermakerBase(simdir=simdir, rng=rng, rng_seed=rng_seed, rng_state=rng_state, **kwargs)
if scaling is not None and isinstance(scaling, Scaling):
if projectile is None:
projectile = scaling.projectile
if target is None:
target = scaling.target
kwargs = {**kwargs, **vars(argproc.common_args)}
target = Target.maker(target, **kwargs)
projectile_args = {
"velocity": args["projectile_velocity"],
"mass": args["projectile_mass"],
"angle": args["projectile_angle"],
"direction": args["orientation"],
"location": args["location"],
"density": args["projectile_density"],
}
projectile = Projectile.maker(
projectile,
target=target,
**{**projectile_args, **kwargs},
).new_projectile(
**projectile_args,
)
scaling = Scaling.maker(
scaling,
target=target,
projectile=projectile,
**kwargs,
)
projectile = scaling.projectile
scaling.recompute()
target = scaling.target
pv = projectile.velocity
pvv = projectile.vertical_velocity
pang = projectile.angle
pdir = projectile.direction
prho = projectile.density
location = projectile.location
# --- Ensure velocity/angle are all set ---
n_set = sum(x is not None for x in [pv, pvv, pang])
if n_set != 3:
raise ValueError("Not enough information to infer a projectile velocity.")
# --- Compute derived quantities ---
fr = radius
fd = diameter
a = args["semimajor_axis"]
b = args["semiminor_axis"]
if a is not None and b is not None:
fr = math.sqrt(a * b)
if fr is not None:
fd = 2 * fr
elif fd is not None:
fr = fd / 2.0
td = args["transient_diameter"]
tr = transient_radius
pd = args["projectile_diameter"]
pr = projectile_radius
pm = projectile_mass
mt = None
if pr is not None:
pd = 2 * pr
if fr is not None:
fd = 2 * fr
if tr is not None:
td = 2 * tr
if fd is not None:
td, mt = scaling.final_to_transient(fd)
pd = scaling.transient_to_projectile(td)
elif td is not None:
fd, mt = scaling.transient_to_final(td)
pd = scaling.transient_to_projectile(td)
elif pd is not None:
td = scaling.projectile_to_transient(pd)
fd, mt = scaling.transient_to_final(td)
elif pm is not None:
pr = ((3.0 * pm) / (4.0 * math.pi * prho)) ** (1.0 / 3.0)
pd = 2.0 * pr
td = scaling.projectile_to_transient(pd)
fd, mt = scaling.transient_to_final(td)
pr = pd / 2
tr = td / 2
fr = fd / 2
pm = (4.0 / 3.0) * math.pi * pr**3 * prho
if a is None:
a = fr
if b is None:
b = fr
if a < 0 or b < 0 or td < 0 or pd < 0:
raise ValueError("Crater and projectile sizes must be non-negative.")
if pv < 0:
raise ValueError("Projectile velocity must be non-negative.")
if prho < 0:
raise ValueError("Projectile density must be non-negative.")
if check_redundant_inputs and measured_radius is not None and measured_diameter is not None:
raise ValueError("Only one of measured_diameter or measured_radius may be set.")
else:
if measured_diameter is not None:
measured_semimajor_axis = measured_diameter / 2.0
measured_semiminor_axis = measured_diameter / 2.0
measured_diameter = None
if measured_radius is not None:
measured_semimajor_axis = measured_radius
measured_semiminor_axis = measured_radius
measured_radius = None
if measured_semimajor_axis is not None or measured_semiminor_axis is not None:
if measured_semimajor_axis is None or measured_semiminor_axis is None:
if check_redundant_inputs and (measured_diameter is not None or measured_radius is not None):
raise ValueError(
"If providing measured properties for an elliptical crater, either both measured_semimajor_axis and measured_semiminor_axis, or one of measured_diameter or measured_radius must be provided."
)
measured_radius = measured_semimajor_axis if measured_semiminor_axis is None else measured_semiminor_axis
measured_semimajor_axis = measured_radius
measured_semiminor_axis = measured_radius
if check_redundant_inputs and (measured_diameter is not None or measured_radius is not None):
raise ValueError(
"measured_diameter or measured_radius cannot be used for elliptical craters; use measured_semimajor_axis and measured_semiminor_axis instead."
)
elif measured_radius is not None:
measured_semimajor_axis = measured_radius
measured_semiminor_axis = measured_radius
else:
measured_semimajor_axis = a if a is not None else None
measured_semiminor_axis = b if b is not None else None
if measured_orientation is None:
measured_orientation = pdir if pdir is not None else None
if measured_location is None:
measured_location = (float(location[0]), float(location[1])) if location is not None else None
# Assemble the final set of arguments to pass to the Crater constructor
args["semimajor_axis"] = float(a) if a is not None else None
args["semiminor_axis"] = float(b) if b is not None else None
args["transient_diameter"] = float(td) if td is not None else None
args["orientation"] = float(pdir) if pdir is not None else None
args["projectile_diameter"] = float(pd) if pd is not None else None
args["projectile_mass"] = float(pm) if pm is not None else None
args["projectile_velocity"] = float(pv) if pv is not None else None
args["projectile_angle"] = float(pang) if pang is not None else None
args["location"] = (float(location[0]), float(location[1]))
args["morphology_type"] = str(mt) if mt is not None else None
if measured_location is not None:
args["measured_location"] = (float(measured_location[0]), float(measured_location[1]))
args["measured_semimajor_axis"] = float(measured_semimajor_axis) if measured_semimajor_axis is not None else None
args["measured_semiminor_axis"] = float(measured_semiminor_axis) if measured_semiminor_axis is not None else None
args["measured_orientation"] = float(measured_orientation) if measured_orientation is not None else None
args["id"] = _set_id(**args)
args["name"] = str(name) if name is not None else None
return cls(**args, **kwargs)
[docs]
def to_geoseries(
self,
n: int = 121,
surface: Surface | None = None,
split_antimeridian: bool = True,
use_measured_properties: bool = True,
**kwargs: Any,
) -> GeoSeries:
"""
Geodesic ellipse on a sphere: for each bearing theta from the center, we shoot a geodesic with distance r(theta) = (a*b)/sqrt((b*ct)^2 + (a*st)^2), then rotate all bearings by the crater's orientation.
Parameters
----------
n : int, optional
Number of points to use for the polygon, by default 121 (every 3 degrees plus the end point).
surface : Surface | None, optional
Surface object providing planetary radius and CRS.
split_antimeridian : bool, optional
If True, split the polygon into a GeometryCollection if it crosses the antimeridian, by default True.
use_measured_properties : bool, optional
If True, use the current measured crater properties (semimajor_axis, semiminor_axis, location, orientation) instead of the initial ones, by default True.
**kwargs : Any
|kwargs|
Returns
-------
A GeoSeries containing a Shapely Polygon in lon/lat degrees.
"""
from pyproj import Geod
from shapely.geometry import GeometryCollection, LineString, Polygon
from shapely.ops import split
from cratermaker.components.surface import Surface
surface = Surface.maker(surface)
geod = Geod(a=surface.target.radius, b=surface.target.radius)
theta = np.linspace(0.0, 360.0, num=n, endpoint=True)
if use_measured_properties:
a = self.measured_semimajor_axis
b = self.measured_semiminor_axis
lon, lat = self.measured_location
phi = self.measured_orientation
radius = self.measured_radius
else:
a = self.semimajor_axis
b = self.semiminor_axis
lon, lat = self.location
phi = self.orientation
radius = self.radius
# Measure the rim height so that the polygon sits on to of the surface rather than underneath
region = surface.extract_region(location=(lon, lat), region_radius=radius, at_least_one_face=True)
rim_height = np.max(region.face_elevation)
z = np.full_like(theta, rim_height)
# Polar radius of an axis-aligned ellipse in a Euclidean tangent plane
ct = np.cos(np.deg2rad(theta))
st = np.sin(np.deg2rad(theta))
r = (a * b) / np.sqrt((b * ct) ** 2 + (a * st) ** 2)
# Bearings (from east, CCW) rotated by azimuth
bearings = theta + phi % 360.0
# Forward geodesic for each bearing/distance
poly_lon, poly_lat, _ = geod.fwd(lon * np.ones_like(bearings), lat * np.ones_like(bearings), bearings, r)
# Correct for potential antimeridian crossing
if split_antimeridian and np.ptp(poly_lon) > 180.0:
center_sign = np.sign(lon)
poly_lon = np.where(np.sign(poly_lon) != center_sign, poly_lon + 360.0 * center_sign, poly_lon)
poly = Polygon(zip(poly_lon, poly_lat, z, strict=True))
merdian_lon = 180.0 * center_sign
meridian = LineString([(merdian_lon, -90.0), (merdian_lon, 90.0)])
poly = split(poly, meridian)
def lon_flip(lon, lat, z=None):
lon = np.where(np.abs(lon) >= 180.0, lon - 360.0 * np.sign(lon), lon)
return lon, lat
new_geoms = []
for p in poly.geoms:
if np.abs(p.centroid.x) > 180.0:
p = transform(lon_flip, p)
new_geoms.append(p)
poly = GeometryCollection(new_geoms)
else:
poly = Polygon(zip(poly_lon, poly_lat, z, strict=True))
return GeoSeries([poly], crs=surface.crs)
[docs]
@classmethod
def from_file(cls, filename: str | Path, **kwargs: Any) -> list[Crater] | None:
"""
Load a list of craters from a file.
Parameters
----------
filename : str | Path
The path to the file to load from.
**kwargs : Any
|kwargs|
Returns
-------
list[Crater] | None
A list of Crater objects loaded from the file, or None if no data.
"""
filename = Path(filename)
if not filename.exists():
raise FileNotFoundError(f"File {filename} does not exist.")
extension = filename.suffix.lower().lstrip(".")
if extension == "nc":
ds = xr.open_dataset(filename)
craters = cls.from_xarray(ds, **kwargs)
elif extension == "csv":
craters = cls.from_csv_file(filename, **kwargs)
elif extension == "scc":
craters = cls.from_scc_file(filename, **kwargs)
return craters
[docs]
@classmethod
def from_csv_file(cls, input_file: Path | str) -> list[Crater]:
"""
Import crater data from a CSV file.
Parameters
----------
input_file : Path | str
The path to the CSV file containing crater data.
Returns
-------
list[Crater]
A list of Crater objects imported from the CSV file.
"""
import csv
craters = []
input_file = Path(input_file)
if not input_file.exists():
raise FileNotFoundError(f"Input file '{input_file}' does not exist.")
with input_file.open(mode="r", newline="") as csvfile:
reader = csv.DictReader(csvfile)
for row in reader:
crater_data = _convert_tuple_vars(input_dict=row, inverse=True)
crater_data = {k: v for k, v in crater_data.items() if v is not None and v != ""}
for key, value in crater_data.items():
if key in ["id", "production_sequence"]:
crater_data[key] = int(value)
elif isinstance(value, tuple | list):
vnew = []
for v in value:
if v is None or v == "":
continue
vnew.append(float(v))
if len(vnew) == len(value):
crater_data[key] = vnew
else:
crater_data[key] = None
else:
try:
crater_data[key] = float(value)
except ValueError:
crater_data[key] = value
crater_data = {k: v for k, v in crater_data.items() if v is not None}
crater = cls.maker(**crater_data, check_redundant_inputs=False)
craters.append(crater)
return craters
[docs]
@classmethod
def from_scc_file(cls, input_file: Path | str) -> list[Crater]:
"""
Import crater data from a Spatial Crater Count file.
Parameters
----------
input_file : Path | str
The path to the SCC file containing crater data.
Returns
-------
list[Crater]
A list of Crater objects imported from the SCC file.
"""
from craterstats import Spatialcount
craters = []
input_file = Path(input_file)
if not input_file.exists():
raise FileNotFoundError(f"Input file '{input_file}' does not exist.")
if input_file.suffix != ".scc":
raise ValueError(f"Input file '{input_file}' is not a .scc file.")
scc = Spatialcount(filename=str(input_file))
for diam, lon, lat in zip(scc.diam, scc.lon, scc.lat, strict=True):
crater = cls.maker(diameter=diam * 1e3, location=(lon, lat))
craters.append(crater)
return craters
[docs]
@classmethod
def from_xarray(cls, dataset: xr.Dataset | dict, interval: int | None = None) -> list[Crater]:
"""
Import crater data from an xarray Dataset.
Parameters
----------
dataset : xr.Dataset | dict
The xarray Dataset containing crater data or a dictionary of xarray Datasets keyed by interval number.
Returns
-------
list[Crater]
A list of Crater objects imported from the xarray Dataset.
"""
craters = []
if type(dataset) is dict:
if interval is None:
dataset = dataset[-1]
elif interval in dataset:
dataset = dataset[interval]
else:
return craters
if "interval" in dataset.coords:
if interval is None:
dataset = dataset.isel(interval=-1)
elif interval in dataset.interval:
dataset = dataset.sel(interval=interval)
else:
return craters
dataset.load()
if len(dataset) == 0:
return craters
for id in tqdm(dataset.id.data, desc="Converting xarray Dataset to Crater objects", unit="crater", position=0, leave=False):
crater_data = dataset.sel(id=id).to_dict()["data_vars"]
crater_data = {k: v["data"] for k, v in crater_data.items()}
if not isinstance(crater_data["semimajor_axis"], (float, int)):
continue
crater_data = _convert_tuple_vars(input_dict=crater_data, inverse=True)
for k, v in crater_data.items():
if v is not None and np.any(np.isreal(v)) and np.any(np.isnan(v)):
crater_data[k] = None
crater = cls.maker(**crater_data, check_redundant_inputs=False)
craters.append(crater)
return craters
@property
def final_diameter(self) -> float | None:
"""Final diameter of the crater in meters."""
return self.diameter
@property
def final_radius(self) -> float | None:
"""Final radius of the crater in meters."""
return self.radius
[docs]
def remove_complex_data(self, **kwargs) -> None:
"""
Remove complex data types from the crater variable properties.
This is useful for storing a lightweight representation, as it removes complex data types that can be recomputed from the fixed properties and the morphology model when needed. The base Crater type does not have any complex data types, but is used by derived types that do, so this method is included here for consistency and to allow for future expansion of the base Crater type without breaking the API.
"""
pass
def _convert_tuple_vars(input_dict: dict, inverse: bool = False) -> dict:
tuple_map = {
"location": ["longitude", "latitude"],
"measured_location": ["measured_longitude", "measured_latitude"],
"production_ND": ["production_D", "production_N", "production_N_stdev"],
"production_time": ["production_time", "production_time_stdev"],
}
new_dict = {k: v for k, v in input_dict.items() if v is not None and v != ""}
if inverse:
# Process aliases
if "production_D" in input_dict and "production_N" in input_dict:
prod_diam = input_dict.pop("production_D")
nval = input_dict.pop("production_N")
nstdev = input_dict.pop("production_N_stdev", 0.0)
if nstdev is None or nstdev == "":
nstdev = 0.0
new_dict["production_ND"] = [prod_diam, nval, nstdev]
if "production_time" in input_dict:
time_mean = input_dict.pop("production_time")
time_stdev = input_dict.pop("production_time_stdev", 0.0)
if time_stdev is None or time_stdev == "":
time_stdev = 0.0
new_dict["production_time"] = [time_mean, time_stdev]
for tup, varlist in tuple_map.items():
if inverse:
if any(var in varlist for var in input_dict):
if tup not in input_dict:
new_dict[tup] = [None] * len(varlist)
for idx, var in enumerate(varlist):
if var != tup and var in input_dict and input_dict[var] is not None:
new_dict[tup][idx] = input_dict.get(var)
else:
if tup in input_dict:
_ = new_dict.pop(tup, None)
for idx, var in enumerate(varlist):
new_dict[var] = input_dict[tup][idx]
return new_dict
import_components(__name__, __path__)
