Tickets/dm 5503

python/lsst/obs/cfht/cfhtIsrTask.py

@@ -1,25 +1,28 @@
+from __future__ import absolute_import, division, print_function
 import numpy as np
 
 import lsst.pex.config as pexConfig
 from lsst.ip.isr import IsrTask
 
-class CfhtIsrTaskConfig(IsrTask.ConfigClass) :
+
+class CfhtIsrTaskConfig(IsrTask.ConfigClass):
     safe = pexConfig.Field(
         dtype = float,
         doc = "Safety margin for CFHT sensors gain determination",
         default = 0.95,
     )
-    
+
     def setDefaults(self):
         IsrTask.ConfigClass.setDefaults(self)
 
-class CfhtIsrTask(IsrTask) :
+
+class CfhtIsrTask(IsrTask):
     ConfigClass = CfhtIsrTaskConfig
-    
+
     def run(self, ccdExposure, bias=None, linearizer=None, dark=None, flat=None, defects=None,
             fringes=None, bfKernel=None):
         """Perform instrument signature removal on an exposure
-        
+
         Steps include:
         - Detect saturation, apply overscan correction, bias, dark and flat
         - Perform CCD assembly
@@ -44,7 +47,7 @@ def run(self, ccdExposure, bias=None, linearizer=None, dark=None, flat=None, def
         ccd = ccdExposure.getDetector()
         floatExposure = self.convertIntToFloat(ccdExposure)
         metadata = floatExposure.getMetadata()
-        
+
         # Detect saturation
         # Saturation values recorded in the fits hader is not reliable, try to estimate it from
         # the pixel vales
@@ -53,13 +56,13 @@ def run(self, ccdExposure, bias=None, linearizer=None, dark=None, flat=None, def
         image = floatExposure.getMaskedImage().getImage()
         imageArray = image.getArray()
         maxValue = np.max(imageArray)
-        if maxValue > 60000.0 :
-            hist, bin_edges = np.histogram(imageArray.ravel(),bins=100,range=(60000.0,maxValue+1.0))
+        if maxValue > 60000.0:
+            hist, bin_edges = np.histogram(imageArray.ravel(), bins=100, range=(60000.0, maxValue+1.0))
             saturate = int(self.config.safe*bin_edges[np.argmax(hist)])
-        else :
+        else:
             saturate = metadata.get("SATURATE")
         self.log.info("Saturation set to %d" % saturate)
-        
+
         for amp in ccd:
             amp.setSaturation(saturate)
             if amp.getName() == "A":
@@ -68,7 +71,7 @@ def run(self, ccdExposure, bias=None, linearizer=None, dark=None, flat=None, def
                 # Check if the noise value is making sense for this amp. If not, replace with value
                 # stored in RDNOISE slot. This change is necessary to process some old CFHT images
                 # (visit : 7xxxxx) where RDNOISEA/B = 65535
-                if rdnA > 60000.0 :
+                if rdnA > 60000.0:
                     rdnA = metadata.get("RDNOISE")
                 amp.setReadNoise(rdnA)
             elif amp.getName() == "B":
@@ -77,18 +80,18 @@ def run(self, ccdExposure, bias=None, linearizer=None, dark=None, flat=None, def
                 # Check if the noise value is making sense for this amp. If not, replace with value
                 # stored in RDNOISE slot. This change is necessary to process some old CFHT images
                 # (visit : 7xxxxx) where RDNOISEA/B = 65535
-                if rdnB > 60000.0 :
+                if rdnB > 60000.0:
                     rdnB = metadata.get("RDNOISE")
                 amp.setReadNoise(rdnB)
-            else :
+            else:
                 raise ValueError("Unexpected amplifier name : %s"%(amp.getName()))
 
         return IsrTask.run(self,
-            ccdExposure = ccdExposure,
-            bias = bias,
-            linearizer = linearizer,
-            dark = dark,
-            flat = flat,
-            defects = defects,
-            fringes = fringes,
-        )
+                           ccdExposure = ccdExposure,
+                           bias = bias,
+                           linearizer = linearizer,
+                           dark = dark,
+                           flat = flat,
+                           defects = defects,
+                           fringes = fringes,
+                           )

python/lsst/obs/cfht/ingest.py

@@ -19,7 +19,7 @@
 # the GNU General Public License along with this program.  If not,
 # see <http://www.lsstcorp.org/LegalNotices/>.
 #
-
+from __future__ import absolute_import, division, print_function
 import re
 
 from lsst.pipe.tasks.ingest import ParseTask
@@ -32,11 +32,13 @@
            'z.MP9801': 'z',
            }
 
+
 class MegacamParseTask(ParseTask):
+
     def translate_ccd(self, md):
         try:
             extname = self.getExtensionName(md)
-            return int(extname[3:]) # chop off "ccd"
+            return int(extname[3:])  # chop off "ccd"
         except:
             # Dummy value, intended for PHU (need something to get filename)
             return 99
@@ -54,15 +56,16 @@ def translate_taiObs(self, md):
         (sec1, sec2) = sec.split('.')
         return "%04d-%02d-%02dT%02d:%02d:%02d.%02d"%(int(yr), int(month), int(day),
                                                      int(hr), int(min), int(sec1), int(sec2))
+
     def translate_defects(self, md):
         maskName = md.get("IMRED_MK").strip()
         maskName, ccd = maskName.split(".fits")
         filter = md.get("FILTER").strip().split('.')[0]
-        if filter == "i" or filter == "i2" or filter == "z" :
+        if filter == "i" or filter == "i2" or filter == "z":
             maskName = maskName+"_enlarged"
         maskFile = maskName+".nn/"+ccd[1:6]+".fits"
         return maskFile
-        
+
     def getInfo(self, filename):
         phuInfo, infoList = super(MegacamParseTask, self).getInfo(filename)
         match = re.search(r"\d+(?P<state>o|p)\.fits.*", filename)

python/lsst/obs/cfht/makeMegacamRawVisitInfo.py

@@ -0,0 +1,72 @@
+#
+# LSST Data Management System
+# Copyright 2016 LSST Corporation.
+#
+# This product includes software developed by the
+# LSST Project (http://www.lsst.org/).
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the LSST License Statement and
+# the GNU General Public License along with this program.  If not,
+# see <http://www.lsstcorp.org/LegalNotices/>.
+#
+from __future__ import absolute_import
+from __future__ import division
+from __future__ import print_function
+
+from lsst.afw.coord import Coord, IcrsCoord, Observatory, Weather
+from lsst.daf.butlerUtils import MakeRawVisitInfo
+
+__all__ = ["MakeMegacamRawVisitInfo"]
+
+
+class MakeMegacamRawVisitInfo(MakeRawVisitInfo):
+
+    """Make a VisitInfo from the FITS header of a raw Megacam image
+    """
+
+    def setArgDict(self, md, argDict):
+        """Set an argument dict for makeVisitInfo and pop associated metadata
+
+        @param[in,out] md  metadata, as an lsst.daf.base.PropertyList or PropertySet
+        @param[in,out] argdict  a dict of arguments
+        """
+        MakeRawVisitInfo.setArgDict(self, md, argDict)
+        argDict["darkTime"] = self.popFloat(md, "DARKTIME")
+        argDict["boresightAzAlt"] = Coord(
+            self.popAngle(md, "TELAZ"),
+            self.popAngle(md, "TELALT"),
+        )
+        argDict["boresightRaDec"] = IcrsCoord(
+            self.popAngle(md, "RA_DEG",),
+            self.popAngle(md, "DEC_DEG"),
+        )
+        argDict["boresightAirmass"] = self.popFloat(md, "AIRMASS")
+        argDict["observatory"] = Observatory(
+            self.popAngle(md, "LONGITUD"),
+            self.popAngle(md, "LATITUDE"),
+            4204,  # from Wikipedia
+        )
+        argDict["weather"] = Weather(
+            self.popFloat(md, "TEMPERAT"),
+            self.popFloat(md, "PRESSURE")*100.0,  # 100 Pascal per millibar
+            self.popFloat(md, "RELHUMID"),
+        )
+
+    def getDateAvg(self, md, exposureTime):
+        """Return date at the middle of the exposure
+
+        @param[in,out] md  FITS metadata; changed in place
+        @param[in] exposureTime  exposure time in sec
+        """
+        dateObs = self.popMjdDate(md, "MJD-OBS")
+        return self.offsetDate(dateObs, 0.5*exposureTime)

python/lsst/obs/cfht/megacamMapper.py

@@ -30,13 +30,17 @@
 
 from lsst.daf.butlerUtils import CameraMapper, exposureFromImage
 import lsst.pex.policy as pexPolicy
+from .makeMegacamRawVisitInfo import MakeMegacamRawVisitInfo
 
 # Solely to get boost serialization registrations for Measurement subclasses
 import lsst.meas.algorithms  # flake8: noqa
 
+
 class MegacamMapper(CameraMapper):
     packageName = "obs_cfht"
 
+    MakeRawVisitInfoClass = MakeMegacamRawVisitInfo
+
     def __init__(self, **kwargs):
         policyFile = pexPolicy.DefaultPolicyFile("obs_cfht", "MegacamMapper.paf", "policy")
         policy = pexPolicy.Policy(policyFile)
@@ -48,12 +52,12 @@ def __init__(self, **kwargs):
 
         self.exposures['raw'].keyDict['ccd'] = int
 
-        afwImageUtils.defineFilter('u',  lambdaEff=350, alias="u.MP9301")
-        afwImageUtils.defineFilter('g',  lambdaEff=450, alias="g.MP9401")
-        afwImageUtils.defineFilter('r',  lambdaEff=600, alias="r.MP9601")
-        afwImageUtils.defineFilter('i',  lambdaEff=750, alias="i.MP9701")
+        afwImageUtils.defineFilter('u', lambdaEff=350, alias="u.MP9301")
+        afwImageUtils.defineFilter('g', lambdaEff=450, alias="g.MP9401")
+        afwImageUtils.defineFilter('r', lambdaEff=600, alias="r.MP9601")
+        afwImageUtils.defineFilter('i', lambdaEff=750, alias="i.MP9701")
         afwImageUtils.defineFilter('i2', lambdaEff=750, alias="i.MP9702")
-        afwImageUtils.defineFilter('z',  lambdaEff=900, alias="z.MP9801")
+        afwImageUtils.defineFilter('z', lambdaEff=900, alias="z.MP9801")
 
         # define filters?
         self.filterIdMap = dict(u=0, g=1, r=2, i=3, z=4, i2=5)
@@ -100,8 +104,6 @@ def bypass_defects(self, datasetType, pythonType, butlerLocation, dataId):
 
         raise RuntimeError("No defects for ccdSerial %s in %s" % (ccdSerial, defectsFitsPath))
 
-
-
     def _defectLookup(self, dataId):
         """Find the defects for a given CCD.
         @param dataId (dict) Dataset identifier
