Coverage for python / lsst / analysis / tools / tasks / coaddDepthSummary.py: 32%

1# This file is part of analysis_tools. 

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/>. 

21from __future__ import annotations 

22 

23__all__ = ( 

24 "CoaddDepthSummaryConfig", 

25 "CoaddDepthSummaryTask", 

26) 

27 

28 

29import numpy as np 

30from astropy.table import Table 

31 

32from lsst.pex.config import ListField 

33from lsst.pipe.base import PipelineTask, PipelineTaskConfig, PipelineTaskConnections, Struct 

34from lsst.pipe.base import connectionTypes as cT 

35 

36 

37class CoaddDepthSummaryConnections( 

38 PipelineTaskConnections, 

39 dimensions=("tract", "skymap"), 

40 defaultTemplates={"coaddName": ""}, # set as either deep or template in the pipeline 

41): 

42 data = cT.Input( 

43 doc="Coadd n_image to load from the butler (pixel values are the number of input images).", 

44 name="{coaddName}_coadd_n_image", 

45 storageClass="ImageU", 

46 multiple=True, 

47 dimensions=("tract", "patch", "band", "skymap"), 

48 deferLoad=True, 

49 ) 

50 

51 statTable = cT.Output( 

52 doc="Table with resulting n_image based depth statistics.", 

53 name="{coaddName}_coadd_depth_table", 

54 storageClass="ArrowAstropy", 

55 dimensions=("tract", "skymap"), 

56 ) 

57 

58 

59class CoaddDepthSummaryConfig(PipelineTaskConfig, pipelineConnections=CoaddDepthSummaryConnections): 

60 threshold_list = ListField( 

61 default=[1, 2, 3, 5, 12], 

62 dtype=int, 

63 doc="The n_image pixel value thresholds, in ascending order.", 

64 ) 

65 

66 quantile_list = ListField( 

67 default=[5, 10, 25, 50, 75, 90, 95], 

68 dtype=int, 

69 doc="The percentiles at which to compute n_image values, in ascending order.", 

70 ) 

71 

72 

73class CoaddDepthSummaryTask(PipelineTask): 

74 ConfigClass = CoaddDepthSummaryConfig 

75 _DefaultName = "coaddDepthSummary" 

76 

77 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

78 inputs = butlerQC.get(inputRefs) 

79 outputs = self.run(inputs) 

80 butlerQC.put(outputs, outputRefs) 

81 

82 def run(self, inputs): 

83 t = Table() 

84 bands = [] 

85 patches = [] 

86 medians = [] 

87 stdevs = [] 

88 stats = [] 

89 quantiles = [] 

90 

91 for n_image_handle in inputs["data"]: 

92 n_image = n_image_handle.get() 

93 data_id = n_image_handle.dataId 

94 band = str(data_id.band.name) 

95 patch = int(data_id.patch.id) 

96 median = np.nanmedian(n_image.array) 

97 stdev = np.nanstd(n_image.array) 

98 

99 bands.append(band) 

100 patches.append(patch) 

101 medians.append(median) 

102 stdevs.append(stdev) 

103 

104 band_patch_stats = [] 

105 for threshold in self.config.threshold_list: 

106 # Calculate the percentage of the image with an image depth 

107 # above the given threshold. 

108 stat = np.sum(n_image.array >= threshold) * 100 / (n_image.getHeight() * n_image.getWidth()) 

109 band_patch_stats.append(stat) 

110 

111 stats.append(band_patch_stats) 

112 

113 # Calculate the quantiles for image depth 

114 # across the whole n_image array. 

115 quantile = list(np.percentile(n_image.array, q=self.config.quantile_list)) 

116 quantiles.append(quantile) 

117 

118 threshold_col_names = [ 

119 f"depth_above_threshold_{threshold}" for threshold in self.config.threshold_list 

120 ] 

121 quantile_col_names = [f"depth_{q}_percentile" for q in self.config.quantile_list] 

122 

123 # Construct the Astropy table 

124 data = [patches, bands, medians, stdevs] + list(zip(*stats)) + list(zip(*quantiles)) 

125 names = ["patch", "band", "medians", "stdevs"] + threshold_col_names + quantile_col_names 

126 dtype = ( 

127 ["int", "str", "float", "float"] 

128 + ["float" for x in range(len(list(zip(*stats))))] 

129 + ["int" for y in range(len(list(zip(*quantiles))))] 

130 ) 

131 t = Table(data=data, names=names, dtype=dtype) 

132 return Struct(statTable=t)