Particle I/O and Catalogue Conversion#

This notebook follows one catalogue through the complete PyHermes input boundary: download the original Quijote group_tab file, read it directly with the native FoF reader, convert the halo catalogue to NPZ and BIN, and read both converted files back through the shared particle interface.

Start with quick_start.ipynb for the shortest complete calculation. Return here when adapting a catalogue of your own.

What this notebook emphasizes#

  • downloading and caching one complete FoF group_tab file

  • reading a Quijote/Pylians FoF file directly

  • constructing and reading an equivalent local NPZ catalogue

  • constructing and reading an equivalent local raw BIN table

  • converging different source formats on one PyHermes particle dictionary

[1]:

from pathlib import Path
import os

import numpy as np

from pyhermes.io import read_particle_data
from pyhermes.param.parambase import read_param

cwd = Path.cwd().resolve()
if cwd.name == "notebooks" and cwd.parent.name == "examples":
    examples_dir = cwd.parent
elif cwd.name == "examples":
    examples_dir = cwd
elif (cwd / "examples").is_dir():
    examples_dir = cwd / "examples"
else:
    raise RuntimeError("Run this notebook from the project root, examples/, or examples/notebooks/.")

os.chdir(examples_dir)
data_dir = Path("./data")
print(f"Working directory: {Path.cwd()}")

Working directory: /Users/xutengpeng/xutp/PycharmProjects/PyHermes/examples

1. Download and read the original FoF catalogue#

The public group_tab_004.0 header reports Nfiles = 1, so it is the complete halo catalogue rather than the first piece of a split dataset. This notebook reuses the exact fin block from the Quick Start configuration. read_particle_data downloads the file into the requested cache path, verifies the SHA256 digest, and then passes it directly to the FoF reader. Later runs reuse the cached file.

[2]:

projection_config = read_param("./configs/param_sfc_projection.yaml")["SFCProjection"]
catalog_fin = projection_config["fin"]
fof_url = catalog_fin["path"]
fof_format = catalog_fin["format"]
fof_download = catalog_fin["download"]
group_tab_path = Path(fof_download["cache_path"])
groups_dir = group_tab_path.parent

print(f"FoF source: {fof_url}")
print(f"Local cache: {group_tab_path}")

11:35:24 - INFO - pyhermes.param.parambase:ParamBase - Reading configure file: './configs/param_sfc_projection.yaml'
11:35:24 - INFO - pyhermes.param.parambase:ParamBase - Input parameter file format is YAML
FoF source: https://pyhermes.astroslacker.com/downloads/group_tab_004.0
Local cache: data/quijote_halos/8000/groups_004/group_tab_004.0

[3]:

halo_data = read_particle_data(
    fof_url,
    data_format=fof_format,
    download=fof_download,
    redshift=0.0,
    fields={
        "vx": "vel_x",
        "vy": "vel_y",
        "vz": "vel_z",
        "mass": "mass",
        "npart": "npart",
    },
)

print(f"haloes: {halo_data['size']:,}")
print(f"positions: {halo_data['pos'].shape} {halo_data['pos'].dtype}")
print(f"fields: {sorted(key for key in halo_data if key not in {'pos', 'size'})}")

11:35:24 - INFO - pyhermes.io.resources:resolve_data_path - Using cached particle data: data/quijote_halos/8000/groups_004/group_tab_004.0
11:35:24 - INFO - pyhermes.io.readers:read_particle_data - Selected input particle format: fof
11:35:24 - INFO - pyhermes.io.readers:read_fof - Reading Quijote FoF halo data from ---> data/quijote_halos/8000/groups_004/group_tab_004.0 <---
haloes: 406,728
positions: (406728, 3) float32
fields: ['mass', 'npart', 'vx', 'vy', 'vz']

2. Construct and read an NPZ catalogue#

NPZ is one of the catalogue formats supported by PyHermes. Here we construct a local NPZ copy from the FoF data to make its array contract explicit:

  • pos: comoving position in \(h^{-1}\mathrm{Mpc}\)

  • vel_x, vel_y, vel_z: peculiar velocity components in \(\mathrm{km}\,\mathrm{s}^{-1}\)

  • mass: halo mass in \(h^{-1}M_\odot\)

  • npart: number of simulation particles in the FoF group

All arrays are stored as float32. This NPZ is a local conversion example, not a second public download. It is placed beside the cached source group_tab file so the provenance and conversion products stay together.

[4]:

