Skip to content

Commit

Permalink
Define and use an Enum for the focal length kind
Browse files Browse the repository at this point in the history
  • Loading branch information
maxnoe committed Jul 8, 2022
1 parent 687f38c commit 9692626
Show file tree
Hide file tree
Showing 20 changed files with 205 additions and 160 deletions.
19 changes: 11 additions & 8 deletions ctapipe/conftest.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,13 @@
from copy import deepcopy

import pytest
from pytest_astropy_header.display import PYTEST_HEADER_MODULES

from ctapipe.instrument import CameraGeometry
from ctapipe.io import SimTelEventSource
from ctapipe.utils import get_dataset_path
from ctapipe.utils.filelock import FileLock

from pytest_astropy_header.display import PYTEST_HEADER_MODULES

PYTEST_HEADER_MODULES.clear()
PYTEST_HEADER_MODULES["eventio"] = "eventio"
PYTEST_HEADER_MODULES["numpy"] = "numpy"
Expand Down Expand Up @@ -54,7 +53,9 @@ def _global_example_event():
print("******************** LOAD TEST EVENT ***********************")

# FIXME: switch to prod5b+ file that contains effective focal length
with SimTelEventSource(input_url=filename, focal_length_choice="nominal") as reader:
with SimTelEventSource(
input_url=filename, focal_length_choice="EQUIVALENT"
) as reader:
event = next(iter(reader))

return event
Expand All @@ -69,7 +70,9 @@ def example_subarray():

print("******************** LOAD TEST EVENT ***********************")

with SimTelEventSource(input_url=filename, focal_length_choice="nominal") as reader:
with SimTelEventSource(
input_url=filename, focal_length_choice="EQUIVALENT"
) as reader:
return reader.subarray


Expand All @@ -94,7 +97,7 @@ def _subarray_and_event_gamma_off_axis_500_gev():

path = get_dataset_path("lst_prod3_calibration_and_mcphotons.simtel.zst")

with SimTelEventSource(path, focal_length_choice="nominal") as source:
with SimTelEventSource(path, focal_length_choice="EQUIVALENT") as source:
it = iter(source)
# we want the second event, first event is a corner clipper
next(it)
Expand Down Expand Up @@ -288,8 +291,8 @@ def dl1_by_type_file(dl1_tmp_path, prod5_gamma_simtel_path):
"""
DL1 file containing both images and parameters from a gamma simulation set.
"""
from ctapipe.tools.process import ProcessorTool
from ctapipe.core import run_tool
from ctapipe.tools.process import ProcessorTool

output = dl1_tmp_path / "gamma_by_type.dl1.h5"

Expand Down Expand Up @@ -384,7 +387,7 @@ def dl1_muon_file(dl1_tmp_path):
"--write-images",
"--DataWriter.write_parameters=False",
"--DataWriter.Contact.name=αℓℓ the äüöß",
"--SimTelEventSource.focal_length_choice=nominal",
"--SimTelEventSource.focal_length_choice=EQUIVALENT",
]
assert run_tool(ProcessorTool(), argv=argv, cwd=dl1_tmp_path) == 0
return output
Expand All @@ -395,8 +398,8 @@ def dl1_proton_file(dl1_tmp_path, prod5_proton_simtel_path):
"""
DL1 file containing images and parameters for a prod5 proton run
"""
from ctapipe.tools.process import ProcessorTool
from ctapipe.core import run_tool
from ctapipe.tools.process import ProcessorTool

output = dl1_tmp_path / "proton.dl1.h5"

Expand Down
21 changes: 11 additions & 10 deletions ctapipe/core/traits.py
Original file line number Diff line number Diff line change
@@ -1,17 +1,17 @@
"""
Traitlet implementations for ctapipe
"""
import copy
import os
import pathlib
from collections import UserList
from fnmatch import fnmatch
from typing import Optional
import copy
from astropy.time import Time
import pathlib
from urllib.parse import urlparse
import os

import traitlets
import traitlets.config
from astropy.time import Time
from traitlets import Undefined

from .component import non_abstract_children
Expand Down Expand Up @@ -65,17 +65,18 @@
Set = traitlets.Set
CRegExp = traitlets.CRegExp
CaselessStrEnum = traitlets.CaselessStrEnum
UseEnum = traitlets.UseEnum
TraitError = traitlets.TraitError
TraitType = traitlets.TraitType
observe = traitlets.observe
flag = traitlets.config.boolean_flag


class AstroTime(TraitType):
""" A trait representing a point in Time, as understood by `astropy.time`"""
"""A trait representing a point in Time, as understood by `astropy.time`"""

def validate(self, obj, value):
""" try to parse and return an ISO time string """
"""try to parse and return an ISO time string"""
try:
the_time = Time(value)
the_time.format = "iso"
Expand Down Expand Up @@ -273,7 +274,7 @@ def __init__(self, *args):

