Skip to content

skais_mapper.simobjects

Tools for manipulating data objects from simulations.

Classes:

Name Description
GasolineGalaxy

A Gasoline Galaxy parser.
SPHGalaxy

A generic base SPH Galaxy simulation parser.
TNGGalaxy

IllustrisTNG Galaxy simulation parser.

ArepoGalaxy 

ArepoGalaxy(
    base_path: str | Path,
    snapshot: int,
    halo_index: int,
    **kwargs,
)

Bases: TNGGalaxy

An Arepo Galaxy parser (alias to TNGGalaxy for now).

Source code in skais_mapper/simobjects.py 
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
def __init__(self, base_path: str | Path, snapshot: int, halo_index: int, **kwargs):
    """Initialize a TNGGalaxy instance.

    Args:
        base_path (str): Base path to the Illustris(TNG) snapshots.
        snapshot (int): Snapshot ID {0-99}.
        halo_index (int): Halo index, index of the subhalo ID list.
        kwargs (dict): Additional keyword arguments
            - verbose (bool): If True, print information to the command line.
    """
    kwargs.setdefault("verbose", False)
    self.base_path = base_path
    self.snapshot = snapshot
    self._halo_index = halo_index
    self.subhalo_ids = self.subhalo_list(base_path, snapshot, verbose=kwargs["verbose"])
    self.subhalo = self.load_subhalo(self.halo_index, verbose=kwargs["verbose"])
    super().__init__(**kwargs)

GasolineGalaxy 

GasolineGalaxy(
    ra: float = 0.0,
    dec: float = 0.0,
    distance: float = 3.0,
    peculiar_v: float = 0.0,
    rotation: Callable | None = None,
    cosmo_pars: dict | None = None,
    units: dict | None = None,
    particle_type: str | None = None,
    as_float32: bool = False,
    verbose: bool = False,
)

Bases: SPHGalaxy

A Gasoline Galaxy parser.

Source code in skais_mapper/simobjects.py 
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
def __init__(
    self,
    ra: float = 0.0,
    dec: float = 0.0,
    distance: float = 3.0,
    peculiar_v: float = 0.0,
    rotation: Callable | None = None,
    cosmo_pars: dict | None = None,
    units: dict | None = None,
    particle_type: str | None = None,
    as_float32: bool = False,
    verbose: bool = False,
):
    """Initialize an SPHGalaxy base class instance.

    Args:
        ra: The right ascension sky coordinate
        dec: The declination sky coordinate
        distance: The distance from observer
        peculiar_v: The galaxy's peculiar velocity
        rotation: Arbitrary rotation of the galaxy
        cosmo_pars: Cosmology settings
        units: unit system in the simulations
        particle_type: the particle type name e.g. 'gas', 'dm', or 'star' (see illustris.util)
        as_float32: If True, load simulation data as 32-bit floats (for efficiency)
        verbose: If True, print information to the command line
    """
    self.as_float32 = as_float32
    self.header = self.load_header()
    cosmo_pars = self.load_cosmology(cosmo_pars)
    self.cosmology = CosmoModel(**cosmo_pars)
    if particle_type is None:
        particle_type = "gas"
    self._p_idx = tng.util.pidx_from_ptype(particle_type)
    self.data = self.load_data(as_float32=self.as_float32, verbose=verbose)
    self.ra = ra * au.deg
    self.dec = dec * au.deg
    self.distance = distance * au.Mpc
    self.peculiar_v = peculiar_v * au.km / au.s
    self.rotation = rotation
    if units:
        self.set_units(**units)
    self.verbose = verbose

SPHGalaxy 

SPHGalaxy(
    ra: float = 0.0,
    dec: float = 0.0,
    distance: float = 3.0,
    peculiar_v: float = 0.0,
    rotation: Callable | None = None,
    cosmo_pars: dict | None = None,
    units: dict | None = None,
    particle_type: str | None = None,
    as_float32: bool = False,
    verbose: bool = False,
)

A generic base SPH Galaxy simulation parser.

Initialize an SPHGalaxy base class instance.

Parameters:

Name Type Description Default
ra float

The right ascension sky coordinate

 0.0
dec float

The declination sky coordinate

 0.0
distance float

The distance from observer

 3.0
peculiar_v float

The galaxy's peculiar velocity

 0.0
rotation Callable | None

Arbitrary rotation of the galaxy

 None
cosmo_pars dict | None

Cosmology settings

 None
units dict | None

unit system in the simulations

 None
particle_type str | None

the particle type name e.g. 'gas', 'dm', or 'star' (see illustris.util)

 None
as_float32 bool

If True, load simulation data as 32-bit floats (for efficiency)

 False
verbose bool

If True, print information to the command line

 False

Methods:

Name Description
angular_distance

Calculate the angular distance of the simulation at the given redshift.
angular_resolution

Calculate the angular scale of the simulation at the given redshift.
kd_tree

Compute the distance of each gas particle to its k nearest neighbors.
load_cosmology

Dummy method, to be overridden in subclasses.
load_data

Dummy method, to be overridden in subclasses.
load_header

Dummy method, to be overridden in subclasses.
set_units

Set units in the header attribute.
units

Get common units from descriptor strings ('l', 'm', 'v', and variants).

Attributes:

Name Type Description
UnitLength au.Quantity

The simulation's length unit as astropy.units.Quantity.
UnitMass au.Quantity

The simulation's mass unit as astropy.units.Quantity.
UnitVelocity au.Quantity

The simulation's velocity unit as astropy.units.Quantity.
boxsize au.Quantity

The virtual boxsize of the simulation.
cell_positions au.Quantity | None

