Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1"""A group of metrics that work together to evaluate season characteristics (length, number, etc). 

2In addition, these supports the time delay metric calculation for strong lensing. 

3""" 

4 

5import numpy as np 

6from .baseMetric import BaseMetric 

7 

8__all__ = ['calcSeason', 'findSeasonEdges', 

9 'SeasonLengthMetric', 'CampaignLengthMetric', 

10 'MeanCampaignFrequencyMetric', 'TdcMetric'] 

11 

12 

13def calcSeason(ra, time): 

14 """Calculate the 'season' in the survey for a series of ra/dec/time values of an observation. 

15 Based only on the RA of the point on the sky, it calculates the 'season' based on when this 

16 point would be overhead .. the season is considered +/- 0.5 years around this time. 

17 

18 Parameters 

19 ---------- 

20 ra : float 

21 The RA (in degrees) of the point on the sky 

22 time : np.ndarray 

23 The times of the observations, in MJD 

24 

25 Returns 

26 ------- 

27 np.ndarray 

28 The season values 

29 """ 

30 # Reference RA and equinox to anchor ra/season reference - RA = 0 is overhead at autumnal equinox 

31 # autumn equinox 2014 happened on september 23 --> equinox MJD 

32 Equinox = 2456923.5 - 2400000.5 

33 # convert ra into 'days' 

34 dayRA = ra / 360 * 365.25 

35 firstSeasonBegan = Equinox + dayRA - 0.5 * 365.25 

36 seasons = (time - firstSeasonBegan) / 365.25 

37 # Set first season to 0 

38 seasons = seasons - np.floor(np.min(seasons)) 

39 return seasons 

40 

41def findSeasonEdges(seasons): 

42 """Given the seasons, return the indexes of each start/end of the season. 

43 

44 Parameters 

45 ---------- 

46 seasons: np.ndarray 

47 Seasons, such as calculated by calcSeason. 

48 Note that seasons should be sorted!! 

49 

50 Returns 

51 ------- 

52 np.ndarray, np.ndarray 

53 The indexes of the first and last date in the season. 

54 """ 

55 intSeasons = np.floor(seasons) 

56 # Get the unique seasons, so that we can separate each one 

57 season_list = np.unique(intSeasons) 

58 # Find the first and last observation of each season. 

59 firstOfSeason = np.searchsorted(intSeasons, season_list) 

60 lastOfSeason = np.searchsorted(intSeasons, season_list, side='right') - 1 

61 return firstOfSeason, lastOfSeason 

62 

63 

64class SeasonLengthMetric(BaseMetric): 

65 """ 

66 Calculate the length of LSST seasons, in days. 

67 

68 Parameters 

69 ---------- 

70 minExpTime: float, opt 

71 Minimum visit exposure time to count for a 'visit', in seconds. Default 20. 

72 reduceFunc : function, optional 

73 Function that can operate on array-like structures. Typically numpy function. 

74 This reduces the season length in each season from 10 separate values to a single value. 

75 Default np.median. 

76 """ 

77 def __init__(self, mjdCol='observationStartMJD', expTimeCol='visitExposureTime', minExpTime=20, 

78 reduceFunc=np.median, metricName='SeasonLength', **kwargs): 

79 units = 'days' 

80 self.mjdCol = mjdCol 

81 self.expTimeCol = expTimeCol 

82 self.minExpTime = minExpTime 

83 self.reduceFunc = reduceFunc 

84 super().__init__(col=[self.mjdCol, self.expTimeCol], 

85 units=units, metricName=metricName, **kwargs) 

86 

87 def run(self, dataSlice, slicePoint): 

88 """Calculate the (reduceFunc) of the length of each season. 

89 Uses the slicePoint RA/Dec to calculate the position in question, then uses the times of the visits 

90 to assign them into seasons (based on where the sun is relative to the slicePoint RA). 

91 

92 Parameters 

93 ---------- 

94 dataSlice : numpy.array 

95 Numpy structured array containing the data related to the visits provided by the slicer. 

96 slicePoint : dict 

97 Dictionary containing information about the slicepoint currently active in the slicer. 

98 

99 Returns 

100 ------- 

101 float 

102 The (reduceFunc) of the length of each season, in days. 

103 """ 

