LSST Applications g013ef56533+d2224463a4,g199a45376c+0ba108daf9,g19c4beb06c+9f335b2115,g1fd858c14a+2459ca3e43,g210f2d0738+2d3d333a78,g262e1987ae+abbb004f04,g2825c19fe3+eedc38578d,g29ae962dfc+0cb55f06ef,g2cef7863aa+aef1011c0b,g35bb328faa+8c5ae1fdc5,g3fd5ace14f+19c3a54948,g47891489e3+501a489530,g4cdb532a89+a047e97985,g511e8cfd20+ce1f47b6d6,g53246c7159+8c5ae1fdc5,g54cd7ddccb+890c8e1e5d,g5fd55ab2c7+951cc3f256,g64539dfbff+2d3d333a78,g67b6fd64d1+501a489530,g67fd3c3899+2d3d333a78,g74acd417e5+0ea5dee12c,g786e29fd12+668abc6043,g87389fa792+8856018cbb,g89139ef638+501a489530,g8d7436a09f+5ea4c44d25,g8ea07a8fe4+81eaaadc04,g90f42f885a+34c0557caf,g9486f8a5af+165c016931,g97be763408+d5e351dcc8,gbf99507273+8c5ae1fdc5,gc2a301910b+2d3d333a78,gca7fc764a6+501a489530,gce8aa8abaa+8c5ae1fdc5,gd7ef33dd92+501a489530,gdab6d2f7ff+0ea5dee12c,ge410e46f29+501a489530,geaed405ab2+e3b4b2a692,gf9a733ac38+8c5ae1fdc5,w.2025.41
LSST Data Management Base Package
Loading...
Searching...
No Matches
lsst.pipe.tasks.computeExposureSummaryStats Namespace Reference

Classes

class  ComputeExposureSummaryStatsConfig
 
class  ComputeExposureSummaryStatsTask
 