The simulation's cell position data in corresponding units.
density au.Quantity | None

The simulation's density data in corresponding units.
internal_energy au.Quantity | None

The simulation's internal energy data in corresponding units.
m_H au.Quantity | None

The simulation's hydrogen mass data in corresponding units.
m_HI au.Quantity | None

The simulation's ionized hydrogen mass data in corresponding units.
masses au.Quantity | None

The simulation's mass data in corresponding units.
n_H au.Quantity | None

The simulation's hydrogen number density data in corresponding units.
n_HI au.Quantity | None

The simulation's ionized hydrogen number density data in corresponding units.
p_idx int

Particle index {0: 'gas', 1: 'dm', 2: 'tracers', 3: 'stars', 4: 'BHs'}.
particle_mass au.Quantity | float | None

The simulation's particle mass (for constant mass particles).
particle_positions au.Quantity | None

The simulation's particle position data in corresponding units.
particle_type str

Particle type {0: 'gas', 1: 'dm', 2: 'tracers', 3: 'stars', 4: 'BHs'}.
velocities au.Quantity | None

The simulation's particle velocity data in corresponding units.
x_H au.Quantity | None

The simulation's neutral hydrogen abundance data in corresponding units.
x_HI au.Quantity | None

The simulation's ionized hydrogen abundance data in corresponding units.
x_e au.Quantity | None

The simulation's electron abundance data in corresponding units.
Source code in skais_mapper/simobjects.py 
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
def __init__(
    self,
    ra: float = 0.0,
    dec: float = 0.0,
    distance: float = 3.0,
    peculiar_v: float = 0.0,
    rotation: Callable | None = None,
    cosmo_pars: dict | None = None,
    units: dict | None = None,
    particle_type: str | None = None,
    as_float32: bool = False,
    verbose: bool = False,
):
    """Initialize an SPHGalaxy base class instance.

    Args:
        ra: The right ascension sky coordinate
        dec: The declination sky coordinate
        distance: The distance from observer
        peculiar_v: The galaxy's peculiar velocity
        rotation: Arbitrary rotation of the galaxy
        cosmo_pars: Cosmology settings
        units: unit system in the simulations
        particle_type: the particle type name e.g. 'gas', 'dm', or 'star' (see illustris.util)
        as_float32: If True, load simulation data as 32-bit floats (for efficiency)
        verbose: If True, print information to the command line
    """
    self.as_float32 = as_float32
    self.header = self.load_header()
    cosmo_pars = self.load_cosmology(cosmo_pars)
    self.cosmology = CosmoModel(**cosmo_pars)
    if particle_type is None:
        particle_type = "gas"
    self._p_idx = tng.util.pidx_from_ptype(particle_type)
    self.data = self.load_data(as_float32=self.as_float32, verbose=verbose)
    self.ra = ra * au.deg
    self.dec = dec * au.deg
    self.distance = distance * au.Mpc
    self.peculiar_v = peculiar_v * au.km / au.s
    self.rotation = rotation
    if units:
        self.set_units(**units)
    self.verbose = verbose

UnitLength property 

UnitLength: au.Quantity

The simulation's length unit as astropy.units.Quantity.

UnitMass property 

UnitMass: au.Quantity

The simulation's mass unit as astropy.units.Quantity.

UnitVelocity property 

UnitVelocity: au.Quantity

The simulation's velocity unit as astropy.units.Quantity.

boxsize property 

boxsize: au.Quantity

The virtual boxsize of the simulation.

cell_positions property 

cell_positions: au.Quantity | None

The simulation's cell position data in corresponding units.

density property 

density: au.Quantity | None

The simulation's density data in corresponding units.

internal_energy property 

internal_energy: au.Quantity | None

The simulation's internal energy data in corresponding units.

m_H property 

m_H: au.Quantity | None

The simulation's hydrogen mass data in corresponding units.

m_HI property 

m_HI: au.Quantity | None

The simulation's ionized hydrogen mass data in corresponding units.

masses property 

masses: au.Quantity | None

The simulation's mass data in corresponding units.

n_H property 

n_H: au.Quantity | None

The simulation's hydrogen number density data in corresponding units.

n_HI property 

n_HI: au.Quantity | None

The simulation's ionized hydrogen number density data in corresponding units.

p_idx property writable 

p_idx: int

Particle index {0: 'gas', 1: 'dm', 2: 'tracers', 3: 'stars', 4: 'BHs'}.

particle_mass property 

particle_mass: au.Quantity | float | None

The simulation's particle mass (for constant mass particles).

particle_positions property 

particle_positions: au.Quantity | None

The simulation's particle position data in corresponding units.

particle_type property writable 

particle_type: str

Particle type {0: 'gas', 1: 'dm', 2: 'tracers', 3: 'stars', 4: 'BHs'}.

velocities property 

velocities: au.Quantity | None

The simulation's particle velocity data in corresponding units.

x_H property 

x_H: au.Quantity | None

The simulation's neutral hydrogen abundance data in corresponding units.

x_HI property 

x_HI: au.Quantity | None

The simulation's ionized hydrogen abundance data in corresponding units.

x_e property 

x_e: au.Quantity | None

The simulation's electron abundance data in corresponding units.

angular_distance 

angular_distance(eps: float = 0.0001) -> float

Calculate the angular distance of the simulation at the given redshift.

Parameters:

Name Type Description Default
eps float

minimum cutoff value for the redshift

 0.0001
