Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
# Clean image using split level cleaning
mask = tailcuts_clean(self.geoms[tel_id], pmt_signal,
picture_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][1],
boundary_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][0])
# Perform Hillas parameterisation
moments = None
try:
moments_cam = hillas_parameters(event.inst.pixel_pos[tel_id][0],
event.inst.pixel_pos[tel_id][1],
pmt_signal * mask)
moments = hillas_parameters(nom_coord.x, nom_coord.y, pmt_signal * mask)
except HillasParameterizationError as e:
print(e)
continue
# Make cut based on Hillas parameters
if self.preselect(moments, tel_id):
# Dialte around edges of image
for i in range(2):
dilate(self.geoms[tel_id], mask)
# Save everything in dicts for reconstruction later
tel_x[tel_id] = tilt_tel.x
tel_y[tel_id] = tilt_tel.y
# cleaning
boundary, picture, min_neighbors = cleaning_level[geom.camera_name]
clean = tailcuts_clean(
geom,
image,
boundary_thresh=boundary,
picture_thresh=picture,
min_number_picture_neighbors=min_neighbors,
)
# ignore images with less than 5 pixels after cleaning
if clean.sum() < 5:
continue
# image parameters
hillas_c = hillas_parameters(geom[clean], image[clean])
leakage_c = leakage(geom, image, clean)
n_islands, island_ids = number_of_islands(geom, clean)
timing_c = timing_parameters(
geom[clean], image[clean], peakpos[clean], hillas_c,
)
# store parameters for stereo reconstruction
hillas_containers[telescope_id] = hillas_c
# store timegradients for plotting
# ASTRI has no timing in PROD3b, so we use skewness instead
if geom.camera_name != "ASTRICam":
time_gradients[telescope_id] = timing_c.slope.value
else:
time_gradients[telescope_id] = hillas_c.skewness
self.imclean = True
if self._counter == self.cleanframes * 2:
plt.suptitle("Image Cleaning OFF")
self.imclean = False
self._counter = 0
disp.clear_overlays()
if self.imclean:
cleanmask = tailcuts_clean(
geom, image, picture_thresh=10.0, boundary_thresh=5.0
)
for ii in range(2):
dilate(geom, cleanmask)
image[cleanmask == 0] = 0 # zero noise pixels
try:
hillas = hillas_parameters(geom, image)
disp.overlay_moments(
hillas,
with_label=False,
color="red",
alpha=0.7,
linewidth=2,
linestyle="dashed",
)
except HillasParameterizationError:
disp.clear_overlays()
pass
self.log.debug(
"Frame=%d image_sum=%.3f max=%.3f",
self._counter,
image.sum(),
length=(0.15 + 0.05 * ii) * u.m,
psi=ii * 20 * u.deg,
)
image, _, _ = model.generate_image(geom, intensity=1500, nsb_level_pe=5,)
mask = tailcuts_clean(
geom,
image,
picture_thresh=6 * image.mean(),
boundary_thresh=4 * image.mean(),
)
cleaned = image.copy()
cleaned[~mask] = 0
hillas = hillas_parameters(geom, cleaned)
disp.image = image
disp.add_colorbar(ax=axs[ii])
disp.set_limits_percent(95)
disp.overlay_moments(hillas, linewidth=3, color="blue")
# Clean image using split level cleaning
mask = tailcuts_clean(self.geoms[tel_id], pmt_signal,
picture_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][1],
boundary_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][0])
# Perform Hillas parameterisation
moments = None
try:
moments_cam = hillas_parameters(event.inst.pixel_pos[tel_id][0],
event.inst.pixel_pos[tel_id][1],
pmt_signal*mask)
moments = hillas_parameters(nom_coord.x, nom_coord.y,pmt_signal*mask)
except HillasParameterizationError as e:
print(e)
continue
# Make cut based on Hillas parameters
if self.preselect(moments, np.sum(mask), tel_id):
# Dialte around edges of image
for i in range(5):
mask = dilate(self.geoms[tel_id], mask)
# Save everything in dicts for reconstruction later
pixel_area[tel_id] = self.geoms[tel_id].pix_area/(fl*fl)
pixel_area[tel_id] *= u.rad*u.rad
pixel_x[tel_id] = nom_coord.x[mask]
rotation=-1 * self.geoms[tel_id].cam_rotation)
nom_coord = camera_coord.transform_to(nom_system)
tx, ty, tz = event.inst.tel_pos[tel_id]
# Clean image using split level cleaning
mask = tailcuts_clean(self.geoms[tel_id], pmt_signal,
picture_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][1],
boundary_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][0])
grd_tel = GroundFrame(x=tx, y=ty, z=tz)
tilt_tel = grd_tel.transform_to(tilted_system)
# Perform Hillas parameterisation
try:
moments = hillas_parameters(nom_coord.x, nom_coord.y, pmt_signal * mask)
# ImPACT reconstruction is performed in the tilted system,
# so we need to transform tel positions
except HillasParameterizationError as e:
print(e)
continue
# Make cut based on Hillas parameters
if self.preselect(moments, tel_id):
mc_tilt = mc_ground.transform_to(tilted_system)
impact_dist = np.sqrt(np.power(tilt_tel.x - mc_tilt.x, 2) +
np.power(tilt_tel.y - mc_tilt.y, 2))
self.output.add_row((event.dl0.event_id, self.geoms[tel_id].cam_id.upper(),
# ImPACT reconstruction is performed in the tilted system,
# so we need to transform tel positions
grd_tel = GroundFrame(x=tx, y=ty, z=tz)
tilt_tel = grd_tel.transform_to(tilted_system)
# Clean image using split level cleaning
mask = tailcuts_clean(self.geoms[tel_id], pmt_signal,
picture_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][1],
boundary_thresh=self.tail_cut[self.geoms[
tel_id].cam_id][0])
# Perform Hillas parameterisation
moments = None
try:
moments_cam = hillas_parameters(event.inst.pixel_pos[tel_id][0],
event.inst.pixel_pos[tel_id][1],
pmt_signal * mask)
moments = hillas_parameters(nom_coord.x, nom_coord.y, pmt_signal * mask)
except HillasParameterizationError as e:
print(e)
continue
# Make cut based on Hillas parameters
if self.preselect(moments, tel_id):
# Dialte around edges of image
for i in range(2):
dilate(self.geoms[tel_id], mask)
camgeom = subarray.tel[tel_id].camera.geometry
# note the [0] is for channel 0 which is high-gain channel
image = event.dl1.tel[tel_id].image
# Cleaning of the image
cleaned_image = image
# create a clean mask of pixels above the threshold
cleanmask = tailcuts_clean(camgeom, image, picture_thresh=10, boundary_thresh=5)
# set all rejected pixels to zero
cleaned_image[~cleanmask] = 0
# Calculate hillas parameters
# It fails for empty pixels
try:
params = hillas_parameters(camgeom, cleaned_image)
except:
continue
if params.width > 0:
hillas_params[tel_id] = params
array_pointing = SkyCoord(
az=event.mcheader.run_array_direction[0],
alt=event.mcheader.run_array_direction[1],
frame=horizon_frame,
)
if len(hillas_params) < 2:
continue
reco_result = reco.predict(