cr_detection_recipe.py

import pdb
import glob
import astrodata

# import dragon packages
from recipe_system import cal_service
from recipe_system.reduction.coreReduce import Reduce

from gempy.utils import logutils
from gempy.adlibrary import dataselect

import gemini_instruments
from recipe_system.mappers.primitiveMapper import PrimitiveMapper

logutils.config(file_name="test_reduce.log")


"""
Approach one, use coreReuce + cal_service + custom recipe

question: this way I could use the reduced bias and flats properly.
        But the reduced CR masks label CR pixels and saturated pxiels
        both as 8.
"""


all_files = glob.glob("./playdata/*.fits")

bias = dataselect.select_data(all_files, ["BIAS"])
flats = dataselect.select_data(all_files, ["FLAT"])

# this doesn't seem to be the right way to filter science...
sciences = dataselect.select_data(all_files, ["IMAGE", "RAW"])

# manual list for now
target = [
    "playdata/N20170614S0201.fits",
    "playdata/N20170614S0202.fits",
    "playdata/N20170614S0203.fits",
    "playdata/N20170614S0204.fits",
    "playdata/N20170614S0205.fits",
]

# init local calibration service
caldb = cal_service.CalibrationService()
caldb.config()
caldb.init()
cal_service.set_calservice()

# reduce and add bias calibration
reduce_bias = Reduce()
reduce_bias.files.extend(bias)
reduce_bias.runr()
caldb.add_cal(reduce_bias.output_filenames[0])

# reduce and add flat calibration
reduce_flats = Reduce()
reduce_flats.files.extend(flats)
reduce_flats.runr()
caldb.add_cal(reduce_flats.output_filenames[0])

# checkout local calibration database
for f in caldb.list_files():
    print(f)


"""
This is what I did but is the proper way to use a custom recipe?
in "recipes_IMAGE.py" I have --

def reduceCR(p):
    p.prepare()
    p.addDQ()
    p.addVAR(read_noise=True)
    p.overscanCorrect()
    p.biasCorrect()
    p.ADUToElectrons()
    p.addVAR(poisson_noise=True)
    p.flatCorrect()
    p.makeFringe()
    p.fringeCorrect()
    p.mosaicDetectors()
    p.detectSources()
    p.adjustWCSToReference()
    p.flagCosmicRaysByStacking()
    p.writeOutputs()
    return

"""


reduce_target = Reduce()
reduce_target.recipename = "./recipes_IMAGE.reduceCR"

vars(reduce_target).items()

reduce_target.files.extend(target)
# reduce_target.uparms.append(("skyCorrect:scale", False))
reduce_target.runr()

reduce_target.output_files


"""
Approach two, use primitives directly

question: I'm not sure how to utilize the cal_serviec this way. 
"""

# ad = astrodata.open("playdata/N20170614S0201.fits")
# ad1 = astrodata.open("playdata/N20170614S0202.fits")
# ad2 = astrodata.open("playdata/N20170614S0203.fits")

# pmapper = PrimitiveMapper([ad, ad1, ad2])
# pset = pmapper.get_applicable_primitives()

# # save out the CR aks
# pset.prepare()
# pset.addDQ()
# pset.addVAR(read_noise=True)
# pset.overscanCorrect()
# pset.biasCorrect()
# pset.ADUToElectrons()
# pset.addVAR(poisson_noise=True)
# pset.flatCorrect()
# pset.mosaicDetectors()
# pset.adjustWCSToReference()
# pset.flagCosmicRaysByStacking()
# pset.writeOutputs()

# # brekapoint
# pdb.set_trace()