Source code in skais_mapper/simobjects.py 
134
135
136
137
138
139
140
141
142
143
144
def angular_distance(self, eps: float = 1e-4) -> float:
    """Calculate the angular distance of the simulation at the given redshift.

    Args:
        eps: minimum cutoff value for the redshift
    """
    z = self.cosmology.z
    z += eps if z < eps else 0
    dz = self.cosmology.d_z(z, cosmo_model=self.cosmology, scaled=False)
    dang = self.cosmology.d_z2kpc(dz, cosmo_model=self.cosmology)
    return dang

angular_resolution 

angular_resolution(
    eps: float = 0.0001,
) -> tuple[float, float]

Calculate the angular scale of the simulation at the given redshift.

Parameters:

Name Type Description Default
eps float

minimum cutoff value for the redshift

 0.0001
Source code in skais_mapper/simobjects.py 
146
147
148
149
150
151
152
153
154
155
156
def angular_resolution(self, eps: float = 1e-4) -> tuple[float, float]:
    """Calculate the angular scale of the simulation at the given redshift.

    Args:
        eps: minimum cutoff value for the redshift
    """
    z = self.cosmology.z
    z += eps if z < eps else 0
    dz = self.cosmology.d_z(z, cosmo_model=self.cosmology, scaled=True)
    arcsec2kpc = self.cosmology.arcsec2kpc(z, dz)
    return 1.0 / arcsec2kpc.to(au.kpc / au.deg), z

kd_tree 

kd_tree(
    k: int = 8,
    threads: int = 1,
    epsilon: float = 0.0001,
    as_float32: bool = False,
    verbose: bool = True,
) -> au.Quantity | NDArray

Compute the distance of each gas particle to its k nearest neighbors.

Parameters:

Name Type Description Default
k int

Number of nearest neighbors.

 8
threads int

Number of computing threads.

 1
epsilon float

Small increase of the box size to avoid segfaults

 0.0001
as_float32 bool

If True, load simulation data as 32-bit floats (for efficiency)

 False
verbose bool

If True, print out results to the command line.

 True

Returns:

Type Description
np.ndarray

The distance of each particle to its farthest neighbor
Source code in skais_mapper/simobjects.py 
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
def kd_tree(
    self,
    k: int = 8,
    threads: int = 1,
    epsilon: float = 1.0e-4,
    as_float32: bool = False,
    verbose: bool = True,
) -> au.Quantity | NDArray:
    """Compute the distance of each gas particle to its k nearest neighbors.

    Args:
        k: Number of nearest neighbors.
        threads: Number of computing threads.
        epsilon: Small increase of the box size to avoid segfaults
        as_float32: If True, load simulation data as 32-bit floats (for efficiency)
        verbose: If True, print out results to the command line.

    Returns:
        (np.ndarray): The distance of each particle to its farthest neighbor
    """
    p = self.particle_positions
    if p is None:
        return None
    p = p.to(au.kpc)
    boxsize = (1 + epsilon) * self.boxsize.to(au.kpc)
    if verbose:
        print("Building KDTree...")
    t_ini = time.time()
    kdtree = sp.spatial.cKDTree(p, leafsize=16, boxsize=boxsize)
    if verbose:
        print(f"Time to build KDTree: {time.time() - t_ini:.3f} secs.")
    if verbose:
        print("Querying KDTree...")
    dist, _ = kdtree.query(p, k, workers=threads)
    if verbose:
        print(f"Time to querying KDTree: {time.time() - t_ini:.3f} secs.")
    if as_float32:
        return (dist[:, -1] * au.kpc).astype(np.float32)
    return dist[:, -1] * au.kpc

load_cosmology 

load_cosmology(
    cosmo_pars: dict | None, in_place: bool = True
) -> dict

Dummy method, to be overridden in subclasses.

Returns:

Type Description
dict

a dictionary with cosmological parameters pulled from a header or set manually. The dictionary will be used for the CosmoModel dataclass.
Source code in skais_mapper/simobjects.py 
120
121
122
123
124
125
126
127
128
129
130
131
132
def load_cosmology(self, cosmo_pars: dict | None, in_place: bool = True) -> dict:
    """Dummy method, to be overridden in subclasses.

    Returns:
        (dict): a dictionary with cosmological parameters pulled from a header
            or set manually. The dictionary will be used for the CosmoModel
            dataclass.
    """
    if cosmo_pars is None:
        cosmo_pars = {}
    if in_place:
        self.cosmology = CosmoModel(**cosmo_pars)
    return cosmo_pars

load_data 

load_data(**kwargs) -> dict

Dummy method, to be overridden in subclasses.

Returns:

Type Description
dict

the data dictionary loaded from the simulations containing: - 'Density', 'Masses', 'Coordinates', 'Velocities', 'InternalEnergy', 'ElectronAbundance' - optionally 'CenterOfMass', 'GFM_Metals' (axis 1 should contain hydrogen gas fractions), 'NeutralHydrogenAbundance'
Source code in skais_mapper/simobjects.py 
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
def load_data(self, **kwargs) -> dict:
    """Dummy method, to be overridden in subclasses.

    Returns:
        (dict): the data dictionary loaded from the simulations containing:
            - 'Density', 'Masses', 'Coordinates', 'Velocities', 'InternalEnergy',
              'ElectronAbundance'
            - optionally 'CenterOfMass', 'GFM_Metals' (axis 1 should contain
              hydrogen gas fractions), 'NeutralHydrogenAbundance'
    """
    kwargs.setdefault(
        "fields",
        self.primary_hdf5_fields[self.p_idx] + self.optional_hdf5_fields[self.p_idx],
    )
    data = {k: None for k in kwargs["fields"]}
    return data

load_header 

load_header() -> dict

Dummy method, to be overridden in subclasses.

Returns:

Type Description
dict

