Coverage for tests / test_generalShapeletPsfModels.py: 22%

1# 

2# LSST Data Management System 

3# 

4# Copyright 2008-2016 AURA/LSST. 

5# 

6# This product includes software developed by the 

7# LSST Project (http://www.lsst.org/). 

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 LSST License Statement and 

20# the GNU General Public License along with this program. If not, 

21# see <https://www.lsstcorp.org/LegalNotices/>. 

22# 

23import unittest 

24import numpy 

25import os 

26 

27import lsst.utils.tests 

28import lsst.shapelet 

29import lsst.geom 

30import lsst.afw.geom.ellipses 

31import lsst.afw.table 

32import lsst.afw.detection 

33import lsst.utils.logging 

34import lsst.meas.modelfit 

35import lsst.meas.base 

36import lsst.meas.algorithms 

37 

38# Set trace to 0-5 to view debug messages. Level 5 enables all traces. 

39lsst.utils.logging.trace_set_at("lsst.meas.modelfit.optimizer.Optimizer", -1) 

40lsst.utils.logging.trace_set_at("lsst.meas.modelfit.optimizer.solveTrustRegion", -1) 

41 

42 

43class GeneralShapeletPsfApproxPluginsTestCase(lsst.utils.tests.TestCase): 

44 

45 def setUp(self): 

46 numpy.random.seed(500) 

47 self.exposure = lsst.afw.image.ExposureF(41, 41) 

48 self.schema = lsst.afw.table.SourceTable.makeMinimalSchema() 

49 self.centroidKey = lsst.afw.table.Point2DKey.addFields(self.schema, "centroid", "centroid", "pixel") 

50 self.schema.getAliasMap().set("slot_Centroid", "centroid") 

51 self.psfDir = os.path.join(os.environ["MEAS_MODELFIT_DIR"], "tests", "data", "psfs") 

52 

53 def tearDown(self): 

54 del self.exposure 

55 del self.schema 

56 del self.centroidKey 

57 del self.psfDir 

58 

59 def makePsf(self, psfname, max=None): 

60 data = lsst.afw.image.ImageF(os.path.join(self.psfDir, psfname)).getArray().astype(numpy.float64) 

61 if max is not None: 

62 trim0 = (data.shape[0] - max)//2 

63 trim1 = (data.shape[1] - max)//2 

64 if trim0 > 0 and trim1 > 0: 

65 data = data[trim0:trim0+max, trim1:trim1+max] 

66 kernel = lsst.afw.math.FixedKernel(lsst.afw.image.ImageD(data)) 

67 return lsst.meas.algorithms.KernelPsf(kernel) 

68 

69 def runTask(self, psftype, sequence): 

70 config = lsst.meas.base.SingleFrameMeasurementTask.ConfigClass() 

71 config.slots.centroid = None 

72 config.slots.shape = None 

73 config.slots.psfFlux = None 

74 config.slots.apFlux = None 

75 config.slots.gaussianFlux = None 

76 config.slots.modelFlux = None 

77 config.slots.calibFlux = None 

78 config.doReplaceWithNoise = False 

79 config.plugins.names = ["modelfit_GeneralShapeletPsfApprox"] 

80 config.plugins["modelfit_GeneralShapeletPsfApprox"].sequence = sequence 

81 task = lsst.meas.base.SingleFrameMeasurementTask(config=config, schema=self.schema) 

82 measCat = lsst.afw.table.SourceCatalog(self.schema) 

83 measRecord = measCat.addNew() 

84 measRecord.set(self.centroidKey, lsst.geom.Point2D(20.0, 20.0)) 

85 measRecord.set(self.centroidKey, lsst.geom.Point2D(20.0, 20.0)) 

86 task.run(measCat, self.exposure) 

87 return measRecord 

88 

89 def testSingleGaussian(self): 

90 sigma1 = 3.0 

91 self.exposure.setPsf(lsst.afw.detection.GaussianPsf(19, 19, sigma1)) 

92 measRecord = self.runTask("Single Gaussian Psf", ["SingleGaussian"]) 

93 keySingleGaussian = lsst.shapelet.MultiShapeletFunctionKey( 

94 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["SingleGaussian"] 

95 ) 

96 msfSingleGaussian = measRecord.get(keySingleGaussian) 

97 self.assertEqual(len(msfSingleGaussian.getComponents()), 1) 

98 comps = msfSingleGaussian.getComponents() 

99 r0 = comps[0].getEllipse().getCore().getDeterminantRadius() 

100 self.assertFloatsAlmostEqual(r0, sigma1, .05) 

101 

102 def testDoubleGaussian(self): 

103 sigma1 = 2.0 

104 sigma2 = 4.0 

105 self.exposure.setPsf(lsst.meas.algorithms.DoubleGaussianPsf(19, 19, sigma1, sigma2, .25)) 

