Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Fix arrays not being masked correctly in flatfield calculator #1149

Merged
merged 6 commits into from
Aug 29, 2023
Merged

Fix arrays not being masked correctly in flatfield calculator #1149

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
ba673a0
Fix arrays not being masked correctly in flatfield calculator
maxnoe Aug 24, 2023
7b8884f
Simplify code a bit
maxnoe Aug 24, 2023
2bdfe4b
Update mask treatment of pixels with NaN median values
FrancaCassol Aug 25, 2023
0933e9b
make all masked pixels as outliers
FrancaCassol Aug 28, 2023
8c23513
Set to zero time corrections of unusable pixels
FrancaCassol Aug 29, 2023
d91e876
Set to Nan values of masked pixels
FrancaCassol Aug 29, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 5 additions & 1 deletion lstchain/calib/camera/calibration_calculator.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -172,7 +172,7 @@ def calculate_calibration_coefficients(self, event):

# find unusable pixel from pedestal and flat-field data
unusable_pixels = np.logical_or(status_data.pedestal_failing_pixels,
status_data.flatfield_failing_pixels)
status_data.flatfield_failing_pixels)

signal = ff_data.charge_median - ped_data.charge_median

Expand Down Expand Up @@ -237,6 +237,10 @@ def calculate_calibration_coefficients(self, event):
fill_array = np.ones((constants.N_GAINS, constants.N_PIXELS)) * median_pedestal_per_sample
calib_data.pedestal_per_sample = np.ma.filled(pedestal_per_sample_masked, fill_array)

# set to zero time corrections of unusable pixels
time_correction_masked = np.ma.array(calib_data.time_correction, mask=calib_data.unusable_pixels)
calib_data.time_correction = time_correction_masked.filled(0)

# in the case FF intensity is not sufficiently high, better to scale low gain calibration from high gain results
if self.use_scaled_low_gain:
calib_data.unusable_pixels[constants.LOW_GAIN] = calib_data.unusable_pixels[constants.HIGH_GAIN]
Expand Down
44 changes: 22 additions & 22 deletions lstchain/calib/camera/flatfield.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -300,11 +300,11 @@ def calculate_time_results(
'sample_time': (time_start + (trigger_time - time_start) / 2).unix * u.s,
'sample_time_min': time_start.unix*u.s,
'sample_time_max': trigger_time.unix*u.s,
'time_mean': np.ma.getdata(pixel_mean)*u.ns,
'time_median': np.ma.getdata(pixel_median)*u.ns,
'time_std': np.ma.getdata(pixel_std)*u.ns,
'relative_time_median': np.ma.getdata(relative_median)*u.ns,
'time_median_outliers': np.ma.getdata(time_median_outliers),
'time_mean': pixel_mean.filled(np.nan)*u.ns,
'time_median': pixel_median.filled(np.nan)*u.ns,
'time_std': pixel_std.filled(np.nan)*u.ns,
'relative_time_median': relative_median.filled(np.nan)*u.ns,
'time_median_outliers': time_median_outliers.filled(True),

}

Expand All @@ -330,13 +330,19 @@ def calculate_relative_gain_results(
axis=0,
)

# mask pixels without defined statistical values (mainly due to hardware problems)
pixel_mean = np.ma.array(pixel_mean, mask=np.isnan(pixel_mean))
pixel_median = np.ma.array(pixel_median, mask=np.isnan(pixel_median))
pixel_std = np.ma.array(pixel_std, mask=np.isnan(pixel_std))

unused_values = np.abs(masked_trace_integral - pixel_mean) > (max_sigma * pixel_std)

# only warn for values discard in the sigma clipping, not those from before
outliers = unused_values & (~masked_trace_integral.mask)
check_outlier_mask(outliers, self.log, "flatfield")

# ignore outliers identified by sigma clipping also for following operations
masked_trace_integral.mask = unused_values
# add outliers identified by sigma clipping for following operations
masked_trace_integral.mask |= unused_values