@property
def tel(self):
""" access the value per telescope_id, e.g. `param.tel[2]`"""
"""access the value per telescope_id, e.g. `param.tel[2]`"""
if self._lookup:
return self._lookup
else:
Expand Down Expand Up @@ -511,23 +512,23 @@ def set(self, obj, value):


class FloatTelescopeParameter(TelescopeParameter):
""" a `~ctapipe.core.traits.TelescopeParameter` with Float trait type"""
"""a `~ctapipe.core.traits.TelescopeParameter` with Float trait type"""

def __init__(self, **kwargs):
"""Create a new IntTelescopeParameter"""
super().__init__(trait=Float(), **kwargs)


class IntTelescopeParameter(TelescopeParameter):
""" a `~ctapipe.core.traits.TelescopeParameter` with Int trait type"""
"""a `~ctapipe.core.traits.TelescopeParameter` with Int trait type"""

def __init__(self, **kwargs):
"""Create a new IntTelescopeParameter"""
super().__init__(trait=Int(), **kwargs)


class BoolTelescopeParameter(TelescopeParameter):
""" a `~ctapipe.core.traits.TelescopeParameter` with Bool trait type"""
"""a `~ctapipe.core.traits.TelescopeParameter` with Bool trait type"""

def __init__(self, **kwargs):
"""Create a new BoolTelescopeParameter"""
Expand Down
10 changes: 5 additions & 5 deletions ctapipe/instrument/__init__.py
Original file line number Diff line number Diff line change
@@ -1,10 +1,9 @@
from .camera import CameraDescription, CameraGeometry, CameraReadout, PixelShape
from .atmosphere import get_atmosphere_profile_functions
from .telescope import TelescopeDescription
from .optics import OpticsDescription
from .subarray import SubarrayDescription, UnknownTelescopeID
from .camera import CameraDescription, CameraGeometry, CameraReadout, PixelShape
from .guess import guess_telescope

from .optics import FocalLengthKind, OpticsDescription
from .subarray import SubarrayDescription, UnknownTelescopeID
from .telescope import TelescopeDescription

__all__ = [
"CameraDescription",
Expand All @@ -17,4 +16,5 @@
"SubarrayDescription",
"TelescopeDescription",
"UnknownTelescopeID",
"FocalLengthKind",
]
22 changes: 20 additions & 2 deletions ctapipe/instrument/optics.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
"""

import logging
from enum import Enum, auto, unique

import astropy.units as u
import numpy as np
Expand All @@ -12,6 +13,23 @@
logger = logging.getLogger(__name__)


@unique
class FocalLengthKind(Enum):
"""
Enumeration for the different kinds.
"""

#: Effective focal length computed from ray tracing a point source
#: and calculating the off-axis center of mass of the light distribution.
#: This focal length should be used in coordinate transforms between camera
#: frame and telescope frame to correct for the mean effect of coma abberation.
EFFECTIVE = auto()
#: Equivalent focal length is the nominal focal length of the main reflector
#: for single mirror telescopes and the thin-lens equivalent for dual mirror
#: telescopes.
EQUIVALENT = auto()


class OpticsDescription:
"""
Describes the optics of a Cherenkov Telescope mirror
Expand All @@ -31,8 +49,8 @@ class OpticsDescription:
mirror telescopes.
effective_focal_length : astropy.units.Quantity[length]
The effective_focal_length is the focal length estimated from
ray tracing to correct for koma aberration. It is thus not automatically
available for all simulations, but only if it was set before hand
ray tracing to correct for coma aberration. It is thus not automatically
available for all simulations, but only if it was set beforehand
in the simtel configuration. This is the focal length that should be
used for transforming from camera frame to telescope frame for all
reconstruction tasks to correct for the mean aberration.
Expand Down
13 changes: 9 additions & 4 deletions ctapipe/instrument/subarray.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@
from .. import __version__ as CTAPIPE_VERSION
from ..coordinates import CameraFrame, GroundFrame
from .camera import CameraDescription, CameraGeometry, CameraReadout
from .optics import OpticsDescription
from .optics import FocalLengthKind, OpticsDescription
from .telescope import TelescopeDescription

__all__ = ["SubarrayDescription"]
Expand Down Expand Up @@ -547,10 +547,13 @@ def to_hdf(self, h5file, overwrite=False, mode="a"):
meta["reference_itrs_z"] = itrs.z.to_value(u.m)

@classmethod
def from_hdf(cls, path, focal_length_choice="effective"):
def from_hdf(cls, path, focal_length_choice=FocalLengthKind.EFFECTIVE):
# here to prevent circular import
from ..io import read_table

if isinstance(focal_length_choice, str):
focal_length_choice = FocalLengthKind[focal_length_choice.upper()]

layout = read_table(
path, "/configuration/instrument/subarray/layout", table_cls=QTable
)
Expand Down Expand Up @@ -617,7 +620,7 @@ def from_hdf(cls, path, focal_length_choice="effective"):
camera = copy(cameras[row["camera_index"]])
optics = optic_descriptions[row["optics_index"]]

