Merge pull request #70 from guaix-ucm/#8

Process Bad-pixels Mask recipe

megaradrp/core/recipe.py

@@ -42,11 +42,15 @@ def __init__(self, version):
                        'ArcCalibrationRecipe':[OverscanCorrector, TrimImage, BiasCorrector],
                        'FiberFlatRecipe':[OverscanCorrector, TrimImage, BiasCorrector],
                        'TraceMapRecipe':[OverscanCorrector, TrimImage, BiasCorrector],
+                       'BadPixelsMaskRecipe':[OverscanCorrector, TrimImage, BiasCorrector],
                        }
         super(MegaraBaseRecipe, self).__init__(version=version)
 
     def __generate_flow(self, params):
-        flow = self.__flow[self.__class__.__name__]
+        import copy
+        ff  = self.__flow[self.__class__.__name__]
+        flow = copy.deepcopy(ff)
+        # flow = self.__flow[self.__class__.__name__]
         try:
             for cont in range(len(flow)):
                 if issubclass(BiasCorrector, flow[cont]):
@@ -58,7 +62,7 @@ def __generate_flow(self, params):
         except Exception as e:
             _logger.error(e)
             raise(e)
-
+        del flow
         return basicflow
 
     def bias_process_common(self, obresult, master_bias):

megaradrp/drp.yaml

@@ -80,11 +80,21 @@ modes:
   summary: Summary of Arc Calibration
   uuid: c1842e00-2005-48bb-9a0e-8a8161f97fdf
   tagger: megaradrp.taggers.tagger_vph
+- date: 2012-05-21
+  description: Bad Pixels
+  key: bpm
+  name: Bad Pixel Mask
+  reference: xxx
+  status: DRAFT
+  summary: Bad Pixel Mask
+  uuid: xxxxxx
+  tagger: null
 pipelines:
   default:
     version: 1
     recipes:
       bias_image: megaradrp.recipes.calibration.bias.BiasRecipe
+      bpm: megaradrp.recipes.calibration.bpm.BadPixelsMaskRecipe
       dark_image: megaradrp.recipes.calibration.dark.DarkRecipe
       fiber_flat_image: megaradrp.recipes.calibration.flat.FiberFlatRecipe
       mos_image: megaradrp.recipes.scientific.FiberMOSRecipe2

megaradrp/recipes/calibration/bpm.py

@@ -19,30 +19,42 @@
 
 """Bad PIxel Mask (BPM) recipe"""
 
-import logging
-
-from numina.core import Product
+import astropy.io.fits as fits
+import numpy as np
+from numina.array.cosmetics import ccdmask
+from numina.core import Product, DataFrameType
 from numina.core.requirements import ObservationResultRequirement
 
 from megaradrp.core import MegaraBaseRecipe
-from megaradrp.products import MasterBias
-
-
-_logger = logging.getLogger('numina.recipes.megara')
+from megaradrp.requirements import MasterBiasRequirement
 
 
 class BadPixelsMaskRecipe(MegaraBaseRecipe):
-
-    '''Process BIAS images and create MASTER_BIAS.'''
-
     obresult = ObservationResultRequirement()
+    master_bias = MasterBiasRequirement()
 
-    biasframe = Product(MasterBias)
+    bpm_image = Product(DataFrameType)
 
     def __init__(self):
-        super(BadPixelsMaskRecipe, self).__init__(
-            version="0.1.0"
-        )
+        super(BadPixelsMaskRecipe, self).__init__(version="0.1.0")
 
     def run(self, rinput):