# median of the median over the camera
median_of_pixel_median = np.ma.median(pixel_median, axis=1)
Expand All @@ -356,26 +362,20 @@ def calculate_relative_gain_results(
charge_median_outliers = (
np.logical_or(charge_deviation < self.charge_median_cut_outliers[0] * median_of_pixel_median[:,np.newaxis],
charge_deviation > self.charge_median_cut_outliers[1] * median_of_pixel_median[:,np.newaxis]))



# outliers from standard deviation
deviation = pixel_std - median_of_pixel_std[:, np.newaxis]
charge_std_outliers = (
np.logical_or(deviation < self.charge_std_cut_outliers[0] * std_of_pixel_std[:, np.newaxis],
deviation > self.charge_std_cut_outliers[1] * std_of_pixel_std[:, np.newaxis]))

# mask pixels with NaN mean, due to missing statistics
pixels_without_stat = np.where(np.isnan(pixel_mean)==True)
charge_median_outliers[pixels_without_stat] = True
charge_std_outliers[pixels_without_stat] = True

return {
'relative_gain_median': np.ma.getdata(np.ma.median(relative_gain_event, axis=0)),
'relative_gain_mean': np.ma.getdata(np.ma.mean(relative_gain_event, axis=0)),
'relative_gain_std': np.ma.getdata(np.ma.std(relative_gain_event, axis=0)),
'charge_median': np.ma.getdata(pixel_median),
'charge_mean': np.ma.getdata(pixel_mean),
'charge_std': np.ma.getdata(pixel_std),
'charge_std_outliers': np.ma.getdata(charge_std_outliers),
'charge_median_outliers': np.ma.getdata(charge_median_outliers),
'relative_gain_median': np.ma.median(relative_gain_event, axis=0).filled(np.nan),
'relative_gain_mean': np.ma.mean(relative_gain_event, axis=0).filled(np.nan),
'relative_gain_std': np.ma.std(relative_gain_event, axis=0).filled(np.nan),
'charge_median': pixel_median.filled(np.nan),
'charge_mean': pixel_mean.filled(np.nan),
'charge_std': pixel_std.filled(np.nan),
'charge_std_outliers': charge_std_outliers.filled(True),
'charge_median_outliers': charge_median_outliers.filled(True)
}
27 changes: 14 additions & 13 deletions lstchain/calib/camera/pedestals.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -268,13 +268,19 @@ def calculate_pedestal_results(self, trace_integral, masked_pixels_of_sample):
axis=0,
)

# mask pixels without defined statistical values (mainly due to hardware problems)
pixel_mean = np.ma.array(pixel_mean, mask=np.isnan(pixel_mean))
pixel_median = np.ma.array(pixel_median, mask=np.isnan(pixel_median))
pixel_std = np.ma.array(pixel_std, mask=np.isnan(pixel_std))

unused_values = np.abs(masked_trace_integral - pixel_mean) > (max_sigma * pixel_std)

# only warn for values discard in the sigma clipping, not those from before
outliers = unused_values & (~masked_trace_integral.mask)
check_outlier_mask(outliers, self.log, "pedestal")

# ignore outliers identified by sigma clipping also for following operations
masked_trace_integral.mask = unused_values
# add outliers identified by sigma clipping for following operations
masked_trace_integral.mask |= unused_values

# median over the camera
median_of_pixel_median = np.ma.median(pixel_median, axis=1)
Expand All @@ -301,18 +307,13 @@ def calculate_pedestal_results(self, trace_integral, masked_pixels_of_sample):
deviation < self.charge_median_cut_outliers[0] * std_of_pixel_median[:,np.newaxis],
deviation > self.charge_median_cut_outliers[1] * std_of_pixel_median[:,np.newaxis],
)

# mask pixels with NaN mean, due to missing statistics
pixels_without_stat = np.where(np.isnan(pixel_mean)==True)
charge_median_outliers[pixels_without_stat] = True
charge_std_outliers[pixels_without_stat] = True


return {
'charge_median': np.ma.getdata(pixel_median),
'charge_mean': np.ma.getdata(pixel_mean),
'charge_std': np.ma.getdata(pixel_std),
'charge_std_outliers': np.ma.getdata(charge_std_outliers),
'charge_median_outliers': np.ma.getdata(charge_median_outliers)
'charge_median': pixel_median.filled(np.nan),
'charge_mean': pixel_mean.filled(np.nan),
'charge_std': pixel_std.filled(np.nan),
'charge_std_outliers': charge_std_outliers.filled(True),
'charge_median_outliers': charge_median_outliers.filled(True)
}


Expand Down
3 changes: 2 additions & 1 deletion lstchain/visualization/plot_calib.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
import numpy as np
from astropy import units as u
from ctapipe.coordinates import EngineeringCameraFrame
from ctapipe.visualization import CameraDisplay
from ctapipe_io_lst import load_camera_geometry
Expand Down Expand Up @@ -139,7 +140,7 @@ def plot_calibration_results(ped_data, ff_data, calib_data, run=0, plot_file=Non
plt.tight_layout()
disp = CameraDisplay(camera)
disp.highlight_pixels(mask[chan], linewidth=2)
disp.image = image[chan]
disp.image = image[chan].to_value(u.ns)
disp.cmap = plt.cm.coolwarm
# disp.axes.text(lposx, 0, f'{channel[chan]} time', rotation=90)
plt.title(f"{channel[chan]} time")
Expand Down
Loading