Functions

 maximum_nearest_psf_distance (image_mask, psf_cat, sampling=8, bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"])
 
 compute_psf_image_deltas (image_mask, image_psf, sampling=96, ap_radius_pix=3.0, bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"])
 
 compute_ap_corr_sigma_scaled_delta (image_mask, image_ap_corr_field, psfSigma, sampling=96, bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"])
 
 compute_magnitude_limit (psfArea, skyBg, zeroPoint, readNoise, gain, snr)
 

Function Documentation

◆ compute_ap_corr_sigma_scaled_delta()

lsst.pipe.tasks.computeExposureSummaryStats.compute_ap_corr_sigma_scaled_delta ( image_mask,
image_ap_corr_field,
psfSigma,
sampling = 96,
bad_mask_bits = ["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"] )
Compute the delta between the maximum and minimum aperture correction
values scaled (divided) by ``psfSigma`` for the given field representation,
``image_ap_corr_field`` evaluated on a grid of points lying in the
unmasked region of the image.

Parameters
----------
image_mask : `lsst.afw.image.Mask`
    The mask plane associated with the exposure.
image_ap_corr_field : `lsst.afw.math.ChebyshevBoundedField`
    The ChebyshevBoundedField representation of the aperture correction
    of interest for the exposure.
psfSigma : `float`
    The PSF model second-moments determinant radius (center of chip)
    in pixels.
sampling : `int`, optional
    Sampling rate in each dimension to create the grid of points at which
    to evaluate ``image_psf``s trace radius value. The tradeoff is between
    adequate sampling versus speed.
bad_mask_bits : `list` [`str`], optional
    Mask bits required to be absent for a pixel to be considered
    "unmasked".

Returns
-------
ap_corr_sigma_scaled_delta : `float`
    The delta between the maximum and minimum of the (multiplicative)
    aperture correction values scaled (divided) by ``psfSigma`` evaluated
    on the x,y-grid subsampled on the unmasked detector pixels by a factor
    of ``sampling``.  If the aperture correction evaluates to NaN on any
    of the grid points, this is set to NaN.

Definition at line 843 of file computeExposureSummaryStats.py.

849):
850 """Compute the delta between the maximum and minimum aperture correction
851 values scaled (divided) by ``psfSigma`` for the given field representation,
852 ``image_ap_corr_field`` evaluated on a grid of points lying in the
853 unmasked region of the image.
854
855 Parameters
856 ----------
857 image_mask : `lsst.afw.image.Mask`
858 The mask plane associated with the exposure.
859 image_ap_corr_field : `lsst.afw.math.ChebyshevBoundedField`
860 The ChebyshevBoundedField representation of the aperture correction
861 of interest for the exposure.
862 psfSigma : `float`
863 The PSF model second-moments determinant radius (center of chip)
864 in pixels.
865 sampling : `int`, optional
866 Sampling rate in each dimension to create the grid of points at which
867 to evaluate ``image_psf``s trace radius value. The tradeoff is between
868 adequate sampling versus speed.
869 bad_mask_bits : `list` [`str`], optional
870 Mask bits required to be absent for a pixel to be considered
871 "unmasked".
872
873 Returns
874 -------
875 ap_corr_sigma_scaled_delta : `float`
876 The delta between the maximum and minimum of the (multiplicative)
877 aperture correction values scaled (divided) by ``psfSigma`` evaluated
878 on the x,y-grid subsampled on the unmasked detector pixels by a factor
879 of ``sampling``. If the aperture correction evaluates to NaN on any
880 of the grid points, this is set to NaN.
881 """
882 mask_arr = image_mask.array[::sampling, ::sampling]
883 bitmask = image_mask.getPlaneBitMask(bad_mask_bits)
884 good = ((mask_arr & bitmask) == 0)
885
886 x = np.arange(good.shape[1], dtype=np.float64) * sampling
887 y = np.arange(good.shape[0], dtype=np.float64) * sampling
888 xx, yy = np.meshgrid(x, y)
889
890 ap_corr = image_ap_corr_field.evaluate(xx.ravel(), yy.ravel()).reshape(xx.shape)
891 ap_corr_good = ap_corr[good]
892 if ~np.isfinite(ap_corr_good).all():
893 ap_corr_sigma_scaled_delta = float("nan")
894 else:
895 ap_corr_sigma_scaled_delta = (np.max(ap_corr_good) - np.min(ap_corr_good))/psfSigma
896
897 return ap_corr_sigma_scaled_delta
898
899

◆ compute_magnitude_limit()

lsst.pipe.tasks.computeExposureSummaryStats.compute_magnitude_limit ( psfArea,
skyBg,
zeroPoint,
readNoise,
gain,
snr )
Compute the expected point-source magnitude limit at a given
signal-to-noise ratio given the exposure-level metadata. Based on
the signal-to-noise formula provided in SMTN-002 (see LSE-40 for
more details on the calculation).

  SNR = C / sqrt( C/g + (B/g + sigma_inst**2) * neff )

where C is the counts from the source, B is counts from the (sky)
background, sigma_inst is the instrumental (read) noise, neff is
the effective size of the PSF, and g is the gain in e-/ADU.  Note
that input values of ``skyBg``, ``zeroPoint``, and ``readNoise``
should all consistently be in electrons or ADU.

Parameters
----------
psfArea : `float`
    The effective area of the PSF [pix].
skyBg : `float`
    The sky background counts for the exposure [ADU or e-].
zeroPoint : `float`
    The zeropoint (includes exposure time) [ADU or e-].
readNoise : `float`
    The instrumental read noise for the exposure [ADU or e-].
gain : `float`
    The instrumental gain for the exposure [e-/ADU]. The gain should
     be 1.0 if the skyBg, zeroPoint, and readNoise are in e-.
snr : `float`
    Signal-to-noise ratio at which magnitude limit is calculated.

Returns
-------
magnitude_limit : `float`
    The expected magnitude limit at the given signal to noise.

Definition at line 900 of file computeExposureSummaryStats.py.

907):
908 """Compute the expected point-source magnitude limit at a given
909 signal-to-noise ratio given the exposure-level metadata. Based on
910 the signal-to-noise formula provided in SMTN-002 (see LSE-40 for
911 more details on the calculation).
912
913 SNR = C / sqrt( C/g + (B/g + sigma_inst**2) * neff )
914
915 where C is the counts from the source, B is counts from the (sky)
916 background, sigma_inst is the instrumental (read) noise, neff is
917 the effective size of the PSF, and g is the gain in e-/ADU. Note
918 that input values of ``skyBg``, ``zeroPoint``, and ``readNoise``
919 should all consistently be in electrons or ADU.
920
921 Parameters
922 ----------
923 psfArea : `float`
924 The effective area of the PSF [pix].
925 skyBg : `float`
926 The sky background counts for the exposure [ADU or e-].
927 zeroPoint : `float`
928 The zeropoint (includes exposure time) [ADU or e-].
929 readNoise : `float`
930 The instrumental read noise for the exposure [ADU or e-].
931 gain : `float`
932 The instrumental gain for the exposure [e-/ADU]. The gain should
933 be 1.0 if the skyBg, zeroPoint, and readNoise are in e-.
934 snr : `float`
935 Signal-to-noise ratio at which magnitude limit is calculated.
936
937 Returns
938 -------
939 magnitude_limit : `float`
940 The expected magnitude limit at the given signal to noise.
941
942 """
943 # Effective number of pixels within the PSF calculated directly
944 # from the PSF model.
945 neff = psfArea
946
947 # Instrumental noise (read noise only)
948 sigma_inst = readNoise
949
950 # Total background counts derived from Eq. 45 of LSE-40
951 background = (skyBg/gain + sigma_inst**2) * neff
952 # Source counts to achieve the desired SNR (from quadratic formula)
953 source = (snr**2)/(2*gain) + np.sqrt((snr**4)/(4*gain) + snr**2 * background)
954
955 # Convert to a magnitude using the zeropoint.
956 # Note: Zeropoint currently includes exposure time
957 magnitude_limit = -2.5*np.log10(source) + zeroPoint
958
959 return magnitude_limit

◆ compute_psf_image_deltas()

lsst.pipe.tasks.computeExposureSummaryStats.compute_psf_image_deltas ( image_mask,
image_psf,
sampling = 96,
ap_radius_pix = 3.0,
bad_mask_bits = ["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"] )
Compute the delta between the maximum and minimum model PSF trace radius
values evaluated on a grid of points lying in the unmasked region of the
image.

Parameters
----------
image_mask : `lsst.afw.image.Mask`
    The mask plane associated with the exposure.
image_psf : `lsst.afw.detection.Psf`
    The PSF model associated with the exposure.
sampling : `int`, optional
    Sampling rate in each dimension to create the grid of points at which
    to evaluate ``image_psf``s trace radius value. The tradeoff is between
    adequate sampling versus speed.
ap_radius_pix : `float`, optional
    Radius in pixels of the aperture on which to measure the flux of the
    PSF model.
bad_mask_bits : `list` [`str`], optional
    Mask bits required to be absent for a pixel to be considered
    "unmasked".

Returns
-------
psf_trace_radius_delta, psf_ap_flux_delta : `float`
    The delta (in pixels) between the maximum and minimum model PSF trace
    radius values and the PSF aperture fluxes (with aperture radius of
    max(2, 3*psfSigma)) evaluated on the x,y-grid subsampled on the
    unmasked detector pixels by a factor of ``sampling``.  If both the
    model PSF trace radius value and aperture flux value on the grid
    evaluate to NaN, then NaNs are returned immediately.

Definition at line 775 of file computeExposureSummaryStats.py.

781):
782 """Compute the delta between the maximum and minimum model PSF trace radius
783 values evaluated on a grid of points lying in the unmasked region of the
784 image.
785
786 Parameters
787 ----------
788 image_mask : `lsst.afw.image.Mask`
789 The mask plane associated with the exposure.
790 image_psf : `lsst.afw.detection.Psf`
791 The PSF model associated with the exposure.
792 sampling : `int`, optional
793 Sampling rate in each dimension to create the grid of points at which
794 to evaluate ``image_psf``s trace radius value. The tradeoff is between
795 adequate sampling versus speed.
796 ap_radius_pix : `float`, optional
797 Radius in pixels of the aperture on which to measure the flux of the
798 PSF model.
799 bad_mask_bits : `list` [`str`], optional
800 Mask bits required to be absent for a pixel to be considered
801 "unmasked".
802
803 Returns
804 -------
805 psf_trace_radius_delta, psf_ap_flux_delta : `float`
806 The delta (in pixels) between the maximum and minimum model PSF trace
807 radius values and the PSF aperture fluxes (with aperture radius of
808 max(2, 3*psfSigma)) evaluated on the x,y-grid subsampled on the
809 unmasked detector pixels by a factor of ``sampling``. If both the
810 model PSF trace radius value and aperture flux value on the grid
811 evaluate to NaN, then NaNs are returned immediately.
812 """
813 psf_trace_radius_list = []
814 psf_ap_flux_list = []
815 mask_arr = image_mask.array[::sampling, ::sampling]
816 bitmask = image_mask.getPlaneBitMask(bad_mask_bits)
817 good = ((mask_arr & bitmask) == 0)
818
819 x = np.arange(good.shape[1]) * sampling
820 y = np.arange(good.shape[0]) * sampling
821 xx, yy = np.meshgrid(x, y)
822
823 for x_mesh, y_mesh, good_mesh in zip(xx, yy, good):
824 for x_point, y_point, is_good in zip(x_mesh, y_mesh, good_mesh):
825 if is_good:
826 position = geom.Point2D(x_point, y_point)
827 psf_trace_radius = image_psf.computeShape(position).getTraceRadius()
828 psf_ap_flux = image_psf.computeApertureFlux(ap_radius_pix, position)
829 if ~np.isfinite(psf_trace_radius) and ~np.isfinite(psf_ap_flux):
830 return float("nan"), float("nan")
831 psf_trace_radius_list.append(psf_trace_radius)
832 psf_ap_flux_list.append(psf_ap_flux)
833
834 psf_trace_radius_delta = np.max(psf_trace_radius_list) - np.min(psf_trace_radius_list)
835 if np.any(np.asarray(psf_ap_flux_list) < 0.0): # Consider any -ve flux entry as "bad".
836 psf_ap_flux_delta = float("nan")
837 else:
838 psf_ap_flux_delta = np.max(psf_ap_flux_list) - np.min(psf_ap_flux_list)
839
840 return psf_trace_radius_delta, psf_ap_flux_delta
841
842

◆ maximum_nearest_psf_distance()

lsst.pipe.tasks.computeExposureSummaryStats.maximum_nearest_psf_distance ( image_mask,
psf_cat,
sampling = 8,
bad_mask_bits = ["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"] )
Compute the maximum distance of an unmasked pixel to its nearest PSF.

Parameters
----------
image_mask : `lsst.afw.image.Mask`
    The mask plane associated with the exposure.
psf_cat : `lsst.afw.table.SourceCatalog` or `astropy.table.Table`
    Catalog containing only the stars used in the PSF modeling.
sampling : `int`
    Sampling rate in each dimension to create the grid of points on which
    to evaluate the distance to the nearest PSF star. The tradeoff is
    between adequate sampling versus speed.
bad_mask_bits : `list` [`str`]
    Mask bits required to be absent for a pixel to be considered
    "unmasked".

Returns
-------
max_dist_to_nearest_psf : `float`
    The maximum distance (in pixels) of an unmasked pixel to its nearest
    PSF model star.

Definition at line 728 of file computeExposureSummaryStats.py.

733):
734 """Compute the maximum distance of an unmasked pixel to its nearest PSF.
735
736 Parameters
737 ----------
738 image_mask : `lsst.afw.image.Mask`
739 The mask plane associated with the exposure.
740 psf_cat : `lsst.afw.table.SourceCatalog` or `astropy.table.Table`
741 Catalog containing only the stars used in the PSF modeling.
742 sampling : `int`
743 Sampling rate in each dimension to create the grid of points on which
744 to evaluate the distance to the nearest PSF star. The tradeoff is
745 between adequate sampling versus speed.
746 bad_mask_bits : `list` [`str`]
747 Mask bits required to be absent for a pixel to be considered
748 "unmasked".
749
750 Returns
751 -------
752 max_dist_to_nearest_psf : `float`
753 The maximum distance (in pixels) of an unmasked pixel to its nearest
754 PSF model star.
755 """
756 mask_arr = image_mask.array[::sampling, ::sampling]
757 bitmask = image_mask.getPlaneBitMask(bad_mask_bits)
758 good = ((mask_arr & bitmask) == 0)
759
760 x = np.arange(good.shape[1]) * sampling
761 y = np.arange(good.shape[0]) * sampling
762 xx, yy = np.meshgrid(x, y)
763
764 dist_to_nearest_psf = np.full(good.shape, np.inf)
765 for psf in psf_cat:
766 x_psf = psf["slot_Centroid_x"]
767 y_psf = psf["slot_Centroid_y"]
768 dist_to_nearest_psf = np.minimum(dist_to_nearest_psf, np.hypot(xx - x_psf, yy - y_psf))
769 unmasked_dists = dist_to_nearest_psf * good
770 max_dist_to_nearest_psf = np.max(unmasked_dists)
771
772 return max_dist_to_nearest_psf
773
774