Package wsatools :: Module QualityControl
[hide private]

Source Code for Module wsatools.QualityControl

   1  #------------------------------------------------------------------------------ 
   2  #$Id: QualityControl.py 9878 2013-06-04 13:16:17Z RossCollins $ 
   3  """ 
   4     Routines to perform automatic quality control on the metadata (via the 
   5     deprecated flag) and the detections (via the ppErrBits flag). 
   6   
   7     @todo: Make bit flags module constants. 
   8     @todo: Seriously consider changing cross-talk to avoid need for detections 
   9            of TWOMASS stars to exist in the WFCAM field for MosaicCrossTalkFlag 
  10            at least - if not others too. Would need to use min/maxRa/Dec from 
  11            CurrentAstrometry table and then convert TWOMASS RA, Decs to x, y 
  12            using pywcs. 
  13     @todo: Add statistics of percentage flagged with each quality issue. 
  14     @todo: Make use of DbSession.update() accepting Join objects as the 
  15            fromTables argument. 
  16     @todo: Additional quality flags (trac tickets: 13, 31). 
  17   
  18     @author: R.S. Collins 
  19     @org:    WFAU, IfA, University of Edinburgh 
  20  """ 
  21  #------------------------------------------------------------------------------ 
  22  from __future__      import division, print_function 
  23  from future_builtins import map, zip 
  24   
  25  from   collections import namedtuple 
  26  from   bisect      import bisect 
  27  import math 
  28  from   operator    import attrgetter 
  29   
  30  import wsatools.Astrometry              as astro 
  31  from   wsatools.DbConnect.CuSession import CuSession 
  32  import wsatools.DbConnect.DbConstants   as dbc 
  33  from   wsatools.DbConnect.DbSession import DbSession, Join, SelectSQL 
  34  import wsatools.FitsUtils               as fits 
  35  from   wsatools.Logger              import Logger, ForLoopMonitor 
  36  from   wsatools.SystemConstants     import DepCodes, SystemConstants 
  37  import wsatools.Utilities               as utils 
  38  #------------------------------------------------------------------------------ 
  39   