npz_path = groups_dir / "quijote_halos_8000_snap004_converted.npz"
np.savez_compressed(
    npz_path,
    pos=np.ascontiguousarray(halo_data["pos"], dtype=np.float32),
    vel_x=np.ascontiguousarray(halo_data["vx"], dtype=np.float32),
    vel_y=np.ascontiguousarray(halo_data["vy"], dtype=np.float32),
    vel_z=np.ascontiguousarray(halo_data["vz"], dtype=np.float32),
    mass=np.ascontiguousarray(halo_data["mass"], dtype=np.float32),
    npart=np.ascontiguousarray(halo_data["npart"], dtype=np.float32),
)
print(f"Wrote {npz_path} ({npz_path.stat().st_size / 2**20:.2f} MiB)")

Wrote data/quijote_halos/8000/groups_004/quijote_halos_8000_snap004_converted.npz (8.55 MiB)

[5]:

npz_fields = {
    "vx": "vel_x",
    "vy": "vel_y",
    "vz": "vel_z",
    "mass": "mass",
    "npart": "npart",
}
npz_data = read_particle_data(
    npz_path,
    data_format="npz",
    fields=npz_fields,
)

for name in ("pos", "vx", "vy", "vz", "mass", "npart"):
    np.testing.assert_allclose(npz_data[name], halo_data[name])

units = {
    "pos": "Mpc/h",
    "vel_x": "km/s",
    "vel_y": "km/s",
    "vel_z": "km/s",
    "mass": "Msun/h",
    "npart": "particle count",
}
with np.load(npz_path) as archive:
    for name in archive.files:
        values = archive[name]
        print(f"{name:>6s}: shape={values.shape}, dtype={values.dtype}, unit={units[name]}")

11:35:25 - INFO - pyhermes.io.readers:read_particle_data - Selected input particle format: npz
11:35:25 - INFO - pyhermes.io.readers:read_npz - Reading particle data from ---> data/quijote_halos/8000/groups_004/quijote_halos_8000_snap004_converted.npz <---
   pos: shape=(406728, 3), dtype=float32, unit=Mpc/h
 vel_x: shape=(406728,), dtype=float32, unit=km/s
 vel_y: shape=(406728,), dtype=float32, unit=km/s
 vel_z: shape=(406728,), dtype=float32, unit=km/s
  mass: shape=(406728,), dtype=float32, unit=Msun/h
 npart: shape=(406728,), dtype=float32, unit=particle count

3. Construct and read a BIN catalogue#

A raw BIN file stores values but no names, shape, dtype, or units. We therefore define its eight-column layout beside both the writer and reader: (x, y, z, vx, vy, vz, mass, npart). Keeping this mapping visible is the binary format contract.

[6]:

bin_path = groups_dir / "quijote_halos_8000_snap004_converted.bin"
bin_table = np.column_stack(
    [
        halo_data["pos"],
        halo_data["vx"],
        halo_data["vy"],
        halo_data["vz"],
        halo_data["mass"],
        halo_data["npart"],
    ]
).astype(np.float32, copy=False)
bin_table.tofile(bin_path)
print(f"Wrote {bin_path} with shape {bin_table.shape}")

Wrote data/quijote_halos/8000/groups_004/quijote_halos_8000_snap004_converted.bin with shape (406728, 8)

[7]:

bin_data = read_particle_data(
    bin_path,
    data_format="bin",
    dtype="float32",
    ncols=8,
    pos_cols=[0, 1, 2],
    fields={
        "vx": 3,
        "vy": 4,
        "vz": 5,
        "mass": 6,
        "npart": 7,
    },
)

for name in ("pos", "vx", "vy", "vz", "mass", "npart"):
    np.testing.assert_allclose(bin_data[name], halo_data[name])
print(f"Reloaded {bin_data['size']:,} haloes from the BIN table")

11:35:25 - INFO - pyhermes.io.readers:read_particle_data - Selected input particle format: bin
11:35:25 - INFO - pyhermes.io.readers:read_bin - Reading particle data from ---> data/quijote_halos/8000/groups_004/quijote_halos_8000_snap004_converted.bin <---
Reloaded 406,728 haloes from the BIN table

4. Other native simulation readers#

The same dispatcher supports format="gadget" for legacy binary snapshots, format="gadget_hdf5" for HDF5 snapshots, and format="gadget-fof" for legacy Gadget FoF catalogues. A direct Quijote group_tab path can describe one complete file, as above, or the .0 member of a local split catalogue; in the latter case all numerically suffixed siblings must be present. Directory input remains available through a catalogue root plus snapnum. Position, velocity, and mass conversion factors belong in reader_params; see the I/O reference for the complete options.

5. Hand the catalogue to SFCProjection#

The FoF, NPZ, and BIN routes now provide the same objects: positions plus optional catalogue weights and field values. Any of them can be passed directly to an in-memory projection.

task.particle_pos = npz_data["pos"]
task.field_value = npz_data["mass"]

Continue to sfc_projection.ipynb to build number-, mass-, and redshift-space SFCField objects.