if focal_length_choice == "effective":
if focal_length_choice is FocalLengthKind.EFFECTIVE:
focal_length = optics.effective_focal_length
if np.isnan(focal_length.value):
raise RuntimeError(
Expand All @@ -626,8 +629,10 @@ def from_hdf(cls, path, focal_length_choice="effective"):
" Use nominal focal length or adapt configuration"
" to include the effective focal length"
)
else:
elif focal_length_choice is FocalLengthKind.EQUIVALENT:
focal_length = optics.equivalent_focal_length
else:
raise ValueError(f"Invalid focal length choice: {focal_length_choice}")

camera.geometry.frame = CameraFrame(focal_length=focal_length)
telescope_descriptions[row["tel_id"]] = TelescopeDescription(
Expand Down
14 changes: 7 additions & 7 deletions ctapipe/instrument/tests/test_subarray.py
Original file line number Diff line number Diff line change
@@ -1,11 +1,11 @@
""" Tests for SubarrayDescriptions """
from copy import deepcopy
from astropy.coordinates.earth import EarthLocation

import numpy as np
import pytest
from astropy import units as u
from astropy.coordinates import SkyCoord
import pytest
from astropy.coordinates.earth import EarthLocation

from ctapipe.coordinates import TelescopeFrame
from ctapipe.instrument import (
Expand Down Expand Up @@ -128,7 +128,7 @@ def test_hdf(example_subarray, tmp_path):
path = tmp_path / "subarray.h5"

example_subarray.to_hdf(path)
read = SubarrayDescription.from_hdf(path, focal_length_choice="nominal")
read = SubarrayDescription.from_hdf(path, focal_length_choice="EQUIVALENT")

assert example_subarray == read

Expand All @@ -147,15 +147,15 @@ def test_hdf(example_subarray, tmp_path):
with tables.open_file(path, "r+") as hdf:
del hdf.root.configuration.instrument.subarray._v_attrs.name

no_name = SubarrayDescription.from_hdf(path, focal_length_choice="nominal")
no_name = SubarrayDescription.from_hdf(path, focal_length_choice="EQUIVALENT")
assert no_name.name == "Unknown"

# Test we can also write and read to an already opened h5file
with tables.open_file(path, "w") as h5file:
example_subarray.to_hdf(h5file)

with tables.open_file(path, "r") as h5file:
read = SubarrayDescription.from_hdf(h5file, focal_length_choice="nominal")
read = SubarrayDescription.from_hdf(h5file, focal_length_choice="EQUIVALENT")
assert read == example_subarray


Expand All @@ -173,7 +173,7 @@ def test_hdf_same_camera(tmp_path):

path = tmp_path / "subarray.h5"
array.to_hdf(path)
read = SubarrayDescription.from_hdf(path, focal_length_choice="nominal")
read = SubarrayDescription.from_hdf(path, focal_length_choice="EQUIVALENT")
assert array == read


Expand Down Expand Up @@ -203,7 +203,7 @@ def test_hdf_duplicate_string_repr(tmp_path):

path = tmp_path / "subarray.h5"
array.to_hdf(path)
read = SubarrayDescription.from_hdf(path, focal_length_choice="nominal")
read = SubarrayDescription.from_hdf(path, focal_length_choice="EQUIVALENT")
assert array == read
assert (
read.tel[1].optics.num_mirror_tiles == read.tel[2].optics.num_mirror_tiles + 1
Expand Down
5 changes: 3 additions & 2 deletions ctapipe/io/eventseeker.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
Handles seeking to a particular event in a `ctapipe.io.EventSource`
"""
from copy import deepcopy

from ctapipe.core import Component

__all__ = ["EventSeeker"]
Expand All @@ -28,7 +29,7 @@ class EventSeeker(Component):
To obtain a particular event in a simtel file:
>>> from ctapipe.io import SimTelEventSource
>>> event_source = SimTelEventSource(input_url="dataset://gamma_test_large.simtel.gz", focal_length_choice="nominal")
>>> event_source = SimTelEventSource(input_url="dataset://gamma_test_large.simtel.gz", focal_length_choice="EQUIVALENT")
>>> seeker = EventSeeker(event_source=event_source)
>>> event = seeker.get_event_index(2)
>>> print(event.count)
Expand All @@ -37,7 +38,7 @@ class EventSeeker(Component):
To obtain a particular event in a simtel file from its event_id:
>>> from ctapipe.io import SimTelEventSource
>>> event_source = SimTelEventSource(input_url="dataset://gamma_test_large.simtel.gz", back_seekable=True, focal_length_choice="nominal")
>>> event_source = SimTelEventSource(input_url="dataset://gamma_test_large.simtel.gz", back_seekable=True, focal_length_choice="EQUIVALENT")
>>> seeker = EventSeeker(event_source=event_source)
>>> event = seeker.get_event_id(31007)
>>> print(event.count)
Expand Down
Loading

0 comments on commit 9692626

Please sign in to comment.