the header dictionary which should contain - a 'snapshot' and 'catalog' subdictionary - the 'snapshot' subdictionary should contain - 'UnitLength_in_cm', 'UnitVelocity_in_cm_per_s', 'UnitMass_in_g' - 'OmegaLambda', 'Omega0', 'Redshift' / 'Time', 'HubbleParam' - optionally: 'OmegaK'
Source code in skais_mapper/simobjects.py 
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
def load_header(self) -> dict:
    """Dummy method, to be overridden in subclasses.

    Returns:
        (dict): the header dictionary which should contain
            - a 'snapshot' and 'catalog' subdictionary
            - the 'snapshot' subdictionary should contain
                - 'UnitLength_in_cm', 'UnitVelocity_in_cm_per_s', 'UnitMass_in_g'
                - 'OmegaLambda', 'Omega0', 'Redshift' / 'Time', 'HubbleParam'
                - optionally: 'OmegaK'
    """
    self.header = {}
    self.header["snapshot"] = {}
    self.header["snapshot"]["UnitLength_in_cm"] = 3.085678e21
    self.header["snapshot"]["UnitVelocity_in_cm_per_s"] = 1e5
    self.header["snapshot"]["UnitMass_in_g"] = 1.989e43
    return self.header

set_units 

set_units(length: float, velocity: float, mass: float)

Set units in the header attribute.

Parameters:

Name Type Description Default
length float

length unit

 required
velocity float

velocity unit (sets the time unit implicitly)

 required
mass float

mass unit

 required
Source code in skais_mapper/simobjects.py 
211
212
213
214
215
216
217
218
219
220
221
222
223
def set_units(self, length: float, velocity: float, mass: float):
    """Set units in the header attribute.

    Args:
        length: length unit
        velocity: velocity unit (sets the time unit implicitly)
        mass: mass unit
    """
    if not hasattr(self, "header"):
        self.header = self.load_header()
    self.header["snapshot"]["UnitLength_in_cm"] = length
    self.header["snapshot"]["UnitVelocity_in_cm_per_s"] = velocity
    self.header["snapshot"]["UnitMass_in_g"] = mass

units 

units(u: str) -> au.Quantity

Get common units from descriptor strings ('l', 'm', 'v', and variants).

Parameters:

Name Type Description Default
u str

the unit string describing the dimensionality; l -> length; m -> mass; v -> velocity; t -> time; vp -> peculiar velocity; sqrtP -> square root of pressure; c[] -> comoving quantities; []/h -> quantities divided by the dimensionless Hubble const.

 required

Returns:

Type Description
astropy.units.Quantity

specified (scaled/composite) unit.
Source code in skais_mapper/simobjects.py 
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
def units(self, u: str) -> au.Quantity:
    """Get common units from descriptor strings ('l', 'm', 'v', and variants).

    Args:
        u: the unit string describing the dimensionality;
            l -> length;
            m -> mass;
            v -> velocity;
            t -> time;
            vp -> peculiar velocity;
            sqrtP -> square root of pressure;
            c[] -> comoving quantities;
            []/h -> quantities divided by the dimensionless Hubble const.

    Returns:
        (astropy.units.Quantity): specified (scaled/composite) unit.
    """
    h = self.cosmology.h
    a = self.cosmology.a
    match u:
        case "l":
            return self.UnitLength
        case "m":
            return self.UnitMass
        case "v":
            return self.UnitVelocity
        case "t":
            return self.UnitLength / self.UnitVelocity
        case "l/h":
            return self.UnitLength * (a / h)
        case "m/h":
            return self.UnitMass / h
        case "vp":
            return self.UnitVelocity * np.sqrt(a)
        case "sqrtP":
            return (
                h
                / a**2
                * (self.UnitMass / self.UnitLength) ** (1.0 / 2)
                / (self.UnitLength / self.UnitVelocity)
            )
        case _:
            return 1

TNGGalaxy 

TNGGalaxy(
    base_path: str | Path,
    snapshot: int,
    halo_index: int,
    **kwargs,
)

Bases: SPHGalaxy

IllustrisTNG Galaxy simulation parser.

Initialize a TNGGalaxy instance.

Parameters:

Name Type Description Default
base_path str

Base path to the Illustris(TNG) snapshots.

 required
snapshot int

Snapshot ID {0-99}.

 required
halo_index int

Halo index, index of the subhalo ID list.

 required
kwargs dict

Additional keyword arguments - verbose (bool): If True, print information to the command line.

 {}

Methods:

Name Description
generate_map

Generate raytracing projection map.
get_mapping_arrays

Fetch data (arrays or scalars) from a TNGGalaxy object given keys.
load_cosmology

Load the cosmological parameters from the header dictionary.
load_data

Load the snapshot's relevant subset data of the specified group ID.
load_header

Load the header (overrides parent class method).
load_subhalo

Load the FOF subhalo metadata of a specified group ID.
p_idx

Setter for particle index.
subhalo_list

Pull a list of subhalo IDs from the snapshot's FOF group catalog.

Attributes:

Name Type Description
N_particles int

Total number of particles in the galaxy.
N_particles_type list[int]

Number of particles per type in the galaxy.
center au.Quantity

Center position getter of the galaxy in units of L/h.
halo_index

The halo index getter.
hsml

Average kernel smoothing length.
magnetic_field

Magnetic vector field in cgs units of Gauss.
magnetic_field_strength

Magnetic field strength in units of SI Gauss.
r_cell

Alias for radii.
radii

Particle radius distances assuming spherical shape.
temperature

