Skip to content

Commit

Permalink
Merge pull request #163 from cta-observatory/feature_add-image-extrac…
Browse files Browse the repository at this point in the history
…tion-status

Add image extraction status for TwoPassWindowSum
  • Loading branch information
HealthyPear committed Feb 1, 2022
2 parents 0a7495f + 9215680 commit 67c5ea6
Show file tree
Hide file tree
Showing 5 changed files with 635 additions and 33 deletions.
7 changes: 5 additions & 2 deletions docs/pipeline/index.rst
Original file line number Diff line number Diff line change
Expand Up @@ -57,8 +57,10 @@ Calibration
The current calibration is performed using:

* automatic gain channel selection (when more than one) above 4000 ADC counts at RO level,
* charge and pulse times extraction via ``ctapipe.image.extractors.TwoPassWindowSum``
* correction for the integration window.
* charge and pulse times extraction via ``ctapipe.image.extractors.TwoPassWindowSum``,
* correction for the integration window,
* each image is labeled with an integer status that describes if its 2nd pass extraction
was successful; the default analysis up to DL2 (CTAMARS-like) will use only such images.

.. figure:: ./double-pass-image-extraction.png
:width: 800
Expand Down Expand Up @@ -86,6 +88,7 @@ but the settings are user-dependent.

**Selection** is performed by the following requirements:

* image extraction status = 1
* at least 50 phe (still biased units),
* image's center of gravity (COG) within 80% of camera radius (radii stored in ``protopipe.pipeline.utils``),
* ellipticity between 0.1 and 0.6,
Expand Down
61 changes: 44 additions & 17 deletions protopipe/pipeline/event_preparer.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,6 @@
# CTAPIPE utilities
from ctapipe.instrument import CameraGeometry
from ctapipe.containers import ReconstructedShowerContainer
from ctapipe.calib import CameraCalibrator
from ctapipe.image.extractor import TwoPassWindowSum
from ctapipe.image import (leakage_parameters,
number_of_islands,
largest_island,
Expand All @@ -34,8 +32,10 @@
from .image_cleaning import ImageCleaner
from .utils import bcolors, effective_focal_lengths, camera_radius, CTAMARS_radii
from .temp import (
MyCameraCalibrator,
TwoPassWindowSum,
HillasParametersTelescopeFrameContainer,
HillasReconstructor,
HillasReconstructor
)

# PiWy utilities
Expand Down Expand Up @@ -70,6 +70,7 @@
"impact_dict",
"good_event",
"good_for_reco",
"image_extraction_status"
],
)

Expand All @@ -85,7 +86,8 @@ def stub(
hillas_dict_reco,
n_tels,
leakage_dict,
concentration_dict
concentration_dict,
passed
):
"""Default container for images that did not survive cleaning."""
return PreparedEvent(
Expand All @@ -108,6 +110,7 @@ def stub(
hillas_dict_reco, np.nan * u.m
), # undefined impact parameter
good_event=False,
image_extraction_status=passed
)


Expand Down Expand Up @@ -216,6 +219,7 @@ def __init__(
OrderedDict(
[
("noCuts", None),
("bad image extraction", lambda p: p == 0),
("min pixel", lambda s: np.count_nonzero(s) < npix_bounds[0]),
("min charge", lambda x: x < charge_bounds[0]),
(
Expand Down Expand Up @@ -244,7 +248,7 @@ def __init__(
# specific to TwoPassWindowSum later on
self.extractorName = list(extractor.get_current_config().items())[0][0]

self.calib = CameraCalibrator(
self.calib = MyCameraCalibrator(
image_extractor=extractor, subarray=subarray,
)

Expand Down Expand Up @@ -291,7 +295,7 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
Dictionary of of MyCameraGeometry objects for each camera in the file
return_stub : bool
If True, yield also images from events that won't be reconstructed.
This feature is not currently available.
This is required for DL1 benchmarking.
save_images : bool
If True, save photoelectron images from reconstructed events.
debug : bool
Expand Down Expand Up @@ -404,7 +408,7 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
+ "Extracting all calibrated images..."
+ bcolors.ENDC
)
self.calib(event) # Calibrate the event
passed = self.calib(event) # Calibrate the event

# =============================================================
# BEGINNING OF LOOP OVER TELESCOPES
Expand Down Expand Up @@ -480,6 +484,28 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
tel_type = str(source.subarray.tel[tel_id])
n_tels[tel_type] += 1

# We now ASSUME that the event will be good at all levels
extracted_image_is_good = True
cleaned_image_is_good = True
good_for_reco[tel_id] = 1
# later we change to 0 at the first condition NOT satisfied

# Apply some selection over image extraction
if self.image_cutflow.cut("bad image extraction", passed[tel_id]):
if debug:
print(
bcolors.WARNING +
"WARNING : bad image extraction!" +
"WARNING : image processed and recorded, but NOT used for DL2." +
bcolors.ENDC
)
# we set immediately this image as BAD at all levels
# for now (and for simplicity) we will process it anyway
# any other failing condition will just confirm this
extracted_image_is_good = False
cleaned_image_is_good = False
good_for_reco[tel_id] = 0

# use ctapipe's functionality to get the calibrated image
# and scale the reconstructed values if required
pmt_signal = event.dl1.tel[tel_id].image
Expand All @@ -490,10 +516,6 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
dl1_phe_image[tel_id] = pmt_signal
mc_phe_image[tel_id] = event.simulation.tel[tel_id].true_image

# We now ASSUME that the event will be good
good_for_reco[tel_id] = 1
# later we change to 0 if any condition is NOT satisfied

if self.cleaner_reco.mode == "tail": # tail uses only ctapipe

# Cleaning used for direction reconstruction
Expand Down Expand Up @@ -626,8 +648,6 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
# PRELIMINARY IMAGE SELECTION
# =============================================================

cleaned_image_is_good = True # we assume this

if self.image_selection_source == "extended":
cleaned_image_to_use = image_extended
elif self.image_selection_source == "biggest":
Expand Down Expand Up @@ -845,7 +865,8 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
hillas_dict_reco,
n_tels,
leakage_dict,
concentration_dict
concentration_dict,
passed
)
continue
else:
Expand Down Expand Up @@ -914,7 +935,8 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
hillas_dict_reco,
n_tels,
leakage_dict,
concentration_dict
concentration_dict,
passed
)
continue
else:
Expand Down Expand Up @@ -1024,8 +1046,10 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
hillas_dict_reco,
n_tels,
leakage_dict,
concentration_dict
concentration_dict,
passed
)
continue
else:
continue

Expand Down Expand Up @@ -1055,8 +1079,10 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
hillas_dict_reco,
n_tels,
leakage_dict,
concentration_dict
concentration_dict,
passed
)
continue
else:
continue

Expand Down Expand Up @@ -1086,4 +1112,5 @@ def prepare_event(self, source, return_stub=True, save_images=False, debug=False
impact_dict=impact_dict_reco,
good_event=True,
good_for_reco=good_for_reco,
image_extraction_status=passed
)
Loading

0 comments on commit 67c5ea6

Please sign in to comment.