104 # Order data Slice/times and exclude visits which are too short. 

105 long = np.where(dataSlice[self.expTimeCol] > self.minExpTime) 

106 data = np.sort(dataSlice[long], order=self.mjdCol) 

107 # SlicePoints ra/dec are always in radians - convert to degrees to calculate season 

108 seasons = calcSeason(np.degrees(slicePoint['ra']), data[self.mjdCol]) 

109 firstOfSeason, lastOfSeason = findSeasonEdges(seasons) 

110 seasonlengths = data[self.mjdCol][lastOfSeason] - data[self.mjdCol][firstOfSeason] 

111 result = self.reduceFunc(seasonlengths) 

112 return result 

113 

114 

115class CampaignLengthMetric(BaseMetric): 

116 """Calculate the number of seasons (roughly, years) a pointing is observed for. 

117 This corresponds to the 'campaign length' for lensed quasar time delays. 

118 """ 

119 def __init__(self, mjdCol='observationStartMJD', expTimeCol='visitExposureTime', minExpTime=20, **kwargs): 

120 units = '' 

121 self.expTimeCol = expTimeCol 

122 self.minExpTime = minExpTime 

123 self.mjdCol = mjdCol 

124 super().__init__(col=[self.mjdCol, self.expTimeCol], units=units, **kwargs) 

125 

126 def run(self, dataSlice, slicePoint): 

127 # Order data Slice/times and exclude visits which are too short. 

128 long = np.where(dataSlice[self.expTimeCol] > self.minExpTime) 

129 data = np.sort(dataSlice[long], order=self.mjdCol) 

130 seasons = calcSeason(np.degrees(slicePoint['ra'], data[self.mjdCol])) 

131 intSeasons = np.floor(seasons) 

132 count = len(np.unique(intSeasons)) 

133 return count 

134 

135 

136class MeanCampaignFrequencyMetric(BaseMetric): 

137 """Calculate the mean separation between nights, within a season - then the mean over the campaign. 

138 Calculate per season, to avoid any influence from season gaps. 

139 """ 

140 def __init__(self, mjdCol='observationStartMJD', expTimeCol='visitExposureTime', minExpTime=20, 

141 nightCol='night', **kwargs): 

142 self.mjdCol = mjdCol 

143 self.expTimeCol = expTimeCol 

144 self.minExpTime = minExpTime 

145 self.nightCol = nightCol 

146 units = 'nights' 

147 super().__init__(col=[self.mjdCol, self.expTimeCol, self.nightCol], units=units, **kwargs) 

148 

149 def run(self, dataSlice, slicePoint): 

150 # Order data Slice/times and exclude visits which are too short. 

151 long = np.where(dataSlice[self.expTimeCol] > self.minExpTime) 

152 data = np.sort(dataSlice[long], order=self.mjdCol) 

153 # SlicePoints ra/dec are always in radians - convert to degrees to calculate season 

154 seasons = calcSeason(np.degrees(slicePoint['ra']), data[self.mjdCol]) 

155 firstOfSeason, lastOfSeason = findSeasonEdges(seasons) 

156 seasonMeans = np.zeros(len(firstOfSeason), float) 

157 for i, (first, last) in enumerate(zip(firstOfSeason, lastOfSeason)): 

158 if first < last: 

159 n = data[self.nightCol][first:last+1] 

160 deltaNights = np.diff(np.unique(n)) 

161 if len(deltaNights) > 0: 

162 seasonMeans[i] = np.mean(deltaNights) 

163 return np.mean(seasonMeans) 

164 

165 

166class TdcMetric(BaseMetric): 