Temperature data getter in units of Kelvin.
Source code in skais_mapper/simobjects.py 
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
def __init__(self, base_path: str | Path, snapshot: int, halo_index: int, **kwargs):
    """Initialize a TNGGalaxy instance.

    Args:
        base_path (str): Base path to the Illustris(TNG) snapshots.
        snapshot (int): Snapshot ID {0-99}.
        halo_index (int): Halo index, index of the subhalo ID list.
        kwargs (dict): Additional keyword arguments
            - verbose (bool): If True, print information to the command line.
    """
    kwargs.setdefault("verbose", False)
    self.base_path = base_path
    self.snapshot = snapshot
    self._halo_index = halo_index
    self.subhalo_ids = self.subhalo_list(base_path, snapshot, verbose=kwargs["verbose"])
    self.subhalo = self.load_subhalo(self.halo_index, verbose=kwargs["verbose"])
    super().__init__(**kwargs)

N_particles property 

N_particles: int

Total number of particles in the galaxy.

N_particles_type property 

N_particles_type: list[int]

Number of particles per type in the galaxy.

center property 

center: au.Quantity

Center position getter of the galaxy in units of L/h.

halo_index property writable 

halo_index

The halo index getter.

hsml property 

hsml

Average kernel smoothing length.

magnetic_field property 

magnetic_field

Magnetic vector field in cgs units of Gauss.

magnetic_field_strength property 

magnetic_field_strength

Magnetic field strength in units of SI Gauss.

r_cell property 

r_cell

Alias for radii.

radii property 

radii

Particle radius distances assuming spherical shape.

temperature property 

temperature

Temperature data getter in units of Kelvin.

generate_map 

generate_map(
    keys: list[str] | None = None,
    factors: list[float] | None = None,
    use_half_mass_rad: bool = True,
    fh: float = 3,
    grid_size: int = 512,
    xaxis: int = 0,
    yaxis: int = 1,
    periodic: bool = True,
    assignment_func: Callable = voronoi_RT_2D,
    tracers: int | None = None,
    divisions: int | None = None,
    rot: list[int] | list[float] | None = None,
    verbose: bool = False,
) -> tuple[au.Quantity, au.Quantity, int] | Any

Generate raytracing projection map.

Parameters:

Name Type Description Default
keys list[str] | None

Keys to fetch data for projecting onto the map.

 None
factors list[float] | None

Factors for modifying the projection data.

 None
use_half_mass_rad bool

If True, the SubhaloHalfmassRad from the subfind catalog is used for selecting relevant particles. Otherwise, a fraction of the entire particle extent is used.

 True
fh float

Expansion factor for the SPH particle radii.

 3
grid_size int

The size of the maps/images. Default: 512.

 512
xaxis int

Projection axis for x.

 0
yaxis int

Projection axis for y.

 1
periodic bool

Use periodic boundary conditions for the projection.

 True
assignment_func Callable

Mass assignment algorithm; one of [voronoi_RT_2D, voronoi_NGP_2D].

 voronoi_RT_2D
tracers int | None

Number of tracer particles to use for the Nearest Grid Point algorithm.

 None
divisions int | None

Number of sphere divisions to use for the Nearest Grid Point algorithm.

 None
rot list[int] | list[float] | None

Angles to rotate the particle positions.

 None
verbose bool

If True, print status updates to command line.

 False

Returns:

Type Description
(np.ndarray, np.ndarray, int)

The projected map, the map extent, and number of particles projected.
Source code in skais_mapper/simobjects.py 
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
def generate_map(
    self,
    keys: list[str] | None = None,
    factors: list[float] | None = None,
    use_half_mass_rad: bool = True,
    fh: float = 3,
    grid_size: int = 512,
    xaxis: int = 0,
    yaxis: int = 1,
    periodic: bool = True,
    assignment_func: Callable = voronoi_RT_2D,
    tracers: int | None = None,
    divisions: int | None = None,
    rot: list[int] | list[float] | None = None,
    verbose: bool = False,
) -> tuple[au.Quantity, au.Quantity, int] | Any:
    """Generate raytracing projection map.

    Args:
        keys: Keys to fetch data for projecting onto the map.
        factors: Factors for modifying the projection data.
        use_half_mass_rad:
          If True, the SubhaloHalfmassRad from the subfind catalog is used for
          selecting relevant particles. Otherwise, a fraction of the entire
          particle extent is used.
        fh: Expansion factor for the SPH particle radii.
        grid_size: The size of the maps/images. Default: 512.
        xaxis: Projection axis for x.
        yaxis: Projection axis for y.
        periodic: Use periodic boundary conditions for the projection.
        assignment_func: Mass assignment algorithm; one of
          [voronoi_RT_2D, voronoi_NGP_2D].
        tracers: Number of tracer particles to use for the Nearest Grid Point algorithm.
        divisions: Number of sphere divisions to use for the Nearest Grid Point algorithm.
        rot:
          Angles to rotate the particle positions.
        verbose: If True, print status updates to command line.

    Returns:
        (np.ndarray, np.ndarray, int): The projected map, the map extent, and
           number of particles projected.
    """
    L, M, T = au.Mpc, au.Msun, au.K
    projected = np.zeros((grid_size, grid_size), dtype=np.float64)
    if keys is None:
        keys = ["particle_positions", "masses", "radii", "center"]
    if factors is None:
        factors = [1, 1, fh, 1]
    if use_half_mass_rad:
        if "SubhaloHalfmassRadType" not in keys:
            keys += ["SubhaloHalfmassRadType"]
            factors += [self.units("l/h")]
        else:
            factors.insert(keys.index("SubhaloHalfmassRadType"), self.units("l/h"))
    if assignment_func not in [voronoi_RT_2D, voronoi_NGP_2D]:
        raise ValueError(f"Assignment function `{assignment_func.__name__}` not compatible.")
    uarrs = self.get_mapping_arrays(keys=keys, factors=factors, verbose=verbose)
    if rot is not None:
        rot_op = R.y(rot[1]) * R.x(rot[0])
        uarrs[0] = rot_op(uarrs[0]).astype(np.float32)
        uarrs[3] = rot_op(uarrs[3]).astype(np.float32)
    # Ngids = uarrs[0].shape[0]
    hmr = uarrs[4][0] if use_half_mass_rad else None  # always use gas (p_idx=0) hmr
    # hmr = (uarrs[4][obj.p_idx] if use_half_mass_rad else None)
    idcs, limits = indices_within_box(uarrs[0], uarrs[3], radius=hmr, verbose=verbose)
    hmr2 = limits[3] - limits[0]
    args = strip_ap_units(*uarrs[:3], *limits[:2], hmr2, mask=idcs, units=[L, M, T * M])
    if tracers is not None:
        args.append(tracers)
    else:
        args.append(xaxis)
    if divisions is not None:
        args.append(divisions)
    else:
        args.append(yaxis)
    if verbose:
        print(f"Raytracing particles with `{assignment_func.__name__}`...")
    try:
        assignment_func(projected, *args, periodic, verbose=True)
    except Exception as e:
        print(e)
        return projected * uarrs[1].unit / L**2
    projected = projected * uarrs[1].unit / L**2
    return projected, hmr2 / 2 * np.array([-1, 1, -1, 1]), idcs.shape[0]

