Rejecting saturated stars based on custom value
I am trying to reduce LCO/MuSCAT3 & 4 data using prose. It has 4 different cameras used simultaneously and each one has different saturation level which depends on the mode of operation. The saturation values saved in the header are somehow wrong so I figured to simply re-define the saturation value by instantiating the telescope class.
saturation = 50_000 #example value
telescope = Telescope(saturation=saturation)
ref = FITSImage(ref_filename, telescope=telescope)
However, is it recommended to instantiate the telescope class inside an aperture photometry sequence? I thought it would be simpler if somehow there is an additional saturation argument in PointSourceDetection block but I couldn't work it out. Alternatively, I created a custom block to reject stars given saturation level and overwrite the sources detected in the previous block. For example
class RejectSaturatedSources(Block):
def __init__(self, saturation_level=None, name=None):
"""Filter out saturated sources from the detected sources.
Parameters
----------
saturation_level : float, optional
The level above which a star is considered saturated.
If None, will use image.saturation_level if available.
name : str, optional
Name of the block, by default None
"""
super().__init__(name=name)
self.saturation_level = saturation_level
def run(self, image):
if len(image.sources) == 0:
return
# Determine saturation level
sat_level = self.saturation_level
if sat_level is None and hasattr(image, 'saturation_level'):
sat_level = image.saturation_level
if sat_level is None:
raise ValueError("No saturation level provided or found in image")
# Filter out sources with peak values above saturation level
non_saturated = [s for s in image.sources if s.peak < sat_level]
image.sources = Sources(non_saturated, type=image.sources.type)
and then it can be inserted after the PointSourceDetection block
calibration = Sequence(
[
blocks.Trim(),
blocks.PointSourceDetection(),
RejectSaturatedSources(saturation),
blocks.Cutouts(),
blocks.MedianEPSF(),
blocks.psf.Moffat2D(),
blocks.CentroidQuadratic(),
blocks.AperturePhotometry(),
blocks.AnnulusBackground(),
]
)
calibration.run(ref, show_progress=False)
What is the correct way to do this? On the other hand, in case I need to define a .telescope file in built-ins, should I create one .telescope file per camera?
Thanks for your question @jpdeleon. I think this way is totally valid. Indeed, you would have to define a .telescope file for each camera.