Coverage for python/lsst/meas/extensions/scarlet/metrics.py: 36%

1# This file is part of meas_extensions_scarlet. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

8# 

9# This program is free software: you can redistribute it and/or modify 

10# it under the terms of the GNU General Public License as published by 

11# the Free Software Foundation, either version 3 of the License, or 

12# (at your option) any later version. 

13# 

14# This program is distributed in the hope that it will be useful, 

15# but WITHOUT ANY WARRANTY; without even the implied warranty of 

16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 

17# GNU General Public License for more details. 

18# 

19# You should have received a copy of the GNU General Public License 

20# along with this program. If not, see <https://www.gnu.org/licenses/>. 

21 

22__all__ = [ 

23 "DeblenderMetrics", 

24 "setDeblenderMetrics", 

25] 

26 

27from dataclasses import dataclass 

28 

29import numpy as np 

30from lsst.scarlet.lite import Blend 

31 

32from .utils import nonzeroBandSupport 

33 

34 

35@dataclass 

36class DeblenderMetrics: 

37 """Metrics and measurements made on single sources. 

38 

39 Store deblender metrics to be added as attributes to a scarlet source 

40 before it is converted into a `SourceRecord`. 

41 TODO: When DM-34414 is finished this class will be eliminated and the 

42 metrics will be added to the schema using a pipeline task that calculates 

43 them from the stored deconvolved models. 

44 

45 All of the parameters are one dimensional numpy arrays, 

46 with an element for each band in the observed images. 

47 

48 `maxOverlap` is useful as a metric for determining how blended a source 

49 is because if it only overlaps with other sources at or below 

50 the noise level, it is likely to be a mostly isolated source 

51 in the deconvolved model frame. 

52 

53 `fluxOverlapFraction` is potentially more useful than the canonical 

54 "blendedness" (or purity) metric because it accounts for potential 

55 biases created during deblending by not weighting the overlapping 

56 flux with the flux of this sources model. 

57 

58 Attributes 

59 ---------- 

60 maxOverlap: 

61 The maximum overlap that the source has with its neighbors in 

62 a single pixel. 

63 fluxOverlap: 

64 The total flux from neighbors overlapping with the current source. 

65 fluxOverlapFraction: 

66 The fraction of `flux from neighbors/source flux` for a 

67 given source within the source's footprint. 

68 blendedness: 

69 The metric for determining how blended a source is using the 

70 Bosch et al. 2018 metric for "blendedness." Note that some 

71 surveys use the term "purity," which is `1-blendedness`. 

72 """ 

73 

74 maxOverlap: np.ndarray 

75 fluxOverlap: np.ndarray 

76 fluxOverlapFraction: np.ndarray 

77 blendedness: np.ndarray 

78 

79 

80def setDeblenderMetrics(blend: Blend): 

81 """Set metrics that can be used to evalute the deblender accuracy 

82 

83 This function calculates the `DeblenderMetrics` for each source in the 

84 blend, and assigns it to that sources `metrics` property in place. 

85 

86 Parameters 

87 ---------- 

88 blend: 

89 The blend containing the sources to measure. 

90 """ 

91 # Store the full model of the scene for comparison 

92 blendModelImage = blend.get_model() 

93 for src in blend.sources: 

94 # Operate over the source's own bbox (intersected with the 

95 # blend's bbox so sources that grew beyond the blend during 

96 # fitting are clipped). The metrics are sums, maxes, and 

97 # pointwise products masked by the source's support; pixels 

98 # outside the source's bbox would contribute zero to every 

99 # one of them, so restricting the math here avoids the 

100 # full-blend-sized per-source allocation done by the previous 

101 # ``project(bbox=blend.bbox)`` call. 

102 srcModelImage = src.get_model() 

103 bbox = srcModelImage.bbox & blendModelImage.bbox 

104 model = srcModelImage[:, bbox].data 

105 blendModel = blendModelImage[:, bbox].data 

106 # The footprint is the 2D array of non-zero pixels in each band 

107 footprint = nonzeroBandSupport(model) 

108 # Calculate the metrics. 

109 # See `DeblenderMetrics` for a description of each metric. 

110 neighborOverlap = (blendModel - model) * footprint[None, :, :] 

111 maxOverlap = np.max(neighborOverlap, axis=(1, 2)) 

112 fluxOverlap = np.sum(neighborOverlap, axis=(1, 2)) 

113 fluxModel = np.sum(model, axis=(1, 2)) 

114 fluxOverlapFraction = np.zeros((len(model),), dtype=float) 

115 isFinite = fluxModel > 0 

116 fluxOverlapFraction[isFinite] = fluxOverlap[isFinite] / fluxModel[isFinite] 

117 blendedness = np.zeros((len(model),), dtype=float) 

118 sqModel = np.sum(model * model, axis=(1, 2)) 

119 crossModel = np.sum(blendModel * model, axis=(1, 2)) 

120 blendedness[isFinite] = 1 - sqModel[isFinite] / crossModel[isFinite] 

121 src.metrics = DeblenderMetrics( 

122 maxOverlap, fluxOverlap, fluxOverlapFraction, blendedness 

123 )