get_mapping_arrays 

get_mapping_arrays(
    keys: list[str] | None = None,
    factors: list[float] = None,
    verbose: bool = False,
) -> list

Fetch data (arrays or scalars) from a TNGGalaxy object given keys.

Parameters:

Name Type Description Default
keys list[str] | None

Keys to fetch data.

 None
factors list[float]

Factors for modifying the fetched data.

 None
verbose bool

If True, print status updates to command line.

 False

Returns:

Type Description
list

List of fetched data arrays or scalars.
Source code in skais_mapper/simobjects.py 
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
def get_mapping_arrays(
    self,
    keys: list[str] | None = None,
    factors: list[float] = None,
    verbose: bool = False,
) -> list:
    """Fetch data (arrays or scalars) from a TNGGalaxy object given keys.

    Args:
        keys: Keys to fetch data.
        factors: Factors for modifying the fetched data.
        verbose: If True, print status updates to command line.

    Returns:
        (list): List of fetched data arrays or scalars.
    """
    if keys is None:
        keys = ["particle_positions", "masses", "radii"]
    if factors is None:
        factors = [1, 1, 3]
    vals = []
    for key, f in zip(keys, factors):
        if isinstance(key, tuple | list):
            try:
                v1 = getattr(self, key[0])
            except Exception:
                v1 = self.subhalo[key[0]]
            try:
                v2 = getattr(self, key[1])
            except Exception:
                v2 = self.subhalo[key[1]]
            v = v1 * v2 * f
        else:
            try:
                v = getattr(self, key) * f
            except Exception:
                v = self.subhalo[key] * f
        vals.append(v)
    if verbose:
        print(f"Loading arrays: {keys}...")
    return [*vals]

load_cosmology 

load_cosmology(
    cosmo_pars: dict | None = None, in_place: bool = True
) -> dict

Load the cosmological parameters from the header dictionary.

Note: This overrides the parent class method.

Parameters:

Name Type Description Default
cosmo_pars dict | None

The default dictionary for the CosmoModel dataclass.

 None
in_place bool

Load CosmoModel into instanace property cosmology directly.

 True

Returns:

Type Description
dict

a dictionary with cosmological parameters pulled the header.
Source code in skais_mapper/simobjects.py 
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
def load_cosmology(self, cosmo_pars: dict | None = None, in_place: bool = True) -> dict:
    """Load the cosmological parameters from the header dictionary.

    Note: This overrides the parent class method.

    Args:
        cosmo_pars: The default dictionary for the CosmoModel dataclass.
        in_place: Load CosmoModel into instanace property `cosmology` directly.

    Returns:
        (dict): a dictionary with cosmological parameters pulled the header.
    """
    cosmo_pars = {} if cosmo_pars is None else cosmo_pars
    cosmo_pars.setdefault("omega_l", self.header["snapshot"].get("OmegaLambda", 0.6911))
    cosmo_pars.setdefault("omega_m", self.header["snapshot"].get("Omega0", 0.3089))
    cosmo_pars.setdefault("omega_k", self.header["snapshot"].get("OmegaK", 0))
    cosmo_pars.setdefault("z", self.header["snapshot"].get("Redshift", None))
    cosmo_pars.setdefault("z", 1.0 / self.header["snapshot"].get("Time", 1) - 1)
    cosmo_pars.setdefault("h", self.header["snapshot"].get("HubbleParam", 0.6774))
    if in_place:
        self.cosmology = CosmoModel(**cosmo_pars)
    return cosmo_pars

load_data 

load_data(
    halo_index: int | None = None,
    primary_fields: list[str] | None = None,
    optional_fields: list[str] | None = None,
    particle_type: str | None = None,
    **kwargs,
)

Load the snapshot's relevant subset data of the specified group ID.

Parameters:

Name Type Description Default
halo_index int | None

Galaxy/Halo ID in the simulation.

 None
primary_fields list[str] | None

The primary fields to load (less fields load faster)

 None
optional_fields list[str] | None

The secondary fields to load (less fields load faster)

 None
