Skip to content

Rejecting saturated stars based on custom value #196

Description

@jpdeleon

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?

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions