40 badMaskPlanes=("NO_DATA", "INTRP", "BAD", "SAT", "EDGE"), detectionThreshold=5):
41 """Make detection mask and set the mask plane.
42
43 Creat a binary image from a masked image by setting all data with signal-to-
44 noise below some threshold to zero, and all data above the threshold to one.
45 If the binning parameter has been set, this procedure will be preceded by a
46 weighted binning of the data in order to smooth the result, after which the
47 result is scaled back to the original dimensions. Set the detection mask
48 plane with this binary image.
49
50 Parameters
51 ----------
52 maskedImage : `lsst.afw.image.maskedImage`
53 Image to be (optionally) binned and converted.
54 forceSlowBin : `bool`, optional
55 Force usage of slower binning method to check that the two methods
56 give the same result.
57 binning : `int`, optional
58 Number of pixels by which to bin image.
59 detectedPlane : `str`, optional
60 Name of mask with pixels that were detected above threshold in image.
61 badMaskPlanes : `set`, optional
62 Names of masks with pixels that are rejected.
63 detectionThreshold : `float`, optional
64 Boundary in signal-to-noise between non-detections and detections for
65 making a binary image from the original input image.
66 """
67 data = maskedImage.image.array
68 weights = 1 / maskedImage.variance.array
69 mask = maskedImage.getMask()
70
71 detectionMask = ((mask.array & mask.getPlaneBitMask(detectedPlane)))
72 badPixelMask = mask.getPlaneBitMask(badMaskPlanes)
73 badMask = (mask.array & badPixelMask) > 0
74 fitMask = detectionMask.astype(bool) & ~badMask
75
76 fitData = np.copy(data)
77 fitData[~fitMask] = 0
78 fitWeights = np.copy(weights)
79 fitWeights[~fitMask] = 0
80
81 if binning:
82
83 ymax, xmax = fitData.shape
84 if (ymax % binning == 0) and (xmax % binning == 0) and (not forceSlowBin):
85
86 binNumeratorReshape = (fitData * fitWeights).reshape(ymax // binning, binning,
87 xmax // binning, binning)
88 binDenominatorReshape = fitWeights.reshape(binNumeratorReshape.shape)
89 binnedNumerator = binNumeratorReshape.sum(axis=3).sum(axis=1)
90 binnedDenominator = binDenominatorReshape.sum(axis=3).sum(axis=1)
91 else:
92
93 warnings.warn('Using slow binning method--consider choosing a binsize that evenly divides '
94 f'into the image size, so that {ymax} mod binning == 0 '
95 f'and {xmax} mod binning == 0', stacklevel=2)
96 xarray = np.arange(xmax)
97 yarray = np.arange(ymax)
98 xmesh, ymesh = np.meshgrid(xarray, yarray)
99 xbins = np.arange(0, xmax + binning, binning)
100 ybins = np.arange(0, ymax + binning, binning)
101 numerator = fitWeights * fitData
102 binnedNumerator, *_ = scipy.stats.binned_statistic_2d(ymesh.ravel(), xmesh.ravel(),
103 numerator.ravel(), statistic='sum',
104 bins=(ybins, xbins))
105 binnedDenominator, *_ = scipy.stats.binned_statistic_2d(ymesh.ravel(), xmesh.ravel(),
106 fitWeights.ravel(), statistic='sum',
107 bins=(ybins, xbins))
108 binnedData = np.zeros(binnedNumerator.shape)
109 ind = binnedDenominator != 0
110 np.divide(binnedNumerator, binnedDenominator, out=binnedData, where=ind)
111 binnedWeight = binnedDenominator
112 binMask = (binnedData * binnedWeight**0.5) > detectionThreshold
113 tmpOutputMask = binMask.repeat(binning, axis=0)[:ymax]
114 outputMask = tmpOutputMask.repeat(binning, axis=1)[:, :xmax]
115 else:
116 outputMask = (fitData * fitWeights**0.5) > detectionThreshold
117
118
119 maskedImage.mask.array &= ~maskedImage.mask.getPlaneBitMask(detectedPlane)
120
121
122 maskedImage.mask.array[outputMask] |= maskedImage.mask.getPlaneBitMask(detectedPlane)
123
124
125@dataclass