particle_type str | None

Manually set the particle type.

 None
kwargs
  • verbose (bool): Print information to the command-line.
 {}

Returns:

Type Description
dict

dictionary containing the relevant snapshot data
Source code in skais_mapper/simobjects.py 
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
def load_data(
    self,
    halo_index: int | None = None,
    primary_fields: list[str] | None = None,
    optional_fields: list[str] | None = None,
    particle_type: str | None = None,
    **kwargs,
):
    """Load the snapshot's relevant subset data of the specified group ID.

    Args:
        halo_index: Galaxy/Halo ID in the simulation.
        primary_fields: The primary fields to load (less fields load faster)
        optional_fields: The secondary fields to load (less fields load faster)
        particle_type: Manually set the particle type.
        kwargs:
          - verbose (bool): Print information to the command-line.

    Returns:
        (dict): dictionary containing the relevant snapshot data
    """
    kwargs.setdefault("verbose", False)
    verbose = kwargs.pop("verbose")
    if particle_type is not None:
        p_idx = tng.util.pidx_from_ptype(particle_type)
    else:
        p_idx = self.p_idx
        particle_type = self.particle_type
    if primary_fields is None:
        primary_fields = self.primary_hdf5_fields[p_idx]
    if optional_fields is None:
        optional_fields = self.optional_hdf5_fields[p_idx]
    if halo_index is None or halo_index < 0 or halo_index >= len(self.subhalo_ids):
        halo_index = self.halo_index
    subset = tng.snapshots.snapshot_offsets(self.base_path, self.snapshot, halo_index, "Group")
    if verbose:
        print(f"Loading primary fields from snapshot [{particle_type}]: {self.snapshot}")
    data = tng.snapshots.load_snapshot(
        self.base_path,
        self.snapshot,
        particle_type,
        subset=subset,
        fields=primary_fields,
        **kwargs,
    )
    try:
        if verbose:
            print(f"Loading optional fields from snapshot [{particle_type}]: {self.snapshot}")
        data.update(
            tng.snapshots.load_snapshot(
                self.base_path,
                self.snapshot,
                particle_type,
                subset=subset,
                fields=optional_fields,
                as_array=False,
                **kwargs,
            )
        )
    except Exception as ex:
        if ("Particle type" in ex.args[0]) and ("does not have field" in ex.args[0]):
            for f in optional_fields:
                data[f] = None
        else:
            raise
    if particle_type == "gas" and ("GFM_Metals" not in data or data["GFM_Metals"] is None):
        data["GFM_Metals"] = np.array([[0.76]])
    if "CenterOfMass" not in data or data["CenterOfMass"] is None:
        data["CenterOfMass"] = data["Coordinates"]
    return data

load_header 

load_header() -> dict

Load the header (overrides parent class method).

Returns:

Type Description
dict

a dictionary containing the headers from a snapshot and the corresponding group catalog.
Source code in skais_mapper/simobjects.py 
531
532
533
534
535
536
537
538
539
540
541
def load_header(self) -> dict:
    """Load the header (overrides parent class method).

    Returns:
        (dict): a dictionary containing the headers from a snapshot and the
          corresponding group catalog.
    """
    header = {}
    header["catalog"] = tng.groupcat.load_header(self.base_path, self.snapshot)
    header["snapshot"] = tng.snapshots.load_header(self.base_path, self.snapshot)
    return header

load_subhalo 

load_subhalo(
    halo_index: int = 0, verbose: bool = False
) -> dict

Load the FOF subhalo metadata of a specified group ID.

Parameters:

Name Type Description Default
halo_index int

the group ID corresponding to a subhalo ID

 0
verbose bool

Print status information to the command-line.

 False

Returns:

Type Description
dict

dictionary containing the FOF subhalo metadata
Source code in skais_mapper/simobjects.py 
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
def load_subhalo(self, halo_index: int = 0, verbose: bool = False) -> dict:
    """Load the FOF subhalo metadata of a specified group ID.

    Args:
        halo_index: the group ID corresponding to a subhalo ID
        verbose: Print status information to the command-line.

    Returns:
        (dict): dictionary containing the FOF subhalo metadata
    """
    if halo_index < 0 or halo_index >= len(self.subhalo_ids):
        halo_index = 0
    sh_id = self.subhalo_ids[halo_index]
    if verbose:
        print(f"Loading group: {sh_id}")
    data = tng.groupcat.load_single(self.base_path, self.snapshot, subhalo_id=sh_id)
    return data

p_idx 

p_idx(p_idx: int | str)

Setter for particle index.

Source code in skais_mapper/simobjects.py 
526
527
528
529
@SPHGalaxy.p_idx.setter
def p_idx(self, p_idx: int | str):
    """Setter for particle index."""
    SPHGalaxy.p_idx.fset(self, p_idx)

subhalo_list staticmethod 

subhalo_list(
    base_path: str,
    snapshot: int,
    filtered: bool = True,
    verbose: bool = False,
) -> NDArray

Pull a list of subhalo IDs from the snapshot's FOF group catalog.

Note: the indices of the list are halo IDs.

Parameters:

Name Type Description Default
base_path str

Base path to the Illustris(TNG) snapshots.

 required
snapshot int

Snapshot ID {0-99}.

 required
filtered bool

If True, negative subhalo IDs are filtered out.

 True
verbose bool

If True, print out results to the command line.

 False

Returns:

Type Description
list | np.ndarray