@@ -110,8 +112,11 @@ def _defectLookup(self, dataId):
         if self.registry is None:
             raise RuntimeError, "No registry for defect lookup"
 
-        rows = self.registry.executeQuery(("defects",), ("raw",),
-                [("visit", "?"),("ccd", "?")], None, (dataId['visit'], dataId['ccd']))
+        rows = self.registry.executeQuery(
+            ("defects",),
+            ("raw",),
+            [("visit", "?"), ("ccd", "?")], None, (dataId['visit'], dataId['ccd']),
+        )
         if len(rows) == 0:
             return None
 
@@ -195,14 +200,15 @@ def bypass_stackExposureId(self, datasetType, pythonType, location, dataId):
 
     def bypass_stackExposureId_bits(self, datasetType, pythonType, location, dataId):
         """Hook to retrieve number of bits in identifier for stack/coadd"""
-        return 32 # not really, but this leaves plenty of space for sources
+        return 32  # not really, but this leaves plenty of space for sources
 
     def _standardizeDetrend(self, detrend, image, dataId, filter=False):
         """Hack up detrend images to remove troublesome keyword"""
         md = image.getMetadata()
-        removeKeyword(md, 'RADECSYS') # Irrelevant, and use of "GAPPT" breaks wcslib
+        removeKeyword(md, 'RADECSYS')  # Irrelevant, and use of "GAPPT" breaks wcslib
         exp = exposureFromImage(image)
-        return self._standardizeExposure(self.calibrations[detrend], exp, dataId, filter=filter, trimmed=False)
+        return self._standardizeExposure(self.calibrations[detrend], exp, dataId, filter=filter,
+                                         trimmed=False)
 
     def std_bias(self, image, dataId):
         return self._standardizeDetrend("bias", image, dataId, filter=False)