106 measRecord = self.runTask("Double Gaussian Psf", ["DoubleGaussian"]) 

107 keyDoubleGaussian = lsst.shapelet.MultiShapeletFunctionKey( 

108 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["DoubleGaussian"] 

109 ) 

110 msf = measRecord.get(keyDoubleGaussian) 

111 comps = msf.getComponents() 

112 self.assertEqual(len(comps), 2) 

113 # amplitudes are equal by construction 

114 A0 = measRecord.get("modelfit_GeneralShapeletPsfApprox_DoubleGaussian_0_0") 

115 A1 = measRecord.get("modelfit_GeneralShapeletPsfApprox_DoubleGaussian_1_0") 

116 self.assertFloatsAlmostEqual(A0, A1, .05) 

117 r0 = comps[0].getEllipse().getCore().getDeterminantRadius() 

118 r1 = comps[1].getEllipse().getCore().getDeterminantRadius() 

119 self.assertFloatsAlmostEqual(r0, sigma1, .05) 

120 self.assertFloatsAlmostEqual(r1, sigma2, .05) 

121 

122 def testDoubleShapelet(self): 

123 self.exposure.setPsf(self.makePsf("galsimPsf_0.5.fits", max=33)) 

124 measRecord = self.runTask("Galsim Psf", ["DoubleShapelet"]) 

125 keyDoubleShapelet = lsst.shapelet.MultiShapeletFunctionKey( 

126 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["DoubleShapelet"] 

127 ) 

128 msf = measRecord.get(keyDoubleShapelet) 

129 comps = msf.getComponents() 

130 self.assertEqual(len(comps), 2) 

131 A0 = measRecord.get("modelfit_GeneralShapeletPsfApprox_DoubleShapelet_0_0") 

132 A1 = measRecord.get("modelfit_GeneralShapeletPsfApprox_DoubleShapelet_1_0") 

133 self.assertGreater(A1, .04) 

134 self.assertGreater(A0, .04) 

135 

136 def testFull(self): 

137 self.exposure.setPsf(self.makePsf("galsimPsf_0.9.fits", max=33)) 

138 measRecord = self.runTask("Galsim Psf", ["Full"]) 

139 keyFull = lsst.shapelet.MultiShapeletFunctionKey( 

140 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["Full"] 

141 ) 

142 msf = measRecord.get(keyFull) 

143 comps = msf.getComponents() 

144 self.assertEqual(len(comps), 4) 

145 A1 = measRecord.get("modelfit_GeneralShapeletPsfApprox_Full_1_0") 

146 A2 = measRecord.get("modelfit_GeneralShapeletPsfApprox_Full_2_0") 

147 # test the primary and wings to be sure we are getting something 

148 self.assertGreater(A2, .04) 

149 self.assertGreater(A1, .04) 

150 

151 def testSequence(self): 

152 sigma1 = 2.0 

153 sigma2 = 4.0 

154 self.exposure.setPsf(lsst.meas.algorithms.DoubleGaussianPsf(19, 19, sigma1, sigma2, .25)) 

155 measRecord = self.runTask("Single Gaussian Psf", ["SingleGaussian", "DoubleGaussian", 

156 "DoubleShapelet"]) 

157 keySingleGaussian = lsst.shapelet.MultiShapeletFunctionKey( 

158 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["SingleGaussian"] 

159 ) 

160 msfSingleGaussian = measRecord.get(keySingleGaussian) 

161 self.assertEqual(len(msfSingleGaussian.getComponents()), 1) 

162 comps = msfSingleGaussian.getComponents() 

163 r0 = comps[0].getEllipse().getCore().getDeterminantRadius() 

164 # don't expect it to be all that close, but the DoubleGaussian should be 

165 self.assertFloatsAlmostEqual(r0, sigma1, .3) 

166 

167 keyDoubleGaussian = lsst.shapelet.MultiShapeletFunctionKey( 

168 self.schema["modelfit"]["GeneralShapeletPsfApprox"]["DoubleGaussian"] 

169 ) 

170 msfDoubleGaussian = measRecord.get(keyDoubleGaussian) 

171 comps = msfDoubleGaussian.getComponents() 

172 r0 = comps[0].getEllipse().getCore().getDeterminantRadius() 

173 r1 = comps[1].getEllipse().getCore().getDeterminantRadius() 

174 self.assertFloatsAlmostEqual(r0, sigma1, .05) 

175 self.assertFloatsAlmostEqual(r1, sigma2, .05) 

176 

177 

178class TestMemory(lsst.utils.tests.MemoryTestCase): 

179 pass 

180 

181 

182def setup_module(module): 

183 lsst.utils.tests.init() 

184 

185 

186if __name__ == "__main__": 186 ↛ 187line 186 didn't jump to line 187 because the condition on line 186 was never true

187 lsst.utils.tests.init() 

188 unittest.main()