Source code for ctapipe.tools.dump_instrument
"""
Dump instrumental descriptions in a monte-carlo (simtelarray) input file to
FITS files that can be loaded independently (e.g. with
CameraGeometry.from_table()). The name of the output files are
automatically generated.
"""
import os
import pathlib
from ..compat import ECSV_FMT
from ..core import Provenance, Tool
from ..core.traits import Enum, Path, Unicode
from ..exceptions import InputMissing
from ..io import EventSource
__all__ = ["DumpInstrumentTool"]
[docs]
class DumpInstrumentTool(Tool):
description = Unicode(__doc__)
name = "ctapipe-dump-instrument"
outdir = Path(
file_ok=False,
directory_ok=True,
allow_none=True,
default_value=None,
help="Output directory. If not given, the current working directory will be used.",
).tag(config=True)
format = Enum(
["fits", "ecsv", "hdf5", "service"],
default_value="fits",
help="Format of output file. 'service' creates CTAO service data format directory structure.",
config=True,
)
aliases = {
("i", "input"): "EventSource.input_url",
("f", "format"): "DumpInstrumentTool.format",
("o", "outdir"): "DumpInstrumentTool.outdir",
}
classes = [EventSource]
[docs]
def setup(self):
try:
with EventSource(parent=self) as source:
self.infile = source.input_url
self.subarray = source.subarray
except InputMissing:
self.log.critical(
"Specifying EventSource.input_url is required (via -i, --input or a config file)."
)
self.exit(1)
[docs]
def start(self):
if self.outdir is None:
self.outdir = pathlib.Path(os.getcwd())
self.outdir.mkdir(exist_ok=True, parents=True)
if self.format == "hdf5":
self.subarray.to_hdf(self.outdir / "subarray.h5")
elif self.format == "service":
self.write_service_data()
else:
self.write_camera_definitions()
self.write_optics_descriptions()
self.write_subarray_description()
@staticmethod
def _get_file_format_info(format_name):
"""returns file extension + dict of required parameters for
Table.write"""
if format_name == "fits":
return "fits.gz", dict()
elif format_name == "ecsv":
return "ecsv", dict()
else:
raise NameError(f"format {format_name} not supported")
[docs]
def write_camera_definitions(self):
"""writes out camgeom and camreadout files for each camera"""
self.subarray.info(printer=self.log.info)
for camera in self.subarray.camera_types:
self.write_single_camera(camera)
[docs]
def write_single_camera(self, camera, outdir=None, name_prefix=None):
"""Write out camera geometry and readout for a single camera.
Parameters
----------
camera : CameraDescription
CameraDescription object to write out
outdir : Path, optional
Directory to write files to. If None, uses self.outdir
name_prefix : str, optional
Prefix for output filenames. If None, uses camera.name
"""
outdir = self.outdir if outdir is None else outdir
name_prefix = camera.name if name_prefix is None else name_prefix
ext, args = self._get_file_format_info(self.format)
self.log.debug("Writing camera %s", camera)
geom = camera.geometry
readout = camera.readout
geom_table = geom.to_table()
geom_table.meta["SOURCE"] = str(self.infile)
geom_filename = outdir / f"{name_prefix}.camgeom.{ext}"
readout_table = readout.to_table()
readout_table.meta["SOURCE"] = str(self.infile)
readout_filename = outdir / f"{name_prefix}.camreadout.{ext}"
try:
geom_table.write(geom_filename, **args)
Provenance().add_output_file(geom_filename, "CameraGeometry")
except OSError as err:
self.log.exception("couldn't write camera geometry because: %s", err)
try:
readout_table.write(readout_filename, **args)
Provenance().add_output_file(readout_filename, "CameraReadout")
except OSError as err:
self.log.exception("couldn't write camera definition because: %s", err)
[docs]
def write_optics_descriptions(self):
"""writes out optics files for each telescope type"""
sub = self.subarray
ext, args = self._get_file_format_info(self.format)
tab = sub.to_table(kind="optics")
tab.meta["SOURCE"] = str(self.infile)
filename = self.outdir / f"{sub.name}.optics.{ext}"
try:
tab.write(filename, **args)
Provenance().add_output_file(filename, "OpticsDescription")
except OSError as err:
self.log.exception(
"couldn't write optics description '%s' because: %s", filename, err
)
[docs]
def write_subarray_description(self):
sub = self.subarray
ext, args = self._get_file_format_info(self.format)
tab = sub.to_table(kind="subarray", meta_convention="fits")
tab.meta["SOURCE"] = str(self.infile)
filename = self.outdir / f"{sub.name}.subarray.{ext}"
try:
tab.write(filename, **args)
Provenance().add_output_file(filename, "SubarrayDescription")
except OSError as err:
self.log.exception(
"couldn't write subarray description '%s' because: %s", filename, err
)
[docs]
def write_service_data(self, subarray_id=1, site=None):
"""
Write SubarrayDescription to service data directory structure.
This creates the directory structure and files, which can later be loaded with
`~ctapipe.instrument.SubarrayDescription.from_service_data`.
Parameters
----------
subarray_id : int, optional
Subarray ID to assign (default: 1)
site : str, optional
Site name (e.g., "CTAO-North", "CTAO-South").
If not provided, it will be inferred from the reference location.
"""
import json
from astropy.table import QTable
from ctapipe.io import metadata as meta
sub = self.subarray
self.outdir.mkdir(exist_ok=True, parents=True)
# Create instrument directory to match CTAO service data structure
instrument_dir = self.outdir / "instrument"
instrument_dir.mkdir(exist_ok=True, parents=True)
self.log.info(
"Writing instrument description in CTAO service data format to %s",
self.outdir,
)
# Infer site from coordinates if not provided
if site is None:
# Simple heuristic based on latitude
lat = sub.reference_location.geodetic.lat.value
if lat > 0:
site = "CTAO-North"
else:
site = "CTAO-South"
# Build reference metadata to embed in instrument.meta.json
activity = Provenance().current_activity
activity_meta = (
meta.Activity.from_provenance(activity.provenance)
if activity is not None
else meta.Activity()
)
instrument_reference = meta.Reference(
contact=meta.Contact(),
product=meta.Product(
description=f"Instrument description for {sub.name}",
data_category="Other",
data_association="Subarray",
data_model_name="CTAO Service Data",
data_model_version=sub.CURRENT_SERVICE_DATA_VERSION,
data_model_url="",
format="json",
),
process=meta.Process(),
activity=activity_meta,
instrument=meta.Instrument(
site=site,
class_="Subarray",
),
)
# Create instrument.meta.json
meta_file = instrument_dir / "instrument.meta.json"
with open(meta_file, "w") as f:
json.dump(instrument_reference.to_dict(), f, indent=2)
Provenance().add_output_file(meta_file, "ServiceDataMeta")
# Create array-element-ids.json
ae_reference = meta.Reference(
contact=meta.Contact(),
product=meta.Product(
description=f"Array element IDs for {sub.name}",
data_category="Other",
data_association="Subarray",
data_model_name="ctao.common.identifiers.array_elements",
data_model_version=sub.CURRENT_ARRAY_ELEMENTS_IDENTIFIERS_VERSION,
data_model_url="https://gitlab.cta-observatory.org/cta-computing/common/identifiers",
format="json",
),
process=meta.Process(),
activity=activity_meta,
instrument=meta.Instrument(site=site),
)
array_element_ids = {
"metadata": ae_reference.to_dict(),
"array_elements": [
{"id": int(tel_id), "name": f"TEL{tel_id:03d}"}
for tel_id, tel in sub.tels.items()
],
}
ae_ids_file = instrument_dir / "array-element-ids.json"
with open(ae_ids_file, "w") as f:
json.dump(array_element_ids, f, indent=2)
Provenance().add_output_file(ae_ids_file, "ServiceDataArrayElements")
# Create subarray-ids.json
subarray_reference = meta.Reference(
contact=meta.Contact(),
product=meta.Product(
description=f"Subarray IDs for {sub.name}",
data_category="Other",
data_association="Subarray",
data_model_name="ctao.common.identifiers.subarrays",
data_model_version=sub.CURRENT_SUBARRAY_IDENTIFIERS_VERSION,
data_model_url="https://gitlab.cta-observatory.org/cta-computing/common/identifiers",
format="json",
),
process=meta.Process(),
activity=activity_meta,
instrument=meta.Instrument(site=site),
)
subarray_ids = {
"metadata": subarray_reference.to_dict(),
"subarrays": [
{
"id": subarray_id,
"name": sub.name,
"site": site,
"array_element_ids": [int(tel_id) for tel_id in sub.tel_ids],
}
],
}
subarray_ids_file = instrument_dir / "subarray-ids.json"
with open(subarray_ids_file, "w") as f:
json.dump(subarray_ids, f, indent=2)
Provenance().add_output_file(subarray_ids_file, "ServiceDataSubarrays")
# Create positions directory and file
positions_dir = instrument_dir / "positions"
positions_dir.mkdir(exist_ok=True)
# Get reference location in ITRS coordinates
itrs = sub.reference_location.itrs
# Create positions table
positions_table = QTable(
{
"ae_id": [int(tel_id) for tel_id in sub.tel_ids],
"name": [tel.name for tel in sub.tels.values()],
"x": [sub.positions[tel_id][0] for tel_id in sub.tel_ids],
"y": [sub.positions[tel_id][1] for tel_id in sub.tel_ids],
"z": [sub.positions[tel_id][2] for tel_id in sub.tel_ids],
}
)
positions_table.meta["reference_x"] = str(itrs.x)
positions_table.meta["reference_y"] = str(itrs.y)
positions_table.meta["reference_z"] = str(itrs.z)
positions_table.meta["site"] = site
positions_reference = meta.Reference(
contact=meta.Contact(),
product=meta.Product(
description=f"Array element positions for {sub.name}",
data_category="Other",
data_association="Subarray",
format="ecsv",
),
process=meta.Process(),
activity=activity_meta,
instrument=meta.Instrument(site=site),
)
positions_table.meta.update(positions_reference.to_dict())
positions_file = (
positions_dir / f"{site.replace(' ', '_')}_ArrayElementPositions.ecsv"
)
positions_table.write(positions_file, format=ECSV_FMT, overwrite=True)
Provenance().add_output_file(positions_file, "ServiceDataPositions")
# Group telescopes by unique TelescopeDescription.
# str(tel_desc) = "{size_type}_{optics_name}_{camera_name}" and is used
# as the shared type-directory name.
array_elements_dir = instrument_dir / "array-elements"
array_elements_dir.mkdir(exist_ok=True, parents=True)
# Write shared files once per unique telescope type
for index, tel_desc in enumerate(sub.telescope_types):
type_key = f"TEL-TYPE-{index:02d}"
type_dir = array_elements_dir / type_key
type_dir.mkdir(exist_ok=True, parents=True)
self.log.debug("Writing shared type directory %s", type_key)
optics_file = type_dir / f"{type_key}.tel_optics.ecsv"
optics_table = tel_desc.optics.to_table()
optics_table.write(optics_file, format=ECSV_FMT, overwrite=True)
Provenance().add_output_file(optics_file, "ServiceDataOptics")
orig_format = self.format
self.format = "fits"
try:
self.write_single_camera(
tel_desc.camera, outdir=type_dir, name_prefix=type_key
)
finally:
self.format = orig_format
# Create a real directory for each array element with per-file symlinks
# pointing to the shared type directory (deduplication).
for tel_id, tel_desc in sub.tel.items():
ae_id_str = f"{tel_id:03d}"
ae_dir = array_elements_dir / ae_id_str
ae_dir.mkdir(exist_ok=True, parents=True)
index = sub.telescope_types.index(tel_desc)
type_key = f"TEL-TYPE-{index:02d}"
self.log.debug(
"Writing array element %s -> %s (symlinks)", ae_id_str, type_key
)
for suffix in ["tel_optics.ecsv", "camgeom.fits.gz", "camreadout.fits.gz"]:
link = ae_dir / f"{ae_id_str}.{suffix}"
link.symlink_to(f"../{type_key}/{type_key}.{suffix}")
self.log.info("Service data written successfully to %s", self.outdir)
def main():
tool = DumpInstrumentTool()
tool.run()
if __name__ == "__main__":
main()