a list of subhalo i.e. galaxy IDs.
Source code in skais_mapper/simobjects.py 
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
@staticmethod
def subhalo_list(
    base_path: str,
    snapshot: int,
    filtered: bool = True,
    verbose: bool = False,
) -> NDArray:
    """Pull a list of subhalo IDs from the snapshot's FOF group catalog.

    Note: the indices of the list are halo IDs.

    Args:
        base_path: Base path to the Illustris(TNG) snapshots.
        snapshot: Snapshot ID {0-99}.
        filtered: If True, negative subhalo IDs are filtered out.
        verbose: If True, print out results to the command line.

    Returns:
        (list | np.ndarray): a list of subhalo i.e. galaxy IDs.
    """
    kwargs = {"fields": ["GroupFirstSub"]}
    if verbose:
        print(f"Searching group catalog: snapshot {snapshot}")
    galaxy_list = tng.groupcat.load_halos(base_path, snapshot, as_array=True, **kwargs)
    galaxy_list = galaxy_list[galaxy_list >= 0] if filtered else galaxy_list
    if verbose:
        N_g = galaxy_list.max()
        N_g += 0 in galaxy_list
        print(f"Found {N_g} groups in catalog: snapshot {snapshot}")
    return galaxy_list

indices_within_box 

indices_within_box(
    pos: np.ndarray | au.Quantity,
    center: list | np.ndarray | au.Quantity,
    radius: float | au.Quantity = None,
    fraction: float = 1.0,
    verbose: bool = False,
) -> tuple[au.Quantity, list[au.Quantity]]

Get particle indices within a box of given radius from the centre, e.g. half-mass radius.

Parameters:

Name Type Description Default
pos np.ndarray | au.Quantity

Particle positions to be filtered.

 required
center list | np.ndarray | au.Quantity

Center position of the cube.

 required
radius float | au.Quantity

Radius, i.e. half-side of the cube.

 None
fraction float

Box fraction for the default radius if not given.

 1.0
verbose bool

If True, print status updates to command line.

 False

Returns:

Type Description
(au.Quantity, list[au.Quantity])

Indices of the filtered particle positions and their 3D extent.
Source code in skais_mapper/simobjects.py 
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
def indices_within_box(
    pos: np.ndarray | au.Quantity,
    center: list | np.ndarray | au.Quantity,
    radius: float | au.Quantity = None,
    fraction: float = 1.0,
    verbose: bool = False,
) -> tuple[au.Quantity, list[au.Quantity]]:
    """Get particle indices within a box of given radius from the centre, e.g. half-mass radius.

    Args:
        pos: Particle positions to be filtered.
        center: Center position of the cube.
        radius: Radius, i.e. half-side of the cube.
        fraction: Box fraction for the default radius if not given.
        verbose: If True, print status updates to command line.

    Returns:
        (au.Quantity, list[au.Quantity]):
          Indices of the filtered particle positions and their 3D extent.
    """
    pos_extent = (
        pos[:, 0].max() - pos[:, 0].min(),
        pos[:, 1].max() - pos[:, 1].min(),
        pos[:, 2].max() - pos[:, 2].min(),
    )
    w = max(pos_extent)
    if verbose:
        print(f"Box extent:  {pos_extent[0]}   {pos_extent[1]}   {pos_extent[2]}")
    if radius is None:
        radius = 0.25 * fraction * w
    else:
        radius *= fraction
    if verbose:
        print(f"Radius: {radius}")
    xmin, xmax = center[0] - radius, center[0] + radius
    ymin, ymax = center[1] - radius, center[1] + radius
    zmin, zmax = center[2] - radius, center[2] + radius
    if verbose:
        print(f"Extent: {xmax - xmin}   {ymax - ymin}   {zmax - zmin}")
    indices = np.where(
        (pos[:, 0] > xmin)
        & (pos[:, 0] < xmax)
        & (pos[:, 1] > ymin)
        & (pos[:, 1] < ymax)
        & (pos[:, 2] > zmin)
        & (pos[:, 2] < zmax)
    )[0]
    if verbose:
        print(f"Selected particles: {indices.shape[0]:,} / {pos.shape[0]:,}")
    return indices, [xmin, ymin, zmin, xmax, ymax, zmax]

strip_ap_units 

strip_ap_units(
    *args,
    mask: list | None = None,
    units: list[au.Unit] | None = None,
    dtype: Any = np.float32,
) -> list

Remove astropy units from data arrays or scalars.

Parameters:

Name Type Description Default
args

Data arrays or scalars with astropy units.

 ()
mask list | None

Mask for the data arrays.

 None
units list[au.Unit] | None

Astropy units to be stripped.

 None
dtype Any

Data type of the stripped data array.

 np.float32

Returns:

Type Description
list[np.ndarray]

Of astropy units stripped data arrays or scalars.
Source code in skais_mapper/simobjects.py 
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
def strip_ap_units(
    *args,
    mask: list | None = None,
    units: list[au.Unit] | None = None,
    dtype: Any = np.float32,
) -> list:
    """Remove astropy units from data arrays or scalars.

    Args:
        args: Data arrays or scalars with astropy units.
        mask: Mask for the data arrays.
        units: Astropy units to be stripped.
        dtype: Data type of the stripped data array.

    Returns:
        (list[np.ndarray]): Of astropy units stripped data arrays or scalars.
    """
    arg_ls = list(args)
    for i, arg in enumerate(arg_ls):
        if units is not None:
            for u in units:
                if arg_ls[i].unit.is_equivalent(u):
                    arg_ls[i] = arg_ls[i].to(u)
        arg_ls[i] = arg_ls[i].value
        if isinstance(arg_ls[i], np.ndarray):
            if mask is not None:
                arg_ls[i] = arg_ls[i][mask].astype(dtype)
    return arg_ls