167 """Calculate the Time Delay Challenge metric, as described in Liao et al 2015 

168 (https://arxiv.org/pdf/1409.1254.pdf). 

169 

170 This combines the MeanCampaignFrequency/MeanNightSeparation, the SeasonLength, and the CampaignLength 

171 metrics above, but rewritten to calculate season information only once. 

172 

173 cadNorm = in units of days 

174 seaNorm = in units of months 

175 campNorm = in units of years 

176 """ 

177 def __init__(self, mjdCol='observationStartMJD', nightCol='night', 

178 expTimeCol='visitExposureTime', minExpTime=20, 

179 metricName = 'TDC', cadNorm=3., seaNorm=4., campNorm=5., badval=-999, **kwargs): 

180 # Save the normalization values. 

181 self.cadNorm = cadNorm 

182 self.seaNorm = seaNorm 

183 self.campNorm = campNorm 

184 self.mjdCol = mjdCol 

185 self.expTimeCol = expTimeCol 

186 self.minExpTime = minExpTime 

187 self.nightCol = nightCol 

188 super().__init__(col=[self.mjdCol, self.expTimeCol, self.nightCol], badval=badval, 

189 metricName = metricName, units = '%s' %('%'), **kwargs) 

190 

191 def run(self, dataSlice, slicePoint): 

192 # Order data Slice/times and exclude visits which are too short. 

193 long = np.where(dataSlice[self.expTimeCol] > self.minExpTime) 

194 data = np.sort(dataSlice[long], order=self.mjdCol) 

195 # SlicePoints ra/dec are always in radians - convert to degrees to calculate season 

196 seasons = calcSeason(np.degrees(slicePoint['ra']), data[self.mjdCol]) 

197 intSeasons = np.floor(seasons) 

198 firstOfSeason, lastOfSeason = findSeasonEdges(seasons) 

199 # Campaign length 

200 camp = len(np.unique(intSeasons)) 

201 # Season length 

202 seasonlengths = data[self.mjdCol][lastOfSeason] - data[self.mjdCol][firstOfSeason] 

203 sea = np.median(seasonlengths) 

204 # Convert to months 

205 sea = sea / 30.0 

206 # Campaign frequency / mean night separation 

207 seasonMeans = np.zeros(len(firstOfSeason), float) 

208 for i, (first, last) in enumerate(zip(firstOfSeason, lastOfSeason)): 

209 n = data[self.nightCol][first:last+1] 

210 deltaNights = np.diff(np.unique(n)) 

211 if len(deltaNights) > 0: 

212 seasonMeans[i] = np.mean(deltaNights) 

213 cad = np.mean(seasonMeans) 

214 # Evaluate precision and accuracy for TDC 

215 if sea == 0 or cad == 0 or camp == 0: 

216 return self.badval 

217 else: 

218 accuracy = 0.06 * (self.seaNorm / sea) * \ 

219 (self.campNorm / camp)**(1.1) 

220 precision = 4.0 * (cad / self.cadNorm)**(0.7) * \ 

221 (self.seaNorm/sea)**(0.3) * \ 

222 (self.campNorm / camp)**(0.6) 

223 rate = 30. * (self.cadNorm / cad)**(0.4) * \ 

224 (sea / self.seaNorm)**(0.8) * \ 

225 (self.campNorm / camp)**(0.2) 

226 return {'accuracy':accuracy, 'precision':precision, 'rate':rate, 

227 'cadence':cad, 'season':sea, 'campaign':camp} 

228 

229 def reduceAccuracy(self, metricValue): 

230 return metricValue['accuracy'] 

231 

232 def reducePrecision(self, metricValue): 

233 return metricValue['precision'] 

234 

235 def reduceRate(self, metricValue): 

236 return metricValue['rate'] 

237 

238 def reduceCadence(self, metricValue): 

239 return metricValue['cadence'] 

240 

241 def reduceSeason(self, metricValue): 

242 return metricValue['season'] 

243 

244 def reduceCampaign(self, metricValue): 

245 return metricValue['campaign']