40 -class CrossTalkFlagger(object):
41 """ 42 Defines the cross-talk artefact flagging algorithm and process. 43 44 @todo: Complete reliability tests (trac ticket:32) 45 46 @group Errors & Exceptions: MissingStarError 47 48 """
49 - class MissingStarError(Exception):
50 """ Exception thrown if a 2MASS star can't be found in a WFCAM image. 51 """
52 - def __init__(self, ra, dec, j_m, curDet):
53 msg = "2MASS star at (%s, %s) of j_m = %s has no WFCAM detection "\ 54 "in %s" % (ra, dec, j_m, curDet) 55 56 super(CrossTalkFlagger.MissingStarError, self).__init__(msg)
57 58 #-------------------------------------------------------------------------- 59 #: Default uncertainty in astrometry due to rounding (degrees). 60 roundingErr = 0.0002 61 #: Search radius limit for non-detected 2MASS stars (degrees). 62 searchErr = 0.003 63 #: 2MASS J-band magnitude limit for star to create cross-talk. 64 magLim = 14 65 #: Field bright star density limit for cross-talk flagging. 66 densityLimit = 150 67 68 #: Only flag fields beneath the density limit? 69 avoidDenseFields = False 70 #: Hexadecimal bit mask string to denote possible cross-talk artefacts. 71 bitMask = None 72 #: Detection table to be flagged. 73 _detTable = None 74 #: Connection to 2MASS database. 75 _twomass = None 76 #: The calling QualityFlagUpdater object. 77 _updater = None 78 79 #-------------------------------------------------------------------------- 80
81 - def __init__(self, bitMask, updater):
82 """ 83 Initialise cross-talk flagger. Prepares connection to 2MASS server. 84 85 @param bitMask: Hexadecimal bit mask string used to denote possible 86 cross-talk artefacts. 87 @type bitMask: str 88 @param updater: The calling QualityFlagUpdater object. 89 @type updater: QualityFlagUpdater 90 91 """ 92 self.bitMask = bitMask 93 self._updater = updater 94 self._detTable = self._updater.programme.getDetectionTable() 95 self.avoidDenseFields = self._updater.programme.isNonSurvey() 96 twomassDbName = self._updater.archive.query("databaseName", 97 "ExternalSurvey", whereStr="surveyName='2MASS'", firstOnly=True) 98 99 self._twomass = DbSession(twomassDbName, 100 userName=dbc.loadServerRwUsername())
101 102 #-------------------------------------------------------------------------- 103
104 - def run(self):
105 """ Flag cross-talk artefacts. 106 """ 107 Logger.addMessage("Determining list of stacks...") 108 multiframes = self._getMultiframeList() 109 Logger.addMessage("%s multiframes will be flagged" % len(multiframes)) 110 111 numArtefacts = 0 112 numFrames = 0 113 progress = ForLoopMonitor(multiframes) 114 for mfID, nuStep in multiframes: 115 pixFactor = math.sqrt(nuStep) 116 artSpacing = 128 * int(pixFactor) 117 numFailures = 0 118 119 # Work through the non-deprecated detector images only, this also 120 # avoids assumption that every multiframe has a fixed number of 121 # detector images and constant xSize, ySize values! 122 detectors = self._updater.archive.query( 123 selectStr="MultiframeDetector.extNum, axis1Length, axis2Length, " 124 "detRows, detCols, xPixSize", 125 fromStr=Join(["MultiframeDetector", "CurrentAstrometry"], 126 ["multiframeID", "extNum"]), 127 whereStr="MultiframeDetector.multiframeID=%s AND %s" 128 % (mfID, DepCodes.selectNonDeprecated)) 129 130 isMfFlagable = True 131 mfArtefacts = [] 132 for extNum, xSize, ySize, detRows, detCols, xPixSize in detectors: 133 curDet = "multiframeID = %s AND extNum = %s" % (mfID, extNum) 134 xOffset = xSize - pixFactor * detRows 135 yOffset = ySize - pixFactor * detCols 136 137 brightStars = self._getBrightStars(curDet) 138 139 isMfFlagable = not self.avoidDenseFields \ 140 or len(brightStars) < self.densityLimit 141 142 if not isMfFlagable: 143 Logger.addMessage("<Info> Too many bright stars in " 144 "multiframeID=%s field to flag cross-talk: %s" % 145 (mfID, len(brightStars)), alwaysLog=False) 146 break 147 148 artefacts = [] 149 for ra, dec, j_m in brightStars: 150 try: 151 # Determine position on WFCAM detector 152 x, y, artErr = \ 153 self._locateBrightStar(ra, dec, j_m, curDet) 154 except CrossTalkFlagger.MissingStarError as error: 155 numFailures += 1 156 157 # Message is distracting to curator but useful for 158 # quality control forensics in the future. 159 Logger.addExceptionWarning(error, echoOff=True) 160 else: 161 artRad = self._calcArtefactRadius(j_m, artErr, 162 xPixSize, pixFactor) 163 artefacts += self._createArtefactList(x, y, artRad, 164 artSpacing, xSize, ySize, xOffset, yOffset) 165 166 mfArtefacts.append(artefacts) 167 168 if isMfFlagable: 169 for detector, artefacts in zip(detectors, mfArtefacts): 170 extNum = detector[0] 171 numArtefacts += self._flagDetector(artefacts, mfID, extNum) 172 173 self._commitMultiframe(mfID) 174 175 if not numFailures: 176 numFrames += 1 177 178 progress.testForOutput("Completed multiframeID %s" % mfID) 179 180 if numArtefacts: 181 Logger.addMessage("done. %s sources affected (these may not be " 182 "unique sources, some may be affected more than " 183 "once)" % numArtefacts) 184 185 Logger.addMessage("%s frames complete." % numFrames) 186 Logger.addMessage("Flagging incomplete for %s frames." 187 % (len(multiframes) - numFrames))
188 189 #-------------------------------------------------------------------------- 190
191 - def _getMultiframeList(self):
192 """ 193 Works through all of the latest multiframe versions that contribute to 194 this detection table, where cross-talk flag has not already been set. 195 196 @return: List of multiframeIDs with number of microsteps. 197 @rtype: list(tuple(int, int)) 198 199 """ 200 attrs = "Multiframe.multiframeID, nuStep" 201 202 return self._updater.archive.query( 203 selectStr=attrs, 204 fromStr=Join(["Multiframe", self._detTable, "ProgrammeFrame"], 205 ["multiframeID"]), 206 whereStr="frameType LIKE '%stack' AND frameType NOT LIKE 'mosaic%'" 207 + " AND Multiframe.%s AND programmeID = %s" 208 " AND ppErrBitsStatus & %s = 0" 209 % (DepCodes.selectNonDeprecated, self._updater.programmeID, 210 self.bitMask), 211 groupBy=attrs)
212 213 #-------------------------------------------------------------------------- 214
215 - def _getBrightStars(self, curDet):
216 """ 217 @return: A list of bright star positions and J-band magnitudes in the 218 current detector field from the 2MASS catalogue. 219 @rtype: list(tuple(float, float, float) 220 221 """ 222 whereStr = curDet + " AND seqNum > 0" 223 224 # Determine range of detector's field for 2MASS search 225 ra_min, ra_max, dec_min, dec_max = \ 226 self._updater.archive.query("min(ra), max(ra), min(dec), max(dec)", 227 self._updater.programme.getAstrometryTable(), whereStr)[0] 228 229 if not ra_min: 230 Logger.addMessage("<Warning> " + curDet + 231 " seems to have no detections! Row count gives %s" % 232 self._updater.archive.queryNumRows( 233 self._updater.programme.getAstrometryTable(), whereStr)) 234 return [] 235 else: 236 # Find bright 2MASS stars within this detector field 237 # beware RA 0h boundary 238 if ra_min < 1 and ra_max > 359: 239 detRange = ("ra BETWEEN %s AND %s " % (0, ra_min) + 240 "AND ra BETWEEN %s AND %s " % (ra_max, 360)) 241 else: 242 detRange = "ra BETWEEN %s AND %s " % (ra_min, ra_max) 243 detRange += " AND dec BETWEEN %s AND %s" % (dec_min, dec_max) 244 return self._twomass.query("ra, dec, j_m", 245 dbc.twomassSourceTableName(), detRange + 246 " AND j_m BETWEEN -10 AND %s" % self.magLim)
247 248 #-------------------------------------------------------------------------- 249
250 - def _locateBrightStar(self, ra, dec, j_m, curDet):
251 """ 252 Finds the location of a bright star on the WFCAM detector given the RA, 253 Dec coordinates, by finding the equivalent WFCAM detection. 254 255 @return: Pixel coordinate of star on WFCAM detector, and error in 256 degrees. 257 @rtype: tuple(float, float, float) 258 259 """ 260 artErr = 0 261 posErr = CrossTalkFlagger.roundingErr 262 starPos = "(ra BETWEEN %s AND %s) AND (dec BETWEEN %s AND %s)" % ( 263 ra - posErr, ra + posErr, dec - posErr, dec + posErr) 264 detJoin = Join([self._detTable, 265 self._updater.programme.getAstrometryTable()], 266 ["multiframeID", "extNum", "seqNum"]) 267 curDetJoin = \ 268 curDet.replace("multiframeID", self._detTable + ".multiframeID")\ 269 .replace("extNum", self._detTable + ".extNum") 270 271 x, y = self._updater.archive.query( 272 selectStr="x, y", 273 fromStr=detJoin, 274 whereStr=curDetJoin + " AND " + starPos, 275 firstOnly=True) 276 277 if x is None: 278 # Widen the search area 279 posErr = CrossTalkFlagger.searchErr 280 starPos = "(ra BETWEEN %s AND %s) AND (dec BETWEEN %s AND %s)" % ( 281 ra - posErr, ra + posErr, dec - posErr, dec + posErr) 282 detList = self._updater.archive.query("x, y, ra, dec", 283 detJoin, curDetJoin + " AND " + starPos) 284 # Find the closest detection using simple maths 285 try: 286 err, x, y = min(((ra - ra1) ** 2 + (dec - dec1) ** 2, x, y) 287 for x, y, ra1, dec1 in detList) 288 except ValueError: 289 raise CrossTalkFlagger.MissingStarError(ra, dec, j_m, curDet) 290 # Determine wider search field for artefacts due to offset in 291 # source position 292 artErr = min(math.sqrt(err), posErr) 293 return x, y, artErr
294 295 #-------------------------------------------------------------------------- 296
297 - def _calcArtefactRadius(self, j_m, artErr, xPixSize, pixFactor):
298 """ 299 Returns the radius in pixels for a first-order cross-talk artefact 300 produced by a star of given brightness (j_m). This radius is increased 301 by artErr to consider the offset between the 2MASS position of the star 302 and the WFCAM position. The radius is increased instead of just 303 shifting centre, because this offset could be due to a 2MASS object 304 becoming decoupled in the WFCAM image into multiple sources. 305 306 @return: Radius of artefact in pixels. 307 @rtype: float 308 309 """ 310 # Convert error from degrees to pixels 311 artErr *= 3600.0 / xPixSize 312 # Calculate radius of artefact limited to 20 pixels 313 return min(40, 10 ** (2.227 - 0.123 * j_m)) * pixFactor / 2.0 + artErr
314 315 #-------------------------------------------------------------------------- 316
317 - def _createArtefactList(self, x, y, artRad, artSpacing, xSize, ySize, 318 xOffset, yOffset):
319 """ 320 Creates a list of artefact positions along each direction from a bright 321 cross-talk producing star. Having first determined which directions 322 artefacts would be produced by a bright-star in this quadrant of the 323 detector, considering the effects of dithering too. 324 325 @return: List of SQL where clause predicates for each strip of 326 cross-talk artefacts for each direction. 327 @rtype: list(str) 328 329 """ 330 # Determine direction of cross-talk 331 # @@TODO: May need to add artRad to x/yOffset... 332 if x < xSize / 2 + xOffset and x > xSize / 2 - xOffset or \ 333 y < ySize / 2 + yOffset and y > ySize / 2 - yOffset: 334 isHoriCrossTalk = True 335 isVertCrossTalk = True 336 else: 337 isHoriCrossTalk = x < xSize / 2 and y < ySize / 2 or \ 338 x >= xSize / 2 and y >= ySize / 2 339 isVertCrossTalk = not isHoriCrossTalk 340 341 # Determine list of cross-talk artefacts regions 342 artefactStrips = [] 343 if isHoriCrossTalk: 344 artefactStrips += self._createStripList(x, y, artRad, artSpacing, 345 xSize, isXY=True) 346 if isVertCrossTalk: 347 artefactStrips += self._createStripList(y, x, artRad, artSpacing, 348 ySize, isXY=False) 349 return artefactStrips
350 351 #-------------------------------------------------------------------------- 352
353 - def _createStripList(self, x1, x2, artRad, artSpacing, size, isXY):
354 """ 355 Determines the list of artefacts in the one-dimension (both directions) 356 from a given position of a bright star on the detector image. The 357 algorithm keeps appending cross-talk artefact positions to the list 358 until the edge of each the detector quadrant is reached. 359 360 @param x1: Position of bright cross-talk producing star in the 361 dimension along-which artefacts will be produced for this 362 position on the detector. 363 @type x1: float 364 @param x2: Position of bright cross-talk producing star in the other 365 dimension. 366 @type x2: float 367 @param artRad: Radius of a cross-talk artefact in pixels. 368 @type artRad: float 369 @param artSpacing: Absolute spacing in pixels between artefacts in this 370 detector image. 371 @type artSpacing: float 372 @param size: Width of the detector image along the x1-dimension. 373 @type size: float 374 @param isXY: If True, then x1=x, x2=y, to find cross-talk artefacts 375 along the x-dimension of the detector, otherwise x1=y, 376 x2=x to find artefacts along the y-dimension. 377 @type isXY: bool 378 379 @return: List containing a single SQL where clause predicate that 380 specifies a strip of artefact positions along one-dimension, 381 with the same radial width. 382 @rtype: list(str) 383 384 """ 385 #@@TODO: These may not be true the limits! 386 minPos, maxPos = ((0, size / 2) if x1 < size / 2 else (size / 2, size)) 387 artefacts = [] 388 predStr = ("x" if isXY else "y") + " BETWEEN %s AND %s" 389 390 for n in range(int(abs(maxPos - x1) / artSpacing)): 391 i = x1 + (n + 1) * artSpacing 392 artefacts.append(predStr % (i - artRad, i + artRad)) 393 394 for n in range(int(abs(x1 - minPos) / artSpacing)): 395 i = x1 - (n + 1) * artSpacing 396 artefacts.append(predStr % (i - artRad, i + artRad)) 397 398 if artefacts: 399 return ["(%s BETWEEN %s AND %s AND (%s))" % ("y" if isXY else "x", 400 x2 - artRad, x2 + artRad, " OR ".join(artefacts))] 401 else: 402 return []
403 404 #-------------------------------------------------------------------------- 405
406 - def _flagDetector(self, artefacts, mfID, extNum):
407 """ 408 Flags possible cross-talk artefacts in the detection table, given a 409 list of artefacts positions as SQL where clause predicates. 410 411 @return: Number of detections flagged. 412 @rtype: int 413 414 """ 415 # Divide long SQL update statements into smaller ones 416 # This number was arbitrarily chosen but fairly optimal 417 numArtefacts = 0 418 for stripUpdates in utils.splitList(artefacts, 256): 419 numArtefacts += self._updater.archive.update(self._detTable, 420 [("ppErrBits", "ppErrBits | " + self.bitMask)], 421 where="multiframeID = %s AND extNum = %s AND (%s)" 422 % (mfID, extNum, " OR ".join(stripUpdates))) 423 return numArtefacts
424 425 #-------------------------------------------------------------------------- 426
427 - def _commitMultiframe(self, mfID):
428 """ 429 Commits all updates to the detection table for a given multiframe and 430 updates the relevant ppErrBitsStatus flag in ProgrammeTable. 431 432 """ 433 self._updater.updateStatus(self.bitMask, mfID) 434 # Don't want log to be too big! 435 self._updater.archive.commitTransaction()
436 437 #------------------------------------------------------------------------------ 438
439 -class MosaicCrossTalkFlagger(CrossTalkFlagger):
440 """ Derivation of the CrossTalkFlagger for mosaic images. 441 """ 442 Artefact = namedtuple("Artefact", "position radius") 443 """ A cross-talk artefact defined by its (x, y) position and radius. 444 """ 445 CrossTalk = namedtuple("CrossTalk", "j_m radii") 446 """ Definition of cross-talk produced by star of brightness j_m and 447 expected radii (in pixels) of artefacts of each order. 448 """ 449 #-------------------------------------------------------------------------- 450 451 artPosErr = 0.2 # arcseconds (1 pixel = 0.134") 452 """ Fuzziness of artefact positions. This is used to overcome rounding 453 errors in artefact positions, so that an artefact at the same position 454 and size from many frames is not duplicated. This reduces the length 455 of the SQL statement required to update all mosaic artefact positions. 456 This shouldn't be too large else the flagging becomes inaccurate, and 457 so is empirically determined by reducing precision until the reduction 458 in the number of duplicate artefacts becomes small and the change in 459 the number of flagged artefacts is minimised. 460 """ 461 magLim = 17 # Overrides CrossTalkFlagger value 462 #: UDS artefact radius grid. Radius in pixels for given j_m. 463 udsArtRads = [CrossTalk(7.0, [50, 50, 50, 50, 50, 50, 50]), 464 CrossTalk(8.5, [30, 25, 25, 25, 25, 25, 25]), 465 CrossTalk(9.5, [20, 15, 15, 15, 15, 15, 15]), 466 CrossTalk(11.5, [15, 10, 10, 6, 6, 6, 6]), 467 CrossTalk(12.5, [12, 6, 6, 6, 6, 3, 3]), 468 CrossTalk(13.0, [10, 6, 6, 3, 3, 3, 0]), 469 CrossTalk(13.5, [10, 6, 3, 3, 3, 0, 0]), 470 CrossTalk(14.0, [8, 6, 3, 3, 0, 0, 0]), 471 CrossTalk(14.5, [6, 3, 3, 0, 0, 0, 0]), 472 CrossTalk(15.0, [6, 3, 0, 0, 0, 0, 0]), 473 CrossTalk(16.0, [5, 0, 0, 0, 0, 0, 0]), 474 CrossTalk(17.0, [3, 0, 0, 0, 0, 0, 0])] 475 476 _mosaicArtefacts = None #: Set of artefacts identified in mosaic. 477 _mosaicMfID = None #: MultiframeID of mosaic currently being flagged. 478 479 #-------------------------------------------------------------------------- 480
481 - def run(self):
482 """ 483 Runs the cross-talk flagging algorithm on every mosaic image that has 484 unflagged detections. 485 486 """ 487 artPrec = 3600 / MosaicCrossTalkFlagger.artPosErr 488 489 detJoin = Join( 490 [self._detTable, self._updater.programme.getAstrometryTable()], 491 ["multiframeID", "extNum", "seqNum"]) 492 493 Logger.addMessage("Determining list of mosaics...") 494 495 # Find all mosaics where cross-talk flag has not already been set. 496 mosaics = self._updater.archive.query( 497 selectStr="Multiframe.multiframeID", 498 fromStr=Join(["Multiframe", "ProgrammeFrame"], ["multiframeID"]), 499 whereStr="frameType LIKE 'mosaicdeep%%stack' AND %s AND " 500 "programmeID = %s AND ppErrBitsStatus & %s = 0" % 501 (DepCodes.selectNonDeprecated, 502 self._updater.programme.getAttr("programmeID"), 503 self.bitMask)) 504 505 for mosaicNum, self._mosaicMfID in enumerate(mosaics): 506 Logger.addMessage( 507 "Flagging mosaic %s of %s" % (mosaicNum + 1, len(mosaics))) 508 509 self._mosaicArtefacts = set() 510 super(MosaicCrossTalkFlagger, self).run() 511 512 # Reduce precision of artefact positions to integer values in 513 # units of artPosErr, so that rounding errors in coordinate 514 # transforms don't lead to too many duplicate identifications 515 # of the same artefact. 516 self._mosaicArtefacts = set( 517 (int(ra * artPrec), int(dec * artPrec), int(artRad * artPrec)) 518 for (ra, dec), artRad in self._mosaicArtefacts) 519 520 Logger.addMessage("Identified %s non-unique artefact positions" % 521 len(self._mosaicArtefacts)) 522 523 Logger.addMessage("Updating mosaic...") 524 numFlagged = 0 525 for artefacts in utils.splitList(list(self._mosaicArtefacts), 250): 526 artefactBoxes = " OR ".join( 527 "ra BETWEEN %s AND %s AND dec BETWEEN %s AND %s" % 528 ((ra - artRad) / artPrec, (ra + artRad) / artPrec, 529 (dec - artRad) / artPrec, (dec + artRad) / artPrec) 530 for ra, dec, artRad in artefacts) 531 numFlagged += self._updater.archive.update(self._detTable, 532 [("ppErrBits", "ppErrBits | " + self.bitMask)], 533 where="%s.multiframeID=%s AND (%s)" % 534 (self._detTable, self._mosaicMfID, artefactBoxes), 535 fromTables=detJoin) 536 537 self._updater.updateStatus(self.bitMask, self._mosaicMfID) 538 self._updater.archive.commitTransaction() 539 Logger.addMessage( 540 "done. %s non-unique sources affected" % numFlagged) 541 542 if mosaics: 543 Logger.addMessage("done. All mosaics flagged.") 544 else: 545 Logger.addMessage("<Info> No non-deprecated mosaics to flag.")
546 547 #-------------------------------------------------------------------------- 548
549 - def _getMultiframeList(self):
550 """ 551 Finds component multiframes of the mosaic currently being flagged. 552 553 @return: List of multiframeIDs with number of microsteps. 554 @rtype: list(tuple(int, int)) 555 556 """ 557 attrs = "Multiframe.multiframeID, nuStep" 558 559 return self._updater.archive.query( 560 selectStr=attrs, 561 fromStr="Multiframe, Provenance AS P", 562 whereStr="P.combiframeID = %s AND " % self._mosaicMfID + 563 "P.multiframeID = Multiframe.multiframeID " 564 "GROUP BY %s" % attrs)
565 566 #-------------------------------------------------------------------------- 567
568 - def _calcArtefactRadius(self, j_m, _artErr, xPixSize, _pixFactor):
569 """ 570 Returns the radii in degrees for each order of cross-talk artefact for 571 a star of given brightness (j_m). The error in the position of the 572 bright star, artErr, is not considered for the UDS mosaic as this can 573 be too large for bright stars, and the image contains no examples of 574 2MASS stars decoupled into multiple point sources by WFCAM. 575 576 @return: List of radii in degrees for each order of artefact. 577 @rtype: list(float) 578 579 """ 580 magIndex = bisect([star.j_m for star in self.udsArtRads], j_m) 581 return [radius * xPixSize / 3600.0 582 for radius in self.udsArtRads[magIndex].radii]
583 584 #-------------------------------------------------------------------------- 585
586 - def _createStripList(self, x1, x2, artRad, artSpacing, size, isXY):
587 """ 588 Determines the list of artefacts in the one-dimension (both directions) 589 from a given position of a bright star on the detector image. The 590 algorithm keeps appending cross-talk artefact positions to the list 591 until the edge of each the detector quadrant is reached. 592 593 @param x1: Position of bright cross-talk producing star in the 594 dimension along-which artefacts will be produced for this 595 position on the detector. 596 @type x1: float 597 @param x2: Position of bright cross-talk producing star in the other 598 dimension. 599 @type x2: float 600 @param artRad: List of cross-talk artefact radii (pixels) of each 601 order. 602 @type artRad: list(float) 603 @param artSpacing: Absolute spacing in pixels between artefacts in this 604 detector image. 605 @type artSpacing: float 606 @param size: Width of the detector image along the x1-dimension. 607 @type size: float 608 @param isXY: If True, then x1=x, x2=y, to find cross-talk artefacts 609 along the x-dimension of the detector, otherwise x1=y, 610 x2=x to find artefacts along the y-dimension. 611 @type isXY: bool 612 613 @return: List of artefact (x, y) positions and radii (degrees) for a 614 single one-dimensional strip. 615 @rtype: list(Artefact(tuple(float, float), float)) 616 617 """ 618 artefacts = [] 619 #@@TODO: This may not be the limits! 620 minPos, maxPos = ((0, size / 2) if x1 < size / 2 else (size / 2, size)) 621 622 for n in range(int(abs(maxPos - x1) / artSpacing)): 623 if artRad[min(6, n)] > 0: 624 i = x1 + (n + 1) * artSpacing 625 x, y = (i, x2) if isXY else (x2, i) 626 artefacts.append( 627 MosaicCrossTalkFlagger.Artefact((x, y), artRad[min(6, n)])) 628 629 for n in range(int(abs(x1 - minPos) / artSpacing)): 630 if artRad[min(6, n)] > 0: 631 i = x1 - (n + 1) * artSpacing 632 x, y = (i, x2) if isXY else (x2, i) 633 artefacts.append( 634 MosaicCrossTalkFlagger.Artefact((x, y), artRad[min(6, n)])) 635 636 return artefacts
637 638 #-------------------------------------------------------------------------- 639
640 - def _flagDetector(self, artefacts, mfID, extNum):
641 """ 642 Detections from individual component frames of the mosaic are not 643 flagged. So this function just maintains a grand total list of artefact 644 positions on mosaic image to be updated. All positions from the same 645 detector are converted from pixel (X, Y) to sky (RA, Dec) coordinates 646 (for that detector's WCS values), which are used to flag the mosaic 647 detections. 648 649 """ 650 self._mosaicArtefacts.update(list(zip( 651 astro.convPixelToSky(map(attrgetter("position"), artefacts), 652 mfID, extNum, db=self._updater.archive), 653 map(attrgetter("radius"), artefacts)))) 654 655 return 0
656 657 #-------------------------------------------------------------------------- 658
659 - def _commitMultiframe(self, mfID):
660 """ No need to commit individual stack frames for the mosaic flagger. 661 """ 662 pass
663 664 #------------------------------------------------------------------------------ 665
666 -class NonSurveyDeprecator(CuSession):
667 """ Deprecates non-survey data according to quality control criteria. 668 If supplied with a specific programmeID then deprecations are just 669 propagated to that programme's detection table, otherwise just 670 frame and multiframe data for all non-surveys are quality controlled. 671 """ 672 #-------------------------------------------------------------------------- 673 # Private class parameters for this procedure - should not be altered 674 675 _autoCommit = True # Overrides CuSession default 676 _isPersistent = True # Overrides CuSession default 677 cuNum = 20 # Overrides CuSession default 678 679 #-------------------------------------------------------------------------- 680 # Define public member variable default values (access as obj.varName) 681 # these need to be set from command-line options 682 683 #: Final date of data to quality control. 684 endDate = CuSession.sysc.obsCal.lateDefaultDateTime.date 685 686 #-------------------------------------------------------------------------- 687
688 - def _onRun(self):
689 """ Perform quality control to deprecate data in the database. 690 """ 691 Logger.addMessage("Applying Quality Control to Non-Survey Data.") 692 693 # Prepare correct list of programmeIDs for ProgrammeCurationHistory 694 if self.programmeID: 695 progIDs = [self.programmeID] 696 else: 697 progIDs = self.programme.getProgIDList(self.onlyNonSurveys) 698 self.programmeID = set(progIDs) 699 700 if any(progID < dbc.calibrationProgrammeID() for progID in progIDs): 701 raise NonSurveyDeprecator.CuError( 702 "Cannot curate programme %s. This script should only be run for" 703 "non-survey programmes." % min(progIDs)) 704 705 if type(self.programmeID) is set: 706 self.comment = "Performing quality control on all non-survey "\ 707 "frames and multiframes" 708 self.deprecateFrames() 709 else: 710 self.comment = "Propagating frame deprecations into %s" % \ 711 self.programme.getDetectionTable() 712 self.deprecateDetections()
713 714 #-------------------------------------------------------------------------- 715
716 - def deprecateDetections(self):
717 """ Propagate deprecations to detection table of current programme. 718 """ 719 detTable = self.programme.getDetectionTable() 720 Logger.addMessage( 721 "Propagating deprecations to detection table %s..." % detTable) 722 723 numRows = self.archive.update(detTable, 724 entryList="deprecated=%s" % DepCodes.propDepDet, 725 where="MFD.%s AND D.deprecated=%s AND " 726 "D.multiframeID=MFD.multiframeID AND D.extNum=MFD.extNum" 727 % (DepCodes.selectDeprecated, DepCodes.nonDep), 728 fromTables="%s AS D, MultiframeDetector AS MFD" % detTable) 729 730 Logger.addMessage("%s detections deprecated." % numRows)
731 732 #-------------------------------------------------------------------------- 733
734 - def deprecateFrames(self):
735 """ Deprecate frames and multiframes for all non-survey data only. 736 """ 737 Logger.addMessage("Deprecating frames and multiframes...") 738 endDate = "'%s'" % self.endDate 739 qcSql = [] 740 sql = [] 741 utils.ParsedFile.commentMarker = '--' 742 for line in utils.ParsedFile(self.sysc.nsQcScript): 743 if not line.startswith('GO'): 744 sql.append(line.lstrip()) 745 else: 746 qcSql.append(' '.join(sql).replace("@endDate", endDate)) 747 sql = [] 748 749 qcSteps = [] 750 if self.sysc.isWSA(): 751 qcSteps += [(DepCodes.ignoredObs, 752 "Frame is flagged as to be ignored in pipeline processing")] 753 754 qcSteps += [ 755 (DepCodes.rejectObs, 756 "Frame is flagged as part of a summit-rejected MSB"), 757 (DepCodes.missingCat, 758 "Stack frames that have no catalogues"), 759 (DepCodes.invalidDet, 760 "Dead detector frames or all bad channels"), 761 (DepCodes.invalidDet, 762 "Dead detector frames or all bad channels: propagating forward 1/2"), 763 (DepCodes.invalidDet, 764 "Dead detector frames or all bad channels: propagating forward 2/2"), 765 (DepCodes.invalidMet, 766 "Undefined and/or nonsensical critical image metadata attributes"), 767 (DepCodes.incompGrp, 768 "Insufficient frame complements within groups"), 769 (DepCodes.invalidCat, 770 "Undefined values of critical catalogue attributes for stacks"), 771 (DepCodes.invalidSee, 772 "Seeing <= 0.0 for stacks"), 773 (DepCodes.incompProv, 774 "Provenance checks: leavs"), 775 (DepCodes.propDepCom, 776 "Leav formed from a deprecated component frame"), 777 (DepCodes.propDepCom, 778 "Stack formed from a deprecated component frame"), 779 (DepCodes.incompProv, 780 "Provenance checks: stacks")] 781 782 if self.sysc.isVSA(): # @TODO: Apply to OSA too? 783 qcSteps += [(DepCodes.poorEsoGrade, 784 "Esograde check for 'C', 'D', 'R': stacks")] 785 786 qcSteps += [ 787 (DepCodes.propDepFra, 788 "Multiframes with all detectors deprecated"), 789 (DepCodes.propDepObs, 790 "Propagate multiframe deprecations to detectors")] 791 792 for step, sql in zip(qcSteps, qcSql): 793 Logger.addMessage("Setting deprecated=%s: %s..." % step) 794 numRows = self.archive._executeScript(sql, wantRowCount=True) 795 Logger.addMessage("%s entries deprecated." % numRows)
796 797 #------------------------------------------------------------------------------ 798
799 -class QualityFlagUpdater(CuSession):
800 """ 801 Fills the ppErrBits attribute of the given programme's detection table. 802 803 Transactions are committed on a per bit flag basis, as there is an option 804 to just run for certain bit flags. The cross-talk flag is the exception, 805 committing on a per multiframe basis, and consequently also updates 806 ppErrBitsStatus on a per multiframe basis to ensure database consistency. 807 808 @note: Only the cross-talk flag skips detections from deprecated frames. 809 This is because skipping deprecated frames for the other, single, 810 updates will probably be slower as there is no index on the 811 deprecated flag. 812 """ 813 cuNum = 24 # Overrides CuSession default 814 815 #: Limit date range of flagging monthly tables to this value. 816 dateRange = SystemConstants().obsCal.dateRange() 817 #: Number of big flag to be updated on this run (-1 = all). 818 doBit = -1 819 #: Only flag detections from deep WFAU-products? 820 isDeepOnly = False 821 #: Skip mosaic provenance completeness check? 822 skipProvenanceCheck = False 823 #: Detection table currently being updated. 824 _detTable = None 825 #: Does this programme have mosaics that need flagging too? 826 _doMosaics = False 827 828 #-------------------------------------------------------------------------- 829
830 - def _onRun(self):
831 """ Update quality bit flags for current programme. 832 """ 833 Logger.addMessage("Determining what needs to be done...") 834 self._doMosaics = self.programme.isMosaicked() 835 836 # Ensure post-CU3 Provenance update has been run! 837 if self._doMosaics and not self.isMosaicProvenanceComplete(): 838 raise CuSession.CuError("ProvenanceFiller.py must be run first!") 839 840 Logger.addMessage("Updating quality bit flags for:") 841 842 detTables = self.programme.getMonthlyDetectionTables(self.dateRange, 843 self.isDeepOnly) 844 progress = ForLoopMonitor(detTables) 845 for self._detTable in detTables: 846 self.updateDeblendedFlag() 847 self.updateBadPixelFlag() 848 self.updateLowConfidenceFlag() 849 self.updateDet16Flag() 850 self.updateSaturatedFlag() 851 self.updateDitherFlag() 852 self.updatePoorFlatFieldFlag() # @GOTCHA: After dither flag! 853 self.updateUncorrectedFlag() 854 self.updateEarFlag() 855 self.updateMissingDetFlag() 856 progress.postMessage("Completed " + self._detTable) 857 858 # Cross-talk flag is last, as is the longest and commits transactions 859 # as it goes along. 860 if self.doBit == 19 or self.doBit == -1 and self.sysc.isWSA() \ 861 and self.programmeID not in (dbc.ptsProgrammeID(), 862 self.sysc.scienceProgs.get("GPS")): 863 864 Logger.addMessage("(Bit 19) Cross-talk artefacts...") 865 self.updateCrossTalkFlag(bitMask="0x00080000")
866 867 #-------------------------------------------------------------------------- 868
869 - def bitMask(self, bitNumber):
870 """ 871 Given a binary bit number, return the corresponding hexadecimal bit 872 mask. 873 874 @param bitNumber: Binary bit number. 875 @type bitNumber: int 876 877 @return: Hexadecimal bit mask. 878 @rtype: str 879 880 """ 881 bits = ['1' if bitNum == bitNumber else '0' for bitNum in range(32)] 882 883 return hex(int('0b' + ''.join(reversed(bits)), 2))
884 885 #-------------------------------------------------------------------------- 886
887 - def bitNumber(self, bitMask):
888 """ 889 Given a hexadecimal bit mask, return the corresponding binary bit 890 number string. 891 892 @param bitMask: Hexadecimal string. 893 @type bitMask: str 894 895 @return: Binary bit number. 896 @rtype: int 897 898 """ 899 return bin(int(bitMask, 16)).count('0') - 1
900 901 #-------------------------------------------------------------------------- 902
903 - def doBitFlag(self, bitNum, msg):
904 """ Checks whether we should apply this bit flag. 905 """ 906 if not (self.doBit == bitNum or self.doBit == -1): 907 return False 908 909 Logger.addMessage("(Bit %s) %s on %s..." 910 % (bitNum, msg, self._detTable)) 911 912 return True
913 914 #-------------------------------------------------------------------------- 915
916 - def isMosaicProvenanceComplete(self):
917 """ @return: True, if the mosaic Provenance is fully up-to-date. 918 @rtype: bool 919 """ 920 Logger.addMessage("Checking that Provenance is complete...") 921 922 # Number of mosaics in Multiframe with detections 923 fromTables = ("Multiframe AS M, %s AS D" % 924 self.programme.getDetectionTable()) 925 where = ("M.multiframeID=D.multiframeID AND " 926 "frameType LIKE '%mosaic%' AND M.deprecated=0") 927 numMosaics = self.archive.queryNumRows(fromTables, where, 928 distinctAttr="M.multiframeID") 929 930 # Number of combiframes with Provenance entries 931 fromTables += ", Provenance AS P" 932 where += " AND P.combiframeID=M.multiframeID" 933 numProvFrames = self.archive.queryNumRows( 934 fromTables, where, distinctAttr="M.multiframeID") 935 936 if numMosaics == numProvFrames and not self.skipProvenanceCheck: 937 # Also check the FITS files 938 catFiles = set(self.archive.query("catName", "Multiframe", 939 "frameType LIKE 'mosaicdeep%stack' AND deprecated=0")) 940 fileNames = [catFile.rpartition(':')[2] for catFile in catFiles] 941 fileNames = [fileName for fileName in fileNames 942 if fileName and fileName.endswith(self.sysc.catType) and 943 not fileName.endswith(self.sysc.emptyFitsCatalogueFileName)] 944 if not fileNames: 945 Logger.addMessage("<ERROR> Can't find any mosaics with " 946 "Multiframe.catName entries") 947 return False 948 949 isComplete = True 950 for fileName in fileNames: 951 mosMfID = self.archive.query("multiframeID", "Multiframe", 952 "frameType LIKE 'mosaicdeep%%stack' AND deprecated=0 AND " 953 "catName LIKE '%%%s'" % fileName, firstOnly=True) 954 if not mosMfID: 955 Logger.addMessage( 956 "<ERROR> Mosaic not in DB: " + fileName) 957 isComplete = False 958 959 numProv = self.archive.queryNumRows("Provenance", 960 "combiFrameID=%s" % mosMfID) 961 962 fitsFile = fits.open(fileName) 963 header = fitsFile[1].header.items() 964 keys = [key for key, _v in header if key.startswith('PROV')] 965 fitsFile.close() 966 967 if len(keys) - 1 != numProv: 968 Logger.addMessage( 969 "<ERROR> Mosaic has incomplete provenance: Expect %s " 970 "entries in Provenance, found %s, for multiframeID=%s, " 971 "catName=%s" % 972 (len(keys) - 1, numProv, mosMfID, fileName)) 973 isComplete = False 974 975 return isComplete 976 977 return numMosaics == numProvFrames
978 979 #-------------------------------------------------------------------------- 980
981 - def requiresUpdate(self, bitMask, frameType=None, isNonDepOnly=False):
982 """ 983 Does the current detection table require quality flagging? 984 985 @param bitMask: Bit mask for the quality issue to check (hexadecimal). 986 @type bitMask: str 987 @param frameType: Optionally only check frames of this type. 988 @type frameType: str 989 @param isNonDepOnly: If True, only check non-deprecated entries. 990 @type isNonDepOnly: bool 991 992 @return: True, if incomplete quality flagging for current detection 993 table. 994 @rtype: bool 995 996 """ 997 table = ("ProgrammeFrame" if not frameType and not isNonDepOnly else 998 Join(["ProgrammeFrame", "Multiframe"], "multiframeID")) 999 1000 whereStr = "ProgrammeFrame." + self._whereProgrammeObs(bitMask) 1001 if frameType: 1002 whereStr += " AND frameType LIKE '%s'" % frameType 1003 1004 if isNonDepOnly: 1005 whereStr += " AND " + DepCodes.selectNonDeprecated 1006 1007 return self.archive.queryEntriesExist(table, whereStr)
1008 1009 #-------------------------------------------------------------------------- 1010
1011 - def updateBadDetFlag(self, bitMask, condition, frameType="tile%stack"):
1012 """ 1013 Updates bit flag for source in a tile region coming from one or more 1014 bad detectors as defined by given condition. Flags sources with a 1015 default aperture that is at least partially in one of these regions. 1016 Also, considers positional offsets in forming deep tiles. 1017 1018 @param bitMask: Bit mask for this quality issue (hexadecimal). 1019 @type bitMask: str 1020 @param condition: Where clause that identifies detector as bad. 1021 @type condition: str 1022 @param frameType: Optionally only flag frames of this type. 1023 @type frameType: str 1024 1025 """ 1026 # First find frames that have condition 1027 safetyMargin = 2.0 1028 badTiles = self.archive.query( 1029 selectStr="Tile.multiframeID, TileAstro.xPixSize" 1030 ", Paw.multiframeID, Paw.extNum", 1031 fromStr="Multiframe AS Tile, CurrentAstrometry AS TileAstro" 1032 ", ProgrammeFrame, Provenance, MultiframeDetector AS Paw", 1033 whereStr="Tile.frameType LIKE '%s'" % frameType 1034 + " AND Tile." + self._whereProgrammeObs(bitMask) 1035 + " AND Tile.%s" % DepCodes.selectNonDeprecated 1036 + " AND TileAstro.multiframeID=Tile.multiframeID" 1037 " AND ProgrammeFrame.multiframeID=Tile.multiframeID" 1038 " AND Provenance.combiframeID=Tile.multiframeID" 1039 " AND Paw.multiframeID=Provenance.multiframeID" 1040 " AND Paw.%s" % condition) 1041 1042 numObjects = 0 1043 for tileMfID, tilePixScale, pawpMfID, pawpExtNum in badTiles: 1044 # Find corners of paw-print and convert to tile X and Y 1045 safMar = safetyMargin / tilePixScale 1046 skyPos = astro.findCorners(pawpMfID, pawpExtNum, self.archive) 1047 tileReg = astro.convSkyToPixel(skyPos, tileMfID, 2, self.archive, 1048 projection="TAN") 1049 1050 minX = min(tileReg, key=attrgetter('x')).x - safMar 1051 maxX = max(tileReg, key=attrgetter('x')).x + safMar 1052 minY = min(tileReg, key=attrgetter('y')).y - safMar 1053 maxY = max(tileReg, key=attrgetter('y')).y + safMar 1054 1055 # Flag all objects within tile region 1056 numObjects += self.archive.update( 1057 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1058 tableName=self._detTable, 1059 where="multiframeID=%s" % tileMfID 1060 + " AND x BETWEEN %s AND %s" % (minX, maxX) 1061 + " AND y BETWEEN %s AND %s" % (minY, maxY)) 1062 1063 # Don't want the log to be too big! 1064 self.archive.commitTransaction() 1065 1066 Logger.addMessage("%s sources affected." % numObjects) 1067 1068 self.updateStatus(bitMask, frameType=frameType, isNonDepOnly=True) 1069 self.archive.commitTransaction() # Don't want the log to be too big!
1070 1071 #-------------------------------------------------------------------------- 1072
1073 - def updateBadPixelFlag(self):
1074 """ Updates bit flag for bad pixel in source condition. 1075 """ 1076 bitNum = 6 1077 bitMask = self.bitMask(bitNum) 1078 if (not self.doBitFlag(bitNum, "Bad pixel in default aperture") 1079 or not self.requiresUpdate(bitMask)): 1080 return 1081 1082 Logger.addMessage("%s sources affected." % 1083 self.archive.update( 1084 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1085 tableName=self._detTable, 1086 fromTables=self._detTable + " AS det, ProgrammeFrame AS pf, " 1087 "Multiframe AS m", 1088 where="det.errBits>0 AND det.multiframeID=pf.multiframeID" 1089 " AND pf.programmeID = %s" % self.programmeID + 1090 " AND pf.ppErrBitsStatus & %s = 0" % bitMask + 1091 " AND m.multiframeID=det.multiframeID" 1092 " AND frameType NOT LIKE 'mosaicdeep%stack'")) 1093 1094 self.updateStatus(bitMask) 1095 self.archive.commitTransaction() # Don't want the log to be too big!
1096 1097 #-------------------------------------------------------------------------- 1098
1099 - def updateBoundaryFlag(self, bitMask, frameTypes, xWidth=None, yWidth=None):
1100 """ 1101 Updates bit flag for source close to boundary condition. Flags sources 1102 with a default aperture that is at least partially in the underexposed 1103 boundary region of a stack (either due to dithering = microstep + 1104 jitter, or for tiles also the "ears"). Also, considers positional 1105 offsets in forming deep stacks. 1106 1107 @param bitMask: Bit mask for this quality issue (hexadecimal). 1108 @type bitMask: str 1109 @param frameTypes: SQL string defining the frame types that this 1110 underexposed region applies to. 1111 @type frameTypes: str 1112 @param xWidth: SQL string defining the fully exposed width on the 1113 x-axis in pixels. 1114 @type xWidth: str 1115 @param yWidth: SQL string defining the fully exposed width on the 1116 y-axis in pixels. 1117 @type yWidth: str 1118 1119 @note: axis1/2Length is more reliable than x/ySize. 1120 1121 """ 1122 # Additional gap from boundary in arcseconds 1123 safetyMargin = 2.0 # i.e. twice default aperture radius 1124 1125 if xWidth: 1126 xWidth += " - %s / xPixSize" % safetyMargin 1127 else: 1128 xWidth = "axis1Length" 1129 1130 if yWidth: 1131 yWidth += " - %s / yPixSize" % safetyMargin 1132 else: 1133 yWidth = "axis2Length" 1134 1135 Logger.addMessage("%s sources affected." % 1136 self.archive.update( 1137 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1138 tableName=self._detTable, 1139 fromTables=Join([("Multiframe", "M"), ("MultiframeDetector", "MD"), 1140 ("CurrentAstrometry", "CA"), ("ProgrammeFrame", "PF"), 1141 (self._detTable, "D")], ["multiframeID"]), 1142 where=('' if self.sysc.isOSA() else "nuStep>0 AND ") + 1143 "(x < axis1Length - (%s) OR y < axis2Length - (%s) OR" 1144 " x > %s OR y > %s) AND " % (2 * (xWidth, yWidth)) + 1145 "MD.extNum=D.extNum AND CA.extNum=D.extNum AND " 1146 "seqNum > 0 AND programmeID = %s AND " 1147 "ppErrBitsStatus & %s = 0" % (self.programmeID, bitMask) 1148 + " AND " + frameTypes))
1149 1150 #-------------------------------------------------------------------------- 1151
1152 - def updateCrossTalkFlag(self, bitMask):
1153 """ 1154 Updates bit flag for possible cross-talk artefact. This method is quite 1155 unlike the other methods where the entire table is updated at once, 1156 here a list of multiframes are updated. 1157 1158 @param bitMask: Bit mask for this quality issue (hexadecimal). 1159 @type bitMask: str 1160 1161 """ 1162 CrossTalkFlagger(bitMask, self).run() 1163 if self._doMosaics: 1164 MosaicCrossTalkFlagger(bitMask, self).run()
1165 1166 #-------------------------------------------------------------------------- 1167
1168 - def updateDeblendedFlag(self):
1169 """ Updates bit flag for deblended source condition. 1170 """ 1171 bitNum = 4 1172 bitMask = self.bitMask(bitNum) 1173 if (not self.doBitFlag(bitNum, "Deblended sources") 1174 or not self.requiresUpdate(bitMask)): 1175 return 1176 1177 Logger.addMessage("%s sources affected." % 1178 self.archive.update( 1179 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1180 tableName=self._detTable, 1181 fromTables=self._detTable + " AS det, ProgrammeFrame AS pf" 1182 ", Multiframe AS m", 1183 where="det.aprof2=-1 AND det.multiframeID=pf.multiframeID" 1184 " AND pf.programmeID = %s" % self.programmeID + 1185 " AND pf.ppErrBitsStatus & %s = 0" % bitMask + 1186 " AND m.multiframeID=det.multiframeID" 1187 " AND frameType NOT LIKE 'mosaicdeep%stack'")) 1188 1189 self.updateStatus(bitMask) 1190 self.archive.commitTransaction() # Don't want the log to be too big!
1191 1192 #-------------------------------------------------------------------------- 1193
1194 - def updateDet16Flag(self):
1195 """ 1196 Updates bit flag for source in a region of a tile coming from 1197 detector 16. Flags sources with a default aperture that is at least 1198 partially in one of these regions. Also, considers positional offsets 1199 in forming deep tiles. 1200 1201 """ 1202 bitNum = 12 1203 bitMask = self.bitMask(bitNum) 1204 frType = "tile%stack" 1205 if (not self.sysc.isVSA() 1206 or not self.doBitFlag(bitNum, "Source in tile region of detector 16") 1207 or not self.requiresUpdate(bitMask, frType, isNonDepOnly=True)): 1208 return 1209 1210 self.updateBadDetFlag(bitMask, condition="extNum=17", frameType=frType)
1211 1212 #-------------------------------------------------------------------------- 1213
1214 - def updateDitherFlag(self):
1215 """ 1216 Updates bit flag for source in a dither region. Flags sources with a 1217 default aperture that is at least partially in the underexposed 1218 boundary region of a stack due to dithering = microstep + jitter. Also, 1219 considers positional offsets in forming deep stacks. 1220 1221 @todo: Boundary flag for UDS mosaic? Probably not worthwhile because 1222 there is only one field with detections consistent in each 1223 band, so users could just make their own x/y cuts. Alternatively 1224 implement a flag based on avConf. 1225 1226 """ 1227 bitNum = 22 1228 bitMask = self.bitMask(bitNum) 1229 if (not self.doBitFlag(bitNum, "Source close to boundary") 1230 or not self.requiresUpdate(bitMask)): 1231 return 1232 1233 nuStepMultiplier = ('' if self.sysc.isOSA() else "sqrt(nuStep) * ") 1234 self.updateBoundaryFlag(bitMask, 1235 frameTypes="frameType LIKE '%stack'" 1236 " AND frameType NOT LIKE '%mosaic%'" 1237 " AND frameType NOT LIKE '%tile%'", 1238 xWidth=nuStepMultiplier + "detRows", 1239 yWidth=nuStepMultiplier + "detCols") 1240 1241 if self.sysc.isVSA(): 1242 Logger.addMessage("...and for tiles...") 1243 self.updateBoundaryFlag(bitMask, 1244 frameTypes="frameType LIKE 'tile%stack'", 1245 xWidth="sqrt(nuStep) * %s" % self.sysc.xWidthTileDith, 1246 yWidth="sqrt(nuStep) * %s" % self.sysc.yWidthTileDith) 1247 1248 self.updateStatus(bitMask) 1249 self.archive.commitTransaction() # Don't want the log to be too big!
1250 1251 #-------------------------------------------------------------------------- 1252
1253 - def updateEarFlag(self):
1254 """ 1255 Updates bit flag for source in a ear region of a tile. Flags sources 1256 with a default aperture that is at least partially in the underexposed 1257 strip of a tile. Also, considers positional offsets in forming deep 1258 tiles. 1259 1260 """ 1261 bitNum = 23 1262 bitMask = self.bitMask(bitNum) 1263 if (not self.sysc.isVSA() 1264 or not self.doBitFlag(bitNum, "Source in tile underexposed strip") 1265 or not self.requiresUpdate(bitMask)): 1266 return 1267 1268 self.updateBoundaryFlag(bitMask, 1269 frameTypes="frameType LIKE 'tile%stack'", 1270 xWidth="sqrt(nuStep) * %s" % self.sysc.xWidthTileEar) 1271 1272 self.updateStatus(bitMask) 1273 self.archive.commitTransaction() # Don't want the log to be too big!
1274 1275 #-------------------------------------------------------------------------- 1276
1277 - def updateLowConfidenceFlag(self):
1278 """ 1279 Updates bit flag for low average confidence in source image condition. 1280 1281 """ 1282 bitNum = 7 1283 bitMask = self.bitMask(bitNum) 1284 if (self.sysc.isWSA() 1285 or not self.doBitFlag(bitNum, "Low confidence in default aperture") 1286 or not self.requiresUpdate(bitMask)): 1287 return 1288 1289 threshold = 90 if self.sysc.isOSA() else 80 1290 Logger.addMessage("%s sources affected." % 1291 self.archive.update( 1292 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1293 tableName=self._detTable, 1294 fromTables=self._detTable + " AS det, ProgrammeFrame AS pf, " 1295 "Multiframe AS m", 1296 where="det.averageConf BETWEEN 0 AND %s" % threshold + 1297 " AND det.multiframeID = pf.multiframeID" 1298 " AND pf.programmeID = %s" % self.programmeID + 1299 " AND pf.ppErrBitsStatus & %s = 0" % bitMask + 1300 " AND m.multiframeID = det.multiframeID" 1301 " AND frameType NOT LIKE 'mosaicdeep%stack'")) 1302 1303 self.updateStatus(bitMask) 1304 self.archive.commitTransaction() # Don't want the log to be too big!
1305 1306 #-------------------------------------------------------------------------- 1307
1308 - def updateMissingDetFlag(self):
1309 """ 1310 Updates bit flag for source in a region of a tile coming from one or 1311 more missing paw-print detectors. Flags sources with a default aperture 1312 that is at least partially in one of these regions. Also, considers 1313 positional offsets in forming deep tiles. 1314 1315 """ 1316 bitNum = 24 1317 bitMask = self.bitMask(bitNum) 1318 t = "tile%stack" 1319 msg = "Source in tile region of missing pawprint detector" 1320 if (not self.sysc.isVSA() or not self.doBitFlag(bitNum, msg) 1321 or not self.requiresUpdate(bitMask, t, isNonDepOnly=True)): 1322 return 1323 1324 self.updateBadDetFlag(bitMask, condition="totalExpTime=0", frameType=t)
1325 1326 #-------------------------------------------------------------------------- 1327
1328 - def updatePoorFlatFieldFlag(self):
1329 """ 1330 Updates bit flag for source in a region where flat fielding was poor. 1331 Flags sources with a default aperture that is at least partially in 1332 this region, accounting for stack dithering (= microstep + jitter) and 1333 positional offsets in forming deep stacks. 1334 1335 """ 1336 bitNum = 15 1337 bitMask = self.bitMask(bitNum) 1338 if (not self.sysc.isWSA() 1339 or not self.doBitFlag(bitNum, "Source in poor flat field region") 1340 or not self.requiresUpdate(bitMask)): 1341 return 1342 1343 self.updateBoundaryFlag(bitMask, 1344 frameTypes="D.extNum IN (2, 3) AND frameType NOT LIKE '%mosaic%'" 1345 " AND ppErrBits & 0x00400000 = 0", 1346 xWidth="SQRT(nuStep) * (detRows - 12.5)", 1347 yWidth="SQRT(nuStep) * (detCols - 12.5)") 1348 1349 Logger.addMessage("...and for the dodgy detector 2 corner...") 1350 1351 Logger.addMessage("%s sources affected." % 1352 self.archive.update( 1353 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1354 tableName=self._detTable, 1355 fromTables=Join([("Multiframe", "M"), 1356 ("CurrentAstrometry", "CA"), ("ProgrammeFrame", "PF"), 1357 (self._detTable, "D")], ["multiframeID"]), 1358 where="SQRT(SQUARE(x - 2270 * SQRT(nuStep))" 1359 " + SQUARE(y - 220 * SQRT(nuStep)))" 1360 " < 300 * SQRT(nuStep) AND CA.extNum=D.extNum AND " 1361 "frameType NOT LIKE '%deep%' AND " 1362 "frameType NOT LIKE '%mosaic%' AND " 1363 + "seqNum > 0 AND programmeID = %s AND D.extNum = 3 AND " 1364 "ppErrBitsStatus & %s = 0" % (self.programmeID, bitMask) 1365 )) 1366 1367 self.updateStatus(bitMask) 1368 self.archive.commitTransaction() # Don't want the log to be too big!
1369 1370 #-------------------------------------------------------------------------- 1371
1372 - def updateSaturatedFlag(self):
1373 """ Updates bit flag for source close to saturation condition. 1374 """ 1375 bitNum = 16 1376 bitMask = self.bitMask(bitNum) 1377 if (not self.doBitFlag(bitNum, "Saturated source image") 1378 or not self.requiresUpdate(bitMask)): 1379 return 1380 1381 # Number of ADU counts where source becomes saturated 1382 correctionLevel = "0.9*avSaturLevel" 1383 saturationLevels = (["40000", correctionLevel] if self.sysc.isWSA() else 1384 [correctionLevel]) 1385 1386 frameTypes = ([" AND frameType NOT LIKE '%deep%'", 1387 " AND frameType LIKE '%deep%' AND frameType NOT LIKE '%mosaic%'"] 1388 if self.sysc.isWSA() else [" AND frameType NOT LIKE '%tile%'" 1389 " AND frameType NOT LIKE '%mosaic%'"]) 1390 1391 for saturationLevel, frameType in zip(saturationLevels, frameTypes): 1392 whereStr = ("pHeight+sky > %s AND programmeID = %s AND " 1393 "ppErrBitsStatus & %s = 0" % (saturationLevel + frameType, 1394 self.programmeID, bitMask)) 1395 1396 joinTables = [self._detTable, "Multiframe", "ProgrammeFrame"] 1397 1398 if saturationLevel is correctionLevel: 1399 if self.sysc.isWSA(): 1400 Logger.addMessage("...and for deeps...") 1401 1402 joinTables.append("MultiframeDetector") 1403 whereStr += " AND %s.extNum = MultiframeDetector.extNum" \ 1404 % self._detTable 1405 1406 Logger.addMessage("%s sources affected." % 1407 self.archive.update( 1408 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1409 tableName=self._detTable, 1410 fromTables=Join(joinTables, ["multiframeID"]), 1411 where=whereStr)) 1412 1413 self.updateStatus(bitMask) 1414 self.archive.commitTransaction() # Don't want the log to be too big!
1415 1416 #-------------------------------------------------------------------------- 1417
1418 - def updateStatus(self, bitMask, mfID=None, frameType=None, 1419 isNonDepOnly=False):
1420 """ 1421 Updates ProgrammeFrame to denote quality flagging done. 1422 1423 @param bitMask: Bit mask of quality flag to mark as done (hexadecimal). 1424 @type bitMask: str 1425 @param mfID: Optionally just set this multiframeID as flagging done. 1426 @type mfID: int 1427 @param frameType: Optionally only flag frames of this type. 1428 @type frameType: str 1429 @param isNonDepOnly: If True, only flag non-deprecated entries. 1430 @type isNonDepOnly: bool 1431 1432 """ 1433 if mfID: 1434 whereStr = self._whereProgrammeObs(mfID=mfID) 1435 else: 1436 Logger.addMessage("Flagging multiframes as done...") 1437 whereStr = "ProgrammeFrame." + self._whereProgrammeObs(bitMask) 1438 1439 # @TODO: These two conditions should really be a part of 1440 # _whereProgrammeObs() and removed from the tile selection 1441 # in updateBadDetFlag() and from requiresUpdate(), but 1442 # requires alias support too! 1443 if frameType: 1444 whereStr += " AND frameType LIKE '%s'" % frameType 1445 1446 if isNonDepOnly: 1447 whereStr += " AND " + DepCodes.selectNonDeprecated 1448 1449 fromTables = ('' if mfID or not frameType and not isNonDepOnly else 1450 Join(["ProgrammeFrame", "Multiframe"], "multiframeID")) 1451 1452 numUpdated = self.archive.update("ProgrammeFrame", 1453 entryList=[("ppErrBitsStatus", "ppErrBitsStatus | " + bitMask)], 1454 where=whereStr, fromTables=fromTables) 1455 1456 if not mfID: 1457 Logger.addMessage("%s multiframes updated." % numUpdated)
1458 1459 #-------------------------------------------------------------------------- 1460
1461 - def updateUncorrectedFlag(self):
1462 """ Updates bit flag for source without illumination-corrected 1463 photometry condition. 1464 """ 1465 bitNum = 21 1466 bitMask = self.bitMask(bitNum) 1467 if (not self.sysc.isOSA() 1468 or not self.doBitFlag(bitNum, "Uncorrected sources") 1469 or not self.requiresUpdate(bitMask)): 1470 return 1471 1472 photTable = self._detTable.replace("Raw", "Photometry") 1473 Logger.addMessage("%s sources affected." % 1474 self.archive.update( 1475 entryList=[("ppErrBits", "ppErrBits | " + bitMask)], 1476 tableName=self._detTable, 1477 fromTables="%s AS det, %s AS phot, ProgrammeFrame AS pf" 1478 ", Multiframe AS m" % (self._detTable, photTable), 1479 where="illumCorr = 0 AND det.seqNum > 0" 1480 " AND det.multiframeID=phot.multiframeID AND " 1481 "det.extNum=phot.extNum AND det.seqNum=phot.seqNum" 1482 " AND det.multiframeID=pf.multiframeID" 1483 " AND pf.programmeID = %s" % self.programmeID + 1484 " AND pf.ppErrBitsStatus & %s = 0" % bitMask + 1485 " AND m.multiframeID=det.multiframeID")) 1486 1487 self.updateStatus(bitMask) 1488 self.archive.commitTransaction() # Don't want the log to be too big!
1489 1490 #-------------------------------------------------------------------------- 1491
1492 - def _whereProgrammeObs(self, bitMask=None, mfID=None):
1493 """ 1494 Selection of the unflagged (or given) observations from the current 1495 programme with existing entries in the current detection table 1496 1497 @param bitMask: Optionally select this bit as unflagged (hexadecimal). 1498 @type bitMask: str 1499 @param mfID: Optionally just select this multiframeID. 1500 @type mfID: int 1501 1502 @return: SQL WHERE string clause to specify just the unflagged 1503 (or given) observations from the current programme in the 1504 current detection table. 1505 @rtype: str 1506 1507 """ 1508 # @@GOTCHA: Need to include the selection on programmeID as well as 1509 # multiframeID to handle the case where a multiframe is shared 1510 # amongst two programmes with different detection tables. Also need to 1511 # handle the case whereby two different programmes share a single 1512 # detection table. 1513 progIDs = self.programme.getProgIDsfromTable(self._detTable) 1514 progSelection = "programmeID IN (%s)" % ', '.join(map(str, progIDs)) 1515 1516 if mfID: 1517 return "multiframeID = %s AND %s" % (mfID, progSelection) 1518 1519 mfSelection = \ 1520 SelectSQL("multiframeID", self._detTable, groupBy="multiframeID") 1521 1522 whereStr = "multiframeID IN (%s) AND %s" % (mfSelection, progSelection) 1523 if bitMask: 1524 whereStr += " AND ppErrBitsStatus & %s = 0" % bitMask 1525 1526 return whereStr
1527 1528 #------------------------------------------------------------------------------ 1529 # Change log: 1530 # 1531 # 20-Nov-2006, RSC: Original version. 1532 # 26-Nov-2008, RSC: Moved classes out of UpdateQualityFlags script and 1533 # combined with NonSurveyDeprecator from NonSurveyRelease 1534 # script to form a new QualityControl module. 1535