-        pass
+        import copy
+
+        N = len(rinput.obresult.frames)
+        obresult1 = copy.copy(rinput.obresult)
+        obresult1.frames = rinput.obresult.frames[:N//2]
+        obresult2 = copy.copy(rinput.obresult)
+        obresult2.frames = rinput.obresult.frames[N//2:]
+
+        reduced1 = self.bias_process_common(obresult1, rinput.master_bias)
+        reduced2 = self.bias_process_common(obresult2, rinput.master_bias)
+
+        mask = np.zeros(reduced1[0].data.shape, dtype='int')
+
+        bpm = ccdmask(reduced1[0].data, reduced2[0].data, mask, mode='full')
+        hdu = fits.PrimaryHDU(bpm)
+
+        reduced = fits.HDUList([hdu])
+
+        return self.create_result(bpm_image=reduced)

megaradrp/recipes/calibration/tests/test_bias.py

@@ -62,4 +62,7 @@ def test_bias(drpmocker):
 
     # In the end, remove the files
     for f in fs:
-        os.remove(f.name)
+        os.remove(f.name)
+
+if __name__ == "__main__":
+    test_bias()

megaradrp/recipes/calibration/tests/test_bpm.py

@@ -0,0 +1,120 @@
+#
+# Copyright 2015 Universidad Complutense de Madrid
+#
+# This file is part of Megara DRP
+#
+# Megara DRP 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.
+#
+# Megara DRP 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 GNU General Public License
+# along with Megara DRP.  If not, see <http://www.gnu.org/licenses/>.
+#
+
+"""Tests for the bpm mode recipe module."""
+import shutil
+from tempfile import mkdtemp
+
+import astropy.io.fits as fits
+import numpy as np
+
+from numina.core import DataFrame, ObservationResult
+
+from megaradrp.tests.simulation import simulate_flat, simulate_bias
+from megaradrp.tests.simulation import ReadParams, MegaraDetectorSat
+from megaradrp.recipes.calibration.bpm import BadPixelsMaskRecipe
+from megaradrp.recipes.calibration.bias import BiasRecipe
+
+
+def generate_bias(detector, number, temporary_path):
+    fs = [simulate_bias(detector) for i in range(number)]
+    for aux in range(len(fs)):
+        fits.writeto('%s/bias_%s.fits' % (temporary_path, aux), fs[aux],
+                     clobber=True)
+
+    fs = ["%s/bias_%s.fits" % (temporary_path, i) for i in range(number)]
+
+    ob = ObservationResult()
+    ob.instrument = 'MEGARA'
+    ob.mode = 'bias_image'
+    ob.frames = [DataFrame(filename=f) for f in fs]
+
+    recipe = BiasRecipe()
+    ri = recipe.create_input(obresult=ob)
+    return recipe.run(ri)
+
+
+def test_bpm():
+    number = 5
+    PSCAN = 50
+    DSHAPE = (2056 * 2, 2048 * 2)
+    OSCAN = 50
+
+    BINR = 1
+    BINC = 1
+
+    SHAPE = DSHAPE[0] // BINR, DSHAPE[1] // BINC
+
+    ron = 2.0
+    gain = 1.0
+    bias = 1000.0
+
+    eq = 0.8 * np.ones(DSHAPE)
+    eq[0:15, 0:170] = 0.0
+
+    temporary_path = mkdtemp()
+
+    fits.writeto('%s/eq.fits' % temporary_path, eq, clobber=True)
+
+    readpars1 = ReadParams(gain=gain, ron=ron, bias=bias)
+    readpars2 = ReadParams(gain=gain, ron=ron, bias=bias)
+
+    detector = MegaraDetectorSat(DSHAPE, OSCAN, PSCAN, eq=eq,
+                                 dark=(3.0 / 3600.0),
+                                 readpars1=readpars1, readpars2=readpars2,
+                                 bins='11')
+
+    source2 = 1.0
+
+    fs = [simulate_flat(detector, exposure=1.0, source=5000 * source2) for i in
+          range(number)]
+    fs2 = [simulate_flat(detector, exposure=1.0, source=40000 * source2) for i
+           in range(number)]
+
+    for aux in range(len(fs)):
+        fits.writeto('%s/flat_%s.fits' % (temporary_path, aux), fs[aux],
+                     clobber=True)
+        fits.writeto('%s/flat_%s.fits' % (temporary_path, aux + number),
+                     fs2[aux], clobber=True)
+
+    master_bias = generate_bias(detector, number, temporary_path)
+    master_bias_data = master_bias.biasframe.frame[0].data
+
+    fits.writeto('%s/master_bias_data0.fits' % temporary_path,
+                 master_bias_data, clobber=True)  # Master Bias
+
+    ob = ObservationResult()
+    ob.instrument = 'MEGARA'
+    ob.mode = 'bias_image'
+    names = []
+
+    for aux in range(number * 2):
+        names.append('%s/flat_%s.fits' % (temporary_path, aux))
+    ob.frames = [DataFrame(filename=open(nombre).name) for nombre in names]
+
+    recipe = BadPixelsMaskRecipe()
+    ri = recipe.create_input(obresult=ob, master_bias=DataFrame(
+        filename=open(temporary_path + '/master_bias_data0.fits').name))
+    recipe.run(ri)
+
+    shutil.rmtree(temporary_path)
+
+
+if __name__ == "__main__":
+    test_bpm()

megaradrp/tests/simulation.py

@@ -159,7 +159,7 @@ def f2(x):
 
         return numpy.piecewise(x, [p1, p2, p3], [f1, f2, f3])
 
-    @classmethod
+    # @classmethod
     def init_regions(cls, detshape, oscan, pscan, bng):
         """Create a image with overscan for testing."""
 
@@ -241,6 +241,21 @@ def create(self, mode, meta, data):
         return hdul
 
 
+class MegaraDetectorSat(MegaraDetector):
+
+    def saturate(self, x):
+        y = super(MegaraDetectorSat, self).saturate(x)
+
+        # Some pixels have a special nonlineary given by
+        # y[100,3000] = self.esp_nonlinearity(x[100,3000])
+        # y[3000,3000] = self.esp_nonlinearity(x[3000,3000])
+        return y
+
+    def esp_nonlinearity(self, x):
+        sat = 12000.0
+        return sat * (1-sat / (sat + x))
+
+
 def simulate_bias(detector):
     """Simulate a BIAS array."""
     detector.expose(source=0.0, time=0.0)
@@ -294,4 +309,4 @@ def simulate_flat_fits(factory, detector, exposure, source):
     data = final
 
     fitsfile = factory.create('flat', meta=meta, data=data)
-    return fitsfile
+    return fitsfile

megaradrp/tests/test_binning.py

@@ -1,7 +1,7 @@
 
 import numpy as np
 
-from .simulation import binning
+from megaradrp.tests.simulation import binning
 
 def test_binning():
     nr = 6
@@ -22,4 +22,8 @@ def test_binning():
 
     res = np.array([[ 18,  26,  34,  42], [ 82,  90,  98, 106], [146, 154, 162, 170]])
 
-    assert np.all(carr.sum(axis=-1) == res)
+    assert np.all(carr.sum(axis=-1) == res)
+
+
+if __name__ == "__main__":
+    test_binning()