From 4669b09763c3cc35b1a48d01cb2c0585fce79a8c Mon Sep 17 00:00:00 2001 From: Steven Crawford Date: Tue, 5 Mar 2013 13:06:34 +0200 Subject: [PATCH 01/24] adding ccdproc --- ccdproc_api.py | 200 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 200 insertions(+) create mode 100644 ccdproc_api.py diff --git a/ccdproc_api.py b/ccdproc_api.py new file mode 100644 index 0000000..8e9d961 --- /dev/null +++ b/ccdproc_api.py @@ -0,0 +1,200 @@ +#!/usr/bin/python +# -*- coding: utf-8 -*- + +from __future__ import print_function + +from astropy.nddata import NDdata +''' +The ccdproc package provides tools for the reduction and +analysis of optical data captured with a CCD. The package +is built around the CCDData class, which has built into +it all of the functions to process data. The CCDData object +contains all the information to describe the 2-D readout +from a single amplifier/detector. + +The CCDData class inherits from the NDData class as its base object +and the object on which all actions will be performed. By +inheriting from the CCD data class, basic array manipulation +and error handling are already built into the object. + +The CCDData task should be able to properly deal with the +propogation of errors and propogate bad pixel frames +through each of the tasks. It should also update the meta +data, units, and WCS information as the data are processed +through each step. + +The following functions are required for performing basic CCD correction: +-creation of variance frame +-overscan subtraction +-bias subtraction +-trimming the data +-gain correction +-xtalk correction +-dark frames correction +-flat field correction +-illumination correction +-fringe correction +-scattered light correction +-cosmic ray cleaning +-distortion correction + +In addition to the CCDData and CCDList class, the ccdproc does +require some additional features in order to properly +reduce CCD data. The following features are required +for basic processing of CCD data: +-fitting data +-combining data +-re-sampling data +-transforming data + +All actions of ccdproc should be logged and recorded. + +Multi-Extension FITS files can be handled by treating +each extension as a CCDData object and + +''' + +# ============ +# Base Objects +# ============ +''' +CCDData is an object that inherits from NDData class and specifically +describes an object created from the single readout of a CCD. + +Users should be able to create a CCDData object from scratch, from +an existing NDData object, or a single extension from a FITS file. + +In the case of the CCDData, the parameter 'uncertainty' will +be mapped to variance as that will be more explicitly describing +the information that will be kept for the processing of the + +''' +data=100+10*np.random.random((110,100)) +ccddata=CCDData(data=data) +ccddata=CCDData(NDData.NDData) +ccddata=CCDData(pyfits.ImageHDU) + +#Setting basic properties of the object +# ---------------------- +ccddata.variance=data**0.5 +ccddata.mask=np.ones(110,100) +ccddata.flags=np.zeros(110,100) +ccddata.wcs=None +ccddata.meta={} +ccddata.units=u.adu #is this valid? + +# Functional Requirements +# ---------------------- +# A number of these different fucntions are convenient functions that +# just outline the process that is needed. The important thing is that +# the process is being logged and that a clear process is being handled +# by each step to make building a pipeline easy. Then again, it might +# not be easy to handle all possible steps which are needed, and the more +# important steps will be things which aren't already handled by NDData. + +#All functions should propogate throught to the variance frame and +#bad pixel mask + +#convenience function based on a given value for +#the readnoise and gain. Units should be specified +#for these values though. +#Question: Do we have an object that is just a scalar +#and a unit? Or a scalar, unit and an error? ie, This +#could actually be handled by the gain and readnoise being +#specified as an NDData object +ccddata.createvariance(gain=1.0, readnoise=5.0) + +#Overscan subtract the data +#Should be able to provide the meta data for +#the keyworkd or provide a section to define and +#possible an axis to specify the oritation of the +#Question: Best way to specify the section? Should it be given +#Error Checks: That the section is within the image +ccddata.subtract_overscan(section='[:,100:110]', function='polynomial', order=3) + +#trim the images--the section gives the part of the image to keep +#That the trim section is within the image +ccddata.trim_image(section='[0:100,0:100]') + +#subtract the master bias. Although this is a convenience function as subtracting +#the two arrays will do the same thing. This should be able to handle logging of +#of subtracting it off (logging should be added to NDData and then this is really +#just a convenience function +#Error checks: the masterbias and image are the same shape +masterbias=NDData.NDData(np.zeros(100,100)) +ccddata.subtract_bias(masterbias) + +#correct for dark frames +#Error checks: the masterbias and image are the same shape +masterdark=NDData.NDData(np.zeros(100,100)) +ccddata.subtract_dark(darkframe) + +#correct for gain--once again gain should have a unit and even an error associated with it. +ccddata.gain_correct(gain=1.0) +#Also the gain may be non-linear +ccddata.gain_correct(gain=np.array([1.0,0.5e-3]) +#although then this step should be apply before any other corrections if it is non-linear +#but that is more up to the person processing their own data. + +#crosstalk corrections--also potential a convenience function, but basically multiples the +#xtalkimage by the coeffient and then subtracts it. It is kept general because this will +#be dependent on the CCD and the set up of the CCD. Not applicable for a single CCD +#situation +#Error checks: the xtalkimage and image are the same shape +xtalkimage=NDData.NDData(np.zeros(100,100)) +ccddata.xtalk_correct(xtalkimage, coef=1e-3) + +#flat field correction--this can either be a dome flat, sky flat, or an +#illumination corrected image. This step should normalize by the value of the +#flatfield after dividing by it. +#Error checks: the flatimage and image are the same shape +#Error checks: check for divive by zero +#Features: If the flat is less than minvalue, minvalue is used +flatimage=NDData.NDData(np.ones(100,100)) +ccddata.flat_correct(flatimage, minvalue=1) + +#fringe correction or any correction that requires subtracting +#off a potentially scaled image +#Error checks: the flatimage and image are the same shape +fringemage=NDData.NDData(np.ones(100,100)) +ccddata.fringe_correct(fringeimage, scale=1) + +#cosmic ray cleaning step--this should have options for different +#ways to do it with their associated steps. We also might want to +#implement this as a slightly different step. The cosmic ray cleaning +#step should update the mask and flags. So the user could have options +#to replace the cosmic rays, only flag the cosmic rays, or flag and +#mask the cosmic rays, or all of the above. +ccddata.cosmicray_clean(method='laplace', args=*kwargs) + +#Apply distortion corrections +#Either update the WCS or transform the frame +ccddata.distortion_correct(distortion) + + +# ================ +# Helper Functions +# ================ + +#fit a 1-D function with iterative rejections and the ability to +#select different functions to fit. +#other options are reject parameters, number of iteractions +#and/or convergernce limit +coef=iterfit(x, y, function='polynomial', order=3) + +#fit a 2-D function with iterative rejections and the ability to +#select different functions to fit. +#other options are reject parameters, number of iteractions +#and/or convergernce limit +coef=iterfit(data, function='polynomial', order=3) + +#combine a set of NDData objects +combine([ccddata, ccddata2], method='average', reject=None, **kwargs) + +#re-sample the data to different binnings (either larger or smaller) +ccddata=rebin(ccddata, binning=(2,2)) + +#tranform the data--ie shift, rotate, etc +#question--add convenience functions for image shifting and rotation? +#should udpate WCS although that would actually be the prefered method +ccddata=transform(ccddata, transform, conserve_flux=True) From 8ecb546cdd83f484e50e45893ce443b093af7e9c Mon Sep 17 00:00:00 2001 From: Steven Crawford Date: Mon, 7 Oct 2013 17:13:25 +0200 Subject: [PATCH 02/24] update ccdproc api --- ccdproc_api.py | 37 ++++++++++++++++++++++++------------- 1 file changed, 24 insertions(+), 13 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 8e9d961..3a899ea 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -83,6 +83,10 @@ ccddata.meta={} ccddata.units=u.adu #is this valid? +#The ccddata class should have a functional form to create a CCDData +#object directory from a fits file +ccddata=fromFITS('img.fits') + # Functional Requirements # ---------------------- # A number of these different fucntions are convenient functions that @@ -102,7 +106,7 @@ #and a unit? Or a scalar, unit and an error? ie, This #could actually be handled by the gain and readnoise being #specified as an NDData object -ccddata.createvariance(gain=1.0, readnoise=5.0) +ccddata=createvariance(ccddata, gain=1.0, readnoise=5.0) #Overscan subtract the data #Should be able to provide the meta data for @@ -110,11 +114,11 @@ #possible an axis to specify the oritation of the #Question: Best way to specify the section? Should it be given #Error Checks: That the section is within the image -ccddata.subtract_overscan(section='[:,100:110]', function='polynomial', order=3) +ccddata=subtract_overscan(ccddata, section='[:,100:110]', function='polynomial', order=3) #trim the images--the section gives the part of the image to keep #That the trim section is within the image -ccddata.trim_image(section='[0:100,0:100]') +ccddata=trim_image(ccddata, section='[0:100,0:100]') #subtract the master bias. Although this is a convenience function as subtracting #the two arrays will do the same thing. This should be able to handle logging of @@ -122,17 +126,17 @@ #just a convenience function #Error checks: the masterbias and image are the same shape masterbias=NDData.NDData(np.zeros(100,100)) -ccddata.subtract_bias(masterbias) +ccddata=subtract_bias(ccddata, masterbias) #correct for dark frames #Error checks: the masterbias and image are the same shape masterdark=NDData.NDData(np.zeros(100,100)) -ccddata.subtract_dark(darkframe) +ccddata=subtract_dark(ccddata,darkframe) #correct for gain--once again gain should have a unit and even an error associated with it. -ccddata.gain_correct(gain=1.0) +ccddata=gain_correct(ccddata, gain=1.0) #Also the gain may be non-linear -ccddata.gain_correct(gain=np.array([1.0,0.5e-3]) +ccddata=gain_correct(ccddata, gain=np.array([1.0,0.5e-3]) #although then this step should be apply before any other corrections if it is non-linear #but that is more up to the person processing their own data. @@ -142,7 +146,7 @@ #situation #Error checks: the xtalkimage and image are the same shape xtalkimage=NDData.NDData(np.zeros(100,100)) -ccddata.xtalk_correct(xtalkimage, coef=1e-3) +ccddata=xtalk_correct(ccddata, xtalkimage, coef=1e-3) #flat field correction--this can either be a dome flat, sky flat, or an #illumination corrected image. This step should normalize by the value of the @@ -151,13 +155,13 @@ #Error checks: check for divive by zero #Features: If the flat is less than minvalue, minvalue is used flatimage=NDData.NDData(np.ones(100,100)) -ccddata.flat_correct(flatimage, minvalue=1) +ccddata=flat_correct(ccddata, flatimage, minvalue=1) #fringe correction or any correction that requires subtracting #off a potentially scaled image #Error checks: the flatimage and image are the same shape fringemage=NDData.NDData(np.ones(100,100)) -ccddata.fringe_correct(fringeimage, scale=1) +ccddata=fringe_correct(ccddata, fringeimage, scale=1, operation='multiple') #cosmic ray cleaning step--this should have options for different #ways to do it with their associated steps. We also might want to @@ -165,11 +169,12 @@ #step should update the mask and flags. So the user could have options #to replace the cosmic rays, only flag the cosmic rays, or flag and #mask the cosmic rays, or all of the above. -ccddata.cosmicray_clean(method='laplace', args=*kwargs) +ccddata=cosmicray_laplace(ccddata, method='laplace', args=*kwargs) +ccddata=cosmicray_median(ccddata, method='laplace', args=*kwargs) #Apply distortion corrections #Either update the WCS or transform the frame -ccddata.distortion_correct(distortion) +ccddata=distortion_correct(ccddata, distortion) # ================ @@ -188,8 +193,14 @@ #and/or convergernce limit coef=iterfit(data, function='polynomial', order=3) +#in addition to these operations, basic addition, subtraction +# multiplication, and division should work for CCDDATA objects +ccddata= ccdata + ccddata +ccddata= ccddata * 2 + + #combine a set of NDData objects -combine([ccddata, ccddata2], method='average', reject=None, **kwargs) +alldata=combine([ccddata, ccddata2], method='average', reject=None, **kwargs) #re-sample the data to different binnings (either larger or smaller) ccddata=rebin(ccddata, binning=(2,2)) From eff09044fa1fa1d2907a46231e0861711f3dd30d Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 31 Jan 2014 16:35:11 -0600 Subject: [PATCH 03/24] Whitespace/formatting changes for PEP8 compliance --- ccdproc_api.py | 157 +++++++++++++++++++++++++------------------------ 1 file changed, 80 insertions(+), 77 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 3a899ea..d0e26c6 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -5,12 +5,12 @@ from astropy.nddata import NDdata ''' -The ccdproc package provides tools for the reduction and +The ccdproc package provides tools for the reduction and analysis of optical data captured with a CCD. The package is built around the CCDData class, which has built into -it all of the functions to process data. The CCDData object +it all of the functions to process data. The CCDData object contains all the information to describe the 2-D readout -from a single amplifier/detector. +from a single amplifier/detector. The CCDData class inherits from the NDData class as its base object and the object on which all actions will be performed. By @@ -26,7 +26,7 @@ The following functions are required for performing basic CCD correction: -creation of variance frame -overscan subtraction --bias subtraction +-bias subtraction -trimming the data -gain correction -xtalk correction @@ -36,10 +36,10 @@ -fringe correction -scattered light correction -cosmic ray cleaning --distortion correction +-distortion correction -In addition to the CCDData and CCDList class, the ccdproc does -require some additional features in order to properly +In addition to the CCDData and CCDList class, the ccdproc does +require some additional features in order to properly reduce CCD data. The following features are required for basic processing of CCD data: -fitting data @@ -49,8 +49,8 @@ All actions of ccdproc should be logged and recorded. -Multi-Extension FITS files can be handled by treating -each extension as a CCDData object and +Multi-Extension FITS files can be handled by treating +each extension as a CCDData object and ''' @@ -59,153 +59,156 @@ # ============ ''' CCDData is an object that inherits from NDData class and specifically -describes an object created from the single readout of a CCD. +describes an object created from the single readout of a CCD. Users should be able to create a CCDData object from scratch, from -an existing NDData object, or a single extension from a FITS file. +an existing NDData object, or a single extension from a FITS file. In the case of the CCDData, the parameter 'uncertainty' will be mapped to variance as that will be more explicitly describing -the information that will be kept for the processing of the +the information that will be kept for the processing of the ''' -data=100+10*np.random.random((110,100)) -ccddata=CCDData(data=data) -ccddata=CCDData(NDData.NDData) -ccddata=CCDData(pyfits.ImageHDU) +data = 100 + 10 * np.random.random((110, 100)) +ccddata = CCDData(data=data) +ccddata = CCDData(NDData.NDData) +ccddata = CCDData(pyfits.ImageHDU) #Setting basic properties of the object # ---------------------- -ccddata.variance=data**0.5 -ccddata.mask=np.ones(110,100) -ccddata.flags=np.zeros(110,100) -ccddata.wcs=None -ccddata.meta={} -ccddata.units=u.adu #is this valid? +ccddata.variance = data**0.5 +ccddata.mask = np.ones(110, 100) +ccddata.flags = np.zeros(110, 100) +ccddata.wcs = None +ccddata.meta = {} +ccddata.units = u.adu # is this valid? #The ccddata class should have a functional form to create a CCDData #object directory from a fits file -ccddata=fromFITS('img.fits') +ccddata = fromFITS('img.fits') # Functional Requirements # ---------------------- # A number of these different fucntions are convenient functions that # just outline the process that is needed. The important thing is that -# the process is being logged and that a clear process is being handled +# the process is being logged and that a clear process is being handled # by each step to make building a pipeline easy. Then again, it might # not be easy to handle all possible steps which are needed, and the more # important steps will be things which aren't already handled by NDData. -#All functions should propogate throught to the variance frame and +#All functions should propogate throught to the variance frame and #bad pixel mask -#convenience function based on a given value for +#convenience function based on a given value for #the readnoise and gain. Units should be specified #for these values though. #Question: Do we have an object that is just a scalar -#and a unit? Or a scalar, unit and an error? ie, This +#and a unit? Or a scalar, unit and an error? ie, This #could actually be handled by the gain and readnoise being -#specified as an NDData object -ccddata=createvariance(ccddata, gain=1.0, readnoise=5.0) +#specified as an NDData object +ccddata = createvariance(ccddata, gain=1.0, readnoise=5.0) #Overscan subtract the data #Should be able to provide the meta data for #the keyworkd or provide a section to define and -#possible an axis to specify the oritation of the -#Question: Best way to specify the section? Should it be given +#possible an axis to specify the oritation of the +#Question: Best way to specify the section? Should it be given #Error Checks: That the section is within the image -ccddata=subtract_overscan(ccddata, section='[:,100:110]', function='polynomial', order=3) +ccddata = subtract_overscan(ccddata, section='[:,100:110]', + function='polynomial', order=3) #trim the images--the section gives the part of the image to keep #That the trim section is within the image -ccddata=trim_image(ccddata, section='[0:100,0:100]') +ccddata = trim_image(ccddata, section='[0:100,0:100]') -#subtract the master bias. Although this is a convenience function as subtracting -#the two arrays will do the same thing. This should be able to handle logging of -#of subtracting it off (logging should be added to NDData and then this is really -#just a convenience function +#subtract the master bias. Although this is a convenience function as +#subtracting the two arrays will do the same thing. This should be able +#to handle logging of subtracting it off (logging should be added to NDData +#and then this is really just a convenience function #Error checks: the masterbias and image are the same shape -masterbias=NDData.NDData(np.zeros(100,100)) -ccddata=subtract_bias(ccddata, masterbias) +masterbias = NDData.NDData(np.zeros(100, 100)) +ccddata = subtract_bias(ccddata, masterbias) #correct for dark frames #Error checks: the masterbias and image are the same shape -masterdark=NDData.NDData(np.zeros(100,100)) -ccddata=subtract_dark(ccddata,darkframe) +masterdark = NDData.NDData(np.zeros(100, 100)) +ccddata = subtract_dark(ccddata, darkframe) -#correct for gain--once again gain should have a unit and even an error associated with it. -ccddata=gain_correct(ccddata, gain=1.0) +#correct for gain--once again gain should have a unit and even an error +#associated with it. +ccddata = gain_correct(ccddata, gain=1.0) #Also the gain may be non-linear -ccddata=gain_correct(ccddata, gain=np.array([1.0,0.5e-3]) -#although then this step should be apply before any other corrections if it is non-linear -#but that is more up to the person processing their own data. - -#crosstalk corrections--also potential a convenience function, but basically multiples the -#xtalkimage by the coeffient and then subtracts it. It is kept general because this will -#be dependent on the CCD and the set up of the CCD. Not applicable for a single CCD -#situation +ccddata = gain_correct(ccddata, gain=np.array([1.0, 0.5e-3])) +#although then this step should be apply before any other corrections +#if it is non-linear, but that is more up to the person processing their +#own data. + +#crosstalk corrections--also potential a convenience function, but basically +#multiples the xtalkimage by the coeffient and then subtracts it. It is kept +#general because this will be dependent on the CCD and the set up of the CCD. +#Not applicable for a single CCD situation #Error checks: the xtalkimage and image are the same shape -xtalkimage=NDData.NDData(np.zeros(100,100)) -ccddata=xtalk_correct(ccddata, xtalkimage, coef=1e-3) +xtalkimage = NDData.NDData(np.zeros(100, 100)) +ccddata = xtalk_correct(ccddata, xtalkimage, coef=1e-3) -#flat field correction--this can either be a dome flat, sky flat, or an +#flat field correction--this can either be a dome flat, sky flat, or an #illumination corrected image. This step should normalize by the value of the #flatfield after dividing by it. #Error checks: the flatimage and image are the same shape #Error checks: check for divive by zero #Features: If the flat is less than minvalue, minvalue is used -flatimage=NDData.NDData(np.ones(100,100)) -ccddata=flat_correct(ccddata, flatimage, minvalue=1) +flatimage = NDData.NDData(np.ones(100, 100)) +ccddata = flat_correct(ccddata, flatimage, minvalue=1) #fringe correction or any correction that requires subtracting #off a potentially scaled image #Error checks: the flatimage and image are the same shape -fringemage=NDData.NDData(np.ones(100,100)) -ccddata=fringe_correct(ccddata, fringeimage, scale=1, operation='multiple') +fringemage = NDData.NDData(np.ones(100, 100)) +ccddata = fringe_correct(ccddata, fringeimage, scale=1, operation='multiple') #cosmic ray cleaning step--this should have options for different -#ways to do it with their associated steps. We also might want to +#ways to do it with their associated steps. We also might want to #implement this as a slightly different step. The cosmic ray cleaning #step should update the mask and flags. So the user could have options -#to replace the cosmic rays, only flag the cosmic rays, or flag and +#to replace the cosmic rays, only flag the cosmic rays, or flag and #mask the cosmic rays, or all of the above. -ccddata=cosmicray_laplace(ccddata, method='laplace', args=*kwargs) -ccddata=cosmicray_median(ccddata, method='laplace', args=*kwargs) +ccddata = cosmicray_laplace(ccddata, method='laplace', **kwargs) +ccddata = cosmicray_median(ccddata, method='laplace', **kwargs) #Apply distortion corrections #Either update the WCS or transform the frame -ccddata=distortion_correct(ccddata, distortion) +ccddata = distortion_correct(ccddata, distortion) # ================ -# Helper Functions +# Helper Functions # ================ -#fit a 1-D function with iterative rejections and the ability to -#select different functions to fit. -#other options are reject parameters, number of iteractions +#fit a 1-D function with iterative rejections and the ability to +#select different functions to fit. +#other options are reject parameters, number of iteractions #and/or convergernce limit -coef=iterfit(x, y, function='polynomial', order=3) +coef = iterfit(x, y, function='polynomial', order=3) -#fit a 2-D function with iterative rejections and the ability to -#select different functions to fit. -#other options are reject parameters, number of iteractions +#fit a 2-D function with iterative rejections and the ability to +#select different functions to fit. +#other options are reject parameters, number of iteractions #and/or convergernce limit -coef=iterfit(data, function='polynomial', order=3) +coef = iterfit(data, function='polynomial', order=3) #in addition to these operations, basic addition, subtraction # multiplication, and division should work for CCDDATA objects -ccddata= ccdata + ccddata -ccddata= ccddata * 2 +ccddata = ccddata + ccddata +ccddata = ccddata * 2 #combine a set of NDData objects -alldata=combine([ccddata, ccddata2], method='average', reject=None, **kwargs) +alldata = combine([ccddata, ccddata2], method='average', reject=None, **kwargs) #re-sample the data to different binnings (either larger or smaller) -ccddata=rebin(ccddata, binning=(2,2)) +ccddata = rebin(ccddata, binning=(2, 2)) #tranform the data--ie shift, rotate, etc #question--add convenience functions for image shifting and rotation? #should udpate WCS although that would actually be the prefered method -ccddata=transform(ccddata, transform, conserve_flux=True) +ccddata = transform(ccddata, transform, conserve_flux=True) From b74352953505e1e5a207900f1b4d14d633d7656d Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 31 Jan 2014 17:48:14 -0600 Subject: [PATCH 04/24] Fix imports/function namespace so that API is closer to ready to execute There are still a few undefined variables (transform, distortion, and x/y in the fit example). --- ccdproc_api.py | 60 ++++++++++++++++++++++++++++---------------------- 1 file changed, 34 insertions(+), 26 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index d0e26c6..5686b10 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -3,7 +3,13 @@ from __future__ import print_function -from astropy.nddata import NDdata +import numpy as np + +from astropy.nddata import NDData +from astropy.io import fits +from astropy import units as u + +import ccdproc ''' The ccdproc package provides tools for the reduction and analysis of optical data captured with a CCD. The package @@ -70,9 +76,9 @@ ''' data = 100 + 10 * np.random.random((110, 100)) -ccddata = CCDData(data=data) -ccddata = CCDData(NDData.NDData) -ccddata = CCDData(pyfits.ImageHDU) +ccddata = ccdproc.CCDData(data=data) +ccddata = ccdproc.CCDData(NDData.NDData(data)) +ccddata = ccdproc.CCDData(fits.ImageHDU(data)) #Setting basic properties of the object # ---------------------- @@ -85,7 +91,7 @@ #The ccddata class should have a functional form to create a CCDData #object directory from a fits file -ccddata = fromFITS('img.fits') +ccddata = ccdproc.CCDData.fromFITS('img.fits') # Functional Requirements # ---------------------- @@ -106,7 +112,7 @@ #and a unit? Or a scalar, unit and an error? ie, This #could actually be handled by the gain and readnoise being #specified as an NDData object -ccddata = createvariance(ccddata, gain=1.0, readnoise=5.0) +ccddata = ccdproc.createvariance(ccddata, gain=1.0, readnoise=5.0) #Overscan subtract the data #Should be able to provide the meta data for @@ -114,12 +120,12 @@ #possible an axis to specify the oritation of the #Question: Best way to specify the section? Should it be given #Error Checks: That the section is within the image -ccddata = subtract_overscan(ccddata, section='[:,100:110]', - function='polynomial', order=3) +ccddata = ccdproc.subtract_overscan(ccddata, section='[:,100:110]', + function='polynomial', order=3) #trim the images--the section gives the part of the image to keep #That the trim section is within the image -ccddata = trim_image(ccddata, section='[0:100,0:100]') +ccddata = ccdproc.trim_image(ccddata, section='[0:100,0:100]') #subtract the master bias. Although this is a convenience function as #subtracting the two arrays will do the same thing. This should be able @@ -127,18 +133,18 @@ #and then this is really just a convenience function #Error checks: the masterbias and image are the same shape masterbias = NDData.NDData(np.zeros(100, 100)) -ccddata = subtract_bias(ccddata, masterbias) +ccddata = ccdproc.subtract_bias(ccddata, masterbias) #correct for dark frames #Error checks: the masterbias and image are the same shape masterdark = NDData.NDData(np.zeros(100, 100)) -ccddata = subtract_dark(ccddata, darkframe) +ccddata = ccdproc.subtract_dark(ccddata, masterdark) #correct for gain--once again gain should have a unit and even an error #associated with it. -ccddata = gain_correct(ccddata, gain=1.0) +ccddata = ccdproc.gain_correct(ccddata, gain=1.0) #Also the gain may be non-linear -ccddata = gain_correct(ccddata, gain=np.array([1.0, 0.5e-3])) +ccddata = ccdproc.gain_correct(ccddata, gain=np.array([1.0, 0.5e-3])) #although then this step should be apply before any other corrections #if it is non-linear, but that is more up to the person processing their #own data. @@ -149,7 +155,7 @@ #Not applicable for a single CCD situation #Error checks: the xtalkimage and image are the same shape xtalkimage = NDData.NDData(np.zeros(100, 100)) -ccddata = xtalk_correct(ccddata, xtalkimage, coef=1e-3) +ccddata = ccdproc.xtalk_correct(ccddata, xtalkimage, coef=1e-3) #flat field correction--this can either be a dome flat, sky flat, or an #illumination corrected image. This step should normalize by the value of the @@ -158,13 +164,14 @@ #Error checks: check for divive by zero #Features: If the flat is less than minvalue, minvalue is used flatimage = NDData.NDData(np.ones(100, 100)) -ccddata = flat_correct(ccddata, flatimage, minvalue=1) +ccddata = ccdproc.flat_correct(ccddata, flatimage, minvalue=1) #fringe correction or any correction that requires subtracting #off a potentially scaled image #Error checks: the flatimage and image are the same shape -fringemage = NDData.NDData(np.ones(100, 100)) -ccddata = fringe_correct(ccddata, fringeimage, scale=1, operation='multiple') +fringeimage = NDData.NDData(np.ones(100, 100)) +ccddata = ccdproc.fringe_correct(ccddata, fringeimage, scale=1, + operation='multiple') #cosmic ray cleaning step--this should have options for different #ways to do it with their associated steps. We also might want to @@ -172,12 +179,12 @@ #step should update the mask and flags. So the user could have options #to replace the cosmic rays, only flag the cosmic rays, or flag and #mask the cosmic rays, or all of the above. -ccddata = cosmicray_laplace(ccddata, method='laplace', **kwargs) -ccddata = cosmicray_median(ccddata, method='laplace', **kwargs) +ccddata = ccdproc.cosmicray_laplace(ccddata, method='laplace') +ccddata = ccdproc.cosmicray_median(ccddata, method='laplace') #Apply distortion corrections #Either update the WCS or transform the frame -ccddata = distortion_correct(ccddata, distortion) +ccddata = ccdproc.distortion_correct(ccddata, distortion) # ================ @@ -188,27 +195,28 @@ #select different functions to fit. #other options are reject parameters, number of iteractions #and/or convergernce limit -coef = iterfit(x, y, function='polynomial', order=3) +coef = ccdproc.iterfit(x, y, function='polynomial', order=3) #fit a 2-D function with iterative rejections and the ability to #select different functions to fit. #other options are reject parameters, number of iteractions #and/or convergernce limit -coef = iterfit(data, function='polynomial', order=3) +coef = ccdproc.iterfit(data, function='polynomial', order=3) #in addition to these operations, basic addition, subtraction # multiplication, and division should work for CCDDATA objects ccddata = ccddata + ccddata -ccddata = ccddata * 2 +ccddata2 = ccddata * 2 #combine a set of NDData objects -alldata = combine([ccddata, ccddata2], method='average', reject=None, **kwargs) +alldata = ccdproc.combine([ccddata, ccddata2], method='average', + reject=None) #re-sample the data to different binnings (either larger or smaller) -ccddata = rebin(ccddata, binning=(2, 2)) +ccddata = ccdproc.rebin(ccddata, binning=(2, 2)) #tranform the data--ie shift, rotate, etc #question--add convenience functions for image shifting and rotation? #should udpate WCS although that would actually be the prefered method -ccddata = transform(ccddata, transform, conserve_flux=True) +ccddata = ccdproc.transform(ccddata, transform, conserve_flux=True) From 62b672c617852858bee97c8805ee7c09da9b0e20 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Thu, 6 Feb 2014 12:23:19 -0600 Subject: [PATCH 05/24] Add API for a Keyword class used to access image metadata and a couple examples of its use --- ccdproc_api.py | 58 +++++++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 55 insertions(+), 3 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 5686b10..bfe7af4 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -87,12 +87,37 @@ ccddata.flags = np.zeros(110, 100) ccddata.wcs = None ccddata.meta = {} + +# making the metadata a Header instance gets us a case-insensitive dictionary... +assert isinstance(ccddata.meta, fits.Header) ccddata.units = u.adu # is this valid? #The ccddata class should have a functional form to create a CCDData #object directory from a fits file ccddata = ccdproc.CCDData.fromFITS('img.fits') +''' +Keyword is an object that represents a key, value pair for use in passing +data between functions in ``ccdproc``. The value is an astropy.units.Quantity, +with the unit specified explicitly when the Keyword instance is created. +The key is case-insensitive, and synonyms can be supplied that will be used +to look for the value in CCDData.meta. +''' +key = ccdproc.Keyword('exposure', unit=u.sec, synonyms=['exptime']) +header = fits.Header() +header['exposure'] = 15.0 +# value matched by keyword name exposure +value = key(header) +assert value == 15 * u.sec +del header['exposure'] +header['exptime'] = 15.0 +# value matched by synonym exptime +value = key(header) +assert value == 15 * u.sec +# inconsistent values in the header raise an error: +header['exposure'] = 25.0 +value = key(header) # raises ValueError + # Functional Requirements # ---------------------- # A number of these different fucntions are convenient functions that @@ -137,12 +162,39 @@ #correct for dark frames #Error checks: the masterbias and image are the same shape -masterdark = NDData.NDData(np.zeros(100, 100)) -ccddata = ccdproc.subtract_dark(ccddata, masterdark) +#Options: Exposure time of the data image and the master dark image can be +# specified as either an astropy.units.Quantity or as a ccdata.Keyword; +# in the second case the exposure time will be extracted from the +# metadata for each image. +masterdark = ccdproc.CCDData(np.zeros(100, 100)) +masterdark.meta['exptime'] = 30.0 +ccddata.meta['EXPOSURE'] = 15.0 + +exposure_time_key = ccdproc.Keyword('exposure', + unit=u.sec, + synonyms=['exptime']) + +# explicitly specify exposure times +ccddata = ccdproc.subtract_dark(ccddata, masterdark, + data_exposure=15 * u.sec, + dark_exposure=30 * u.sec, + scale=True + ) + +# get exposure times from metadata +ccddata = ccdproc.subtract_dark(ccddata, masterdark, + exposure_time=exposure_time_key, + scale=True) #correct for gain--once again gain should have a unit and even an error #associated with it. -ccddata = ccdproc.gain_correct(ccddata, gain=1.0) + +# gain can be specified as a Quantity... +ccddata = ccdproc.gain_correct(ccddata, gain=1.0 * u.ph / u.adu) +# ...or the gain can be specified as a ccdproc.Keyword: +gain_key = ccdproc.Keyword('gain', unit=u.ph / u.adu) +ccddata = ccdproc.gain_correct(ccddata, gain=gain_key) + #Also the gain may be non-linear ccddata = ccdproc.gain_correct(ccddata, gain=np.array([1.0, 0.5e-3])) #although then this step should be apply before any other corrections From 23a872104915445fafea195546cabe5a8dbd5329 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Tue, 11 Feb 2014 01:41:57 -0600 Subject: [PATCH 06/24] Add logging proposals --- ccdproc_api.py | 53 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 53 insertions(+) diff --git a/ccdproc_api.py b/ccdproc_api.py index bfe7af4..de5f8a1 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -239,6 +239,59 @@ ccddata = ccdproc.distortion_correct(ccddata, distortion) +# ======= +# Logging +# ======= + +# By logging we mean simply keeping track of what has been to each image in +# its as opposed to logging in the sense of the python logging module. Logging +# at that level is expected to be done by pipelines using the functions in +# ccdproc. + +# for the purposes of illustration this document describes how logging would +# be handled for subtract_bias; handling for other functions would be similar. + +# OPTION: One entry is added to the metadata for each processing step and the +# key added is the __name__ of the processing step. + +# Subtracting bias like this: + +ccddata = ccdproc.subtract_bias(ccddata, masterbias) + +# adds a keyword to the metadata: + +assert ccddata['subtract_bias'] # name is the __name__ of the + # processing step + +# this allows fairly easy checking of whether the processing step is being +# repeated. + +# OPTION: One entry is added to the metadata for each processing step and the +# key added is more human-friendly. + +# Subtracting bias like this: + +ccddata = ccdproc.subtract_bias(ccddata, masterbias) + +# adds a keyword to the metadata: + +assert ccddata.meta['bias_subtracted'] # name reads more naturally than previous + # option + +# OPTION: There is a single keyword called, e.g. calibration_status, whose +# value is a logical OR of flags defined for each potential operation. + +# in this option several package-level constants would be defined: +ccdproc.BIAS_DONE +ccdproc.DARK_DONE +# and so on. + +# Then logging has this effect: + +assert not (ccddata.meta['calibration_status'] & ccdproc.BIAS_DONE) +ccddata = ccdproc.subtract_bias(ccddata, masterbias) +assert (ccddata.meta['calibration_status'] & ccdproc.BIAS_DONE) + # ================ # Helper Functions # ================ From 7e77d23d3b4fb95b58dee377299d998f2dae136e Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 14 Feb 2014 21:24:30 -0600 Subject: [PATCH 07/24] Incorporate suggestions by @crawfordsm --- ccdproc_api.py | 26 ++++++++++++++------------ 1 file changed, 14 insertions(+), 12 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index de5f8a1..3b3bdb0 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -260,8 +260,8 @@ # adds a keyword to the metadata: -assert ccddata['subtract_bias'] # name is the __name__ of the - # processing step +assert ccddata.meta['subtract_bias'] # name is the __name__ of the + # processing step # this allows fairly easy checking of whether the processing step is being # repeated. @@ -278,19 +278,21 @@ assert ccddata.meta['bias_subtracted'] # name reads more naturally than previous # option -# OPTION: There is a single keyword called, e.g. calibration_status, whose -# value is a logical OR of flags defined for each potential operation. +# OPTION: Each of the processing steps allows the user to specify a keyword +# that is added to the metadata. An optional value can also be given. -# in this option several package-level constants would be defined: -ccdproc.BIAS_DONE -ccdproc.DARK_DONE -# and so on. +# SUB-OPTION: Allow the Keyword class to have string values? Then could use +# those as an argument for add_keywords -# Then logging has this effect: +# add keyword named SUBBIAS but don't set a value for it: -assert not (ccddata.meta['calibration_status'] & ccdproc.BIAS_DONE) -ccddata = ccdproc.subtract_bias(ccddata, masterbias) -assert (ccddata.meta['calibration_status'] & ccdproc.BIAS_DONE) +ccddata = ccdproc.subtract_bias(ccddata, masterbias, add_keyword='SUBBIAS') + +# add keyword CALSTAT with value 'B': + +ccddata = ccdproc.subtract_bias(ccddata, masterbias, + add_keyword='CALSTAT', + value='B') # ================ # Helper Functions From fe818ebf5706ded27067ebdc952c24c3e325facf Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Mon, 17 Feb 2014 20:40:57 -0600 Subject: [PATCH 08/24] Update user-specified keyword option to allow string or ccdproc.Keyword This also updates the API for ccdproc.Keyword to allow string values. --- ccdproc_api.py | 27 ++++++++++++++++----------- 1 file changed, 16 insertions(+), 11 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 3b3bdb0..95127ae 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -118,6 +118,14 @@ header['exposure'] = 25.0 value = key(header) # raises ValueError +# the value of a Keyword can also be set directly: +key.value = 20 * u.sec + +# String values are accommodated by setting the unit to the python type str: + +string_key = ccdproc.Keyword('filter', unit=str) +string_key.value = 'V' + # Functional Requirements # ---------------------- # A number of these different fucntions are convenient functions that @@ -279,20 +287,17 @@ # option # OPTION: Each of the processing steps allows the user to specify a keyword -# that is added to the metadata. An optional value can also be given. - -# SUB-OPTION: Allow the Keyword class to have string values? Then could use -# those as an argument for add_keywords - -# add keyword named SUBBIAS but don't set a value for it: +# that is added to the metadata. The keyword can either be a string or a +# ccdproc.Keyword instance +# add keyword as string: ccddata = ccdproc.subtract_bias(ccddata, masterbias, add_keyword='SUBBIAS') -# add keyword CALSTAT with value 'B': - -ccddata = ccdproc.subtract_bias(ccddata, masterbias, - add_keyword='CALSTAT', - value='B') +# add keyword/value using a ccdproc.Keyword object: +key = ccdproc.Keyword('calstat', unit=str) +key.value = 'B' +ccddata = ccdproc.subtract_bias(ccddata, masterbias, + add_keyword=key) # ================ # Helper Functions From ac603c5648d33d0572ef22bd13068ed7f815f408 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Mon, 17 Feb 2014 20:41:29 -0600 Subject: [PATCH 09/24] Fix assert statements for first two options and logging so they will work --- ccdproc_api.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 95127ae..9b6bfd8 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -268,8 +268,8 @@ # adds a keyword to the metadata: -assert ccddata.meta['subtract_bias'] # name is the __name__ of the - # processing step +assert 'subtract_bias' in ccddata.meta # name is the __name__ of the + # processing step # this allows fairly easy checking of whether the processing step is being # repeated. @@ -283,8 +283,8 @@ # adds a keyword to the metadata: -assert ccddata.meta['bias_subtracted'] # name reads more naturally than previous - # option +assert 'bias_subtracted' in ccddata.meta # name reads more naturally than + # previous option # OPTION: Each of the processing steps allows the user to specify a keyword # that is added to the metadata. The keyword can either be a string or a From 839022f714311c1c951466ce041a7677306835ef Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Tue, 11 Feb 2014 22:03:52 -0600 Subject: [PATCH 10/24] Initial attempt at image combination API --- ccdproc_api.py | 80 +++++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 79 insertions(+), 1 deletion(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 9b6bfd8..cba8bf4 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -8,6 +8,7 @@ from astropy.nddata import NDData from astropy.io import fits from astropy import units as u +from astropy.stats.funcs import sigma_clip import ccdproc ''' @@ -246,7 +247,6 @@ #Either update the WCS or transform the frame ccddata = ccdproc.distortion_correct(ccddata, distortion) - # ======= # Logging # ======= @@ -299,6 +299,84 @@ ccddata = ccdproc.subtract_bias(ccddata, masterbias, add_keyword=key) +# ================= +# Image combination +# ================= + +# The ``combine`` task from IRAF performs several functions: +# 1. Selection of images to be combined by image type with optional grouping +# into subsets. +# 2. Offsetting of images based on either user-specified shifts or on WCS +# information in the image metadata. +# 3. Rejection of pixels from inclusion in the combination based on masking, +# threshold rejection prior to any image scaling or zero offsets, and +# automatic rejection through a variety of algorithms (minmax, sigmaclip, +# ccdclip, etc) that allow for scaling, zero offset and in some cases +# weighting of the images being combined. +# 4. Scaling and/or zero offset of images before combining based on metadata +# (e.g. image exposure) or image statistics (e.g image median, mode or average +# determined by either an IRAF-selected subset of points or a region of the +# image supplied by the user). +# 5. Combination of remaining pixels by either median or average. +# 6. Logging of the choices made by IRAF in carrying out the operation (e.g. +# recording what zero offset was used for each image). + +# As much as is practical, the ccdproc API separates these functions, discussed +# in detail below. + +# 1. Image selection: this will not be provided by ccdproc (or at least not +# considered part of image combination). We assume that the user will have +# selected a set of images prior to beginning combination. + +# 2. Position offsets: offsets and other transforms are handled by +# ccdproc.transform, described below under "Helper Function" + +# 3. (One option: build up a list of rejected pixels) +# Masking: ccdpro.CCDData objects are already masked arrays, allowing +# automatic exclusion of masked pixels from all operations. +# +# Threshold rejection of all pixels with data value over 30000 or under -100: + +rejected_pixels = ccdproc.threshold_reject(ccddata1, ccddata2, ccddata3, + max=30000, min=-100) + +# automatic rejection by min/max, sigmaclip, ccdclip, etc. provided through +# one interface with separate helper functions + +# min/max +rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, + method=ccdproc.minmax) + +# sigmaclip (relies on astropy.stats.funcs) +rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, + method=sigma_clip, + sigma_high=3.0, sigma_low=2.0, + centerfunc=np.mean, + exclude_extrema=True) # are min/max pixels excluded? + +# ccdclip +rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, + method=ccdproc.ccdclip, + sigma_high=3.0, sigma_low=2.0, + gain=1.0, read_noise=5.0, + centerfunc=np.mean + exclude_extrema=True) + +# 4. Image scaling/zero offset with scaling determined by the mode of each +# image and the offset by the median + +scales, zero_offset = ccdproc.calc_weights(cddata1, ccddata2, ccddata3, + scale_by=np.mode, + offset_by=np.median) + +# 5. The actual combination + +combined_image = ccdproc.combine(cddata1, ccddata2, ccddata3, + reject=rejected_pixels, + scale=scales, + offset=zero_offset, + method=np.median) + # ================ # Helper Functions # ================ From 3e2b0dd4c5fa8f543be97d76d9f33cb2b5d8e52b Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 14 Feb 2014 21:43:45 -0600 Subject: [PATCH 11/24] Clarify handling of offsets/transforms --- ccdproc_api.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index cba8bf4..6c2b23e 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -328,8 +328,9 @@ # considered part of image combination). We assume that the user will have # selected a set of images prior to beginning combination. -# 2. Position offsets: offsets and other transforms are handled by -# ccdproc.transform, described below under "Helper Function" +# 2. Position offsets: In ccdproc this is not part of combine. Instead, +# offsets and other transforms are handled by ccdproc.transform, described +# below under "Helper Function" # 3. (One option: build up a list of rejected pixels) # Masking: ccdpro.CCDData objects are already masked arrays, allowing From f288e32a33e5468217ef2dfe2af2cc617029b1cc Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Mon, 17 Feb 2014 21:19:58 -0600 Subject: [PATCH 12/24] Perform image combination with a Combiner class --- ccdproc_api.py | 70 ++++++++++++++++++++++++++++---------------------- 1 file changed, 39 insertions(+), 31 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 6c2b23e..6b7a92d 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -127,6 +127,12 @@ string_key = ccdproc.Keyword('filter', unit=str) string_key.value = 'V' +''' +Combiner is an object to facilitate image combination. Its methods perform +the steps that are bundled into the IRAF task combine + +It is discussed in detail below, in the section "Image combination" +''' # Functional Requirements # ---------------------- # A number of these different fucntions are convenient functions that @@ -306,7 +312,7 @@ # The ``combine`` task from IRAF performs several functions: # 1. Selection of images to be combined by image type with optional grouping # into subsets. -# 2. Offsetting of images based on either user-specified shifts or on WCS +# 2. Offsetting of images based on either user-specified shifts or on WCS # information in the image metadata. # 3. Rejection of pixels from inclusion in the combination based on masking, # threshold rejection prior to any image scaling or zero offsets, and @@ -314,9 +320,9 @@ # ccdclip, etc) that allow for scaling, zero offset and in some cases # weighting of the images being combined. # 4. Scaling and/or zero offset of images before combining based on metadata -# (e.g. image exposure) or image statistics (e.g image median, mode or average -# determined by either an IRAF-selected subset of points or a region of the -# image supplied by the user). +# (e.g. image exposure) or image statistics (e.g image median, mode or +# average determined by either an IRAF-selected subset of points or a +# region of the image supplied by the user). # 5. Combination of remaining pixels by either median or average. # 6. Logging of the choices made by IRAF in carrying out the operation (e.g. # recording what zero offset was used for each image). @@ -332,51 +338,53 @@ # offsets and other transforms are handled by ccdproc.transform, described # below under "Helper Function" -# 3. (One option: build up a list of rejected pixels) +# The combination process begins with the creation of a Combiner object, +# initialized with the images to be combined + +combine = ccdproc.Combiner(ccddata1, ccddata2, ccddata3) + +# 3. # Masking: ccdpro.CCDData objects are already masked arrays, allowing # automatic exclusion of masked pixels from all operations. # # Threshold rejection of all pixels with data value over 30000 or under -100: -rejected_pixels = ccdproc.threshold_reject(ccddata1, ccddata2, ccddata3, - max=30000, min=-100) +combine.threshold_reject(max=30000, min=-100) # automatic rejection by min/max, sigmaclip, ccdclip, etc. provided through -# one interface with separate helper functions +# one method, with different helper functions # min/max -rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, - method=ccdproc.minmax) +combine.clip(method=ccdproc.minmax) # sigmaclip (relies on astropy.stats.funcs) -rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, - method=sigma_clip, - sigma_high=3.0, sigma_low=2.0, - centerfunc=np.mean, - exclude_extrema=True) # are min/max pixels excluded? +combine.clip(method=sigma_clip, + sigma_high=3.0, sigma_low=2.0, + centerfunc=np.mean, + exclude_extrema=True) # min/max pixels can be excluded in the + # sigma_clip. IRAF's clip excludes them # ccdclip -rejected_pixels = ccdproc.clip(ccddata1, ccddata2, ccddata3, - method=ccdproc.ccdclip, - sigma_high=3.0, sigma_low=2.0, - gain=1.0, read_noise=5.0, - centerfunc=np.mean - exclude_extrema=True) +combine.clip(method=ccdproc.ccdclip, + sigma_high=3.0, sigma_low=2.0, + gain=1.0, read_noise=5.0, + centerfunc=np.mean, + exclude_extrema=True) # 4. Image scaling/zero offset with scaling determined by the mode of each -# image and the offset by the median +# image and the offset by the median. This method calculates what are +# effectively weights for each image + +combine.calc_weights(scale_by=np.mode, offset_by=np.median) -scales, zero_offset = ccdproc.calc_weights(cddata1, ccddata2, ccddata3, - scale_by=np.mode, - offset_by=np.median) +# 5. The actual combination -- a couple of ways images can be combined -# 5. The actual combination +# median; the excluded pixels based on the individual image masks, threshold +# rejection, clipping, etc, are wrapped into a single mask for each image +combined_image = combine(method=np.ma.median) -combined_image = ccdproc.combine(cddata1, ccddata2, ccddata3, - reject=rejected_pixels, - scale=scales, - offset=zero_offset, - method=np.median) +# average; in this case image weights can also be specified +combined_image = combine(method=np.ma.mean, weights=[0.5, 1, 2]) # ================ # Helper Functions From d6fc09bf784ba5f8ca740ed27e31585e5805152a Mon Sep 17 00:00:00 2001 From: Steven Crawford Date: Tue, 11 Mar 2014 12:23:37 +0200 Subject: [PATCH 13/24] Update the reader information to use the astropy.io.registry --- ccdproc_api.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 6b7a92d..483cd26 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -95,7 +95,13 @@ #The ccddata class should have a functional form to create a CCDData #object directory from a fits file -ccddata = ccdproc.CCDData.fromFITS('img.fits') +ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits') + +#This function should then be registered with astropy.io.registry so +# the standard way for reading in a fits image will be +# the following: +ccddata = ccdproc.CCDData.read('img.fits', format='fits') + ''' Keyword is an object that represents a key, value pair for use in passing From 16b18789e46071c28b378ca9c77398923f3d77f3 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 4 Apr 2014 10:09:55 -0500 Subject: [PATCH 14/24] Simplify behavior of Keyword class In response to comments by @embray --- ccdproc_api.py | 15 +++------------ 1 file changed, 3 insertions(+), 12 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 483cd26..dbd377f 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -107,23 +107,14 @@ Keyword is an object that represents a key, value pair for use in passing data between functions in ``ccdproc``. The value is an astropy.units.Quantity, with the unit specified explicitly when the Keyword instance is created. -The key is case-insensitive, and synonyms can be supplied that will be used -to look for the value in CCDData.meta. +The key is case-insensitive. ''' -key = ccdproc.Keyword('exposure', unit=u.sec, synonyms=['exptime']) +key = ccdproc.Keyword('exposure', unit=u.sec) header = fits.Header() header['exposure'] = 15.0 # value matched by keyword name exposure -value = key(header) +value = key.value_from(header) assert value == 15 * u.sec -del header['exposure'] -header['exptime'] = 15.0 -# value matched by synonym exptime -value = key(header) -assert value == 15 * u.sec -# inconsistent values in the header raise an error: -header['exposure'] = 25.0 -value = key(header) # raises ValueError # the value of a Keyword can also be set directly: key.value = 20 * u.sec From 3b5b776bf628cfc9cee26e18cba5cd2bd6899371 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Fri, 4 Apr 2014 10:10:12 -0500 Subject: [PATCH 15/24] Clarify handling of units --- ccdproc_api.py | 22 +++++++++++++++++++++- 1 file changed, 21 insertions(+), 1 deletion(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index dbd377f..335e911 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -102,6 +102,8 @@ # the following: ccddata = ccdproc.CCDData.read('img.fits', format='fits') +# CCDData should raise an error if no unit is provided and a unit cannot be +# extracted from the FITS header. ''' Keyword is an object that represents a key, value pair for use in passing @@ -149,7 +151,25 @@ #and a unit? Or a scalar, unit and an error? ie, This #could actually be handled by the gain and readnoise being #specified as an NDData object -ccddata = ccdproc.createvariance(ccddata, gain=1.0, readnoise=5.0) + +ccddata.unit = u.adu +# The call below should raise an error because gain and readnoise are provided +# without units. +ccddata = ccdproc.create_variance(ccddata, gain=1.0, readnoise=5.0) + +# The electron unit is provided by ccddata +gain = 1.0 ccddata.electron / u.adu +readnoise = 5.0 * u.electron +# This succeeds because the units are consistent +ccddata = ccdproc.create_variance(ccddata, gain-gain, readnoise=readnoise) + +# with the gain units below there is a mismatch between the units of the +# image after scaling by the gain and the units of the readnoise... +gain = 1.0 u.photon / u.adu + +# ...so this should fail with an error. +ccddata = ccdproc.create_variance(ccddata, gain-gain, readnoise=readnoise) + #Overscan subtract the data #Should be able to provide the meta data for From 2a5fb642d00a63747ef3973e345d61aebf2e9d22 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Wed, 9 Apr 2014 10:41:54 -0500 Subject: [PATCH 16/24] Explicitly state that units are required for initialization --- ccdproc_api.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 335e911..0299802 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -77,9 +77,11 @@ ''' data = 100 + 10 * np.random.random((110, 100)) -ccddata = ccdproc.CCDData(data=data) -ccddata = ccdproc.CCDData(NDData.NDData(data)) -ccddata = ccdproc.CCDData(fits.ImageHDU(data)) +# initializing without a unit raises an error +ccddata = ccdproc.CCDData(data=data) # ValueError +ccddata = ccdproc.CCDData(data=data, unit=u.adu) +ccddata = ccdproc.CCDData(NDData.NDData(data), unit=u.photon) +ccddata = ccdproc.CCDData(fits.ImageHDU(data), unit=u.adu) #Setting basic properties of the object # ---------------------- From d2292c3cd3396b0be98189a62ae1a1e37380575e Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Wed, 9 Apr 2014 10:42:20 -0500 Subject: [PATCH 17/24] Add detail to i/o portion of api --- ccdproc_api.py | 44 ++++++++++++++++++++++++++++++++++++++------ 1 file changed, 38 insertions(+), 6 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 0299802..d5289b4 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -97,15 +97,47 @@ #The ccddata class should have a functional form to create a CCDData #object directory from a fits file -ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits') +ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits', image_unit=u.adu) -#This function should then be registered with astropy.io.registry so -# the standard way for reading in a fits image will be +# omitting a unit causes an error for now -- in the future an attempt should +# be made to extract the unit from the FITS header. +ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits') # raises ValueError + +# If a file has multiple extensions the desired hdu should be specified. It +# defaults to zero. +ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', + image_unit=u.adu, + hdu=2) + + +# This function should then be registered with astropy.io.registry, and the +# FITS format auto-identified using the fits.connect.is_fits so +# the standard way for reading in a fits image will be # the following: -ccddata = ccdproc.CCDData.read('img.fits', format='fits') +ccddata = ccdproc.CCDData.read('img.fits', image_unit=u.adu) + + +# CCDData raises an error if no unit is provided; eventually an attempt to +# extract the unit from the FITS header should be made. + +# Writing is handled in a similar fashion; the image ccddata is written to +# the file img2.fits with: +ccdproc.CCDData.fits_ccddata_writer(ccddata, 'img2.fits') + +# all additional keywords are passed on to the underlying FITS writer, e.g. +ccdproc.CCDData.fits_ccddata_writer(ccddata, 'img2.fits', clobber=True) + +# The writer is registered with unified io so that in practice a user will do +ccddata.write('img2.fits') + +# NOTE: for now any flag, mask and/or unit for ccddata is discarded when +# writing. If you want all or some of that information preserved you must +# create the FITS files manually. -# CCDData should raise an error if no unit is provided and a unit cannot be -# extracted from the FITS header. +# To be completely explicit about this: +ccddata2 = ccdproc.CCDData.read('img2.fits', image_unit=u.adu) +assert ccddata2.mask is None # even though we set ccddata.mask before saving +assert ccddata2.flag is None # even though we set ccddata.flag before saving ''' Keyword is an object that represents a key, value pair for use in passing From 0bdcb50aee83c93c489d9b6bade37d20fc260714 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Wed, 9 Apr 2014 11:03:24 -0500 Subject: [PATCH 18/24] Pass read keywords to FITS reader --- ccdproc_api.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index d5289b4..98d6e06 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -109,7 +109,10 @@ image_unit=u.adu, hdu=2) - +# any additional keywords are passed through to fits.open, e.g. +ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', + image_unit=u.adu, + do_not_scale_image_data=True) # This function should then be registered with astropy.io.registry, and the # FITS format auto-identified using the fits.connect.is_fits so # the standard way for reading in a fits image will be @@ -135,6 +138,10 @@ # create the FITS files manually. # To be completely explicit about this: +ccddata.mask = np.ones(110, 100) +ccddata.flags = np.zeros(110, 100) +ccddata.write('img2.fits') + ccddata2 = ccdproc.CCDData.read('img2.fits', image_unit=u.adu) assert ccddata2.mask is None # even though we set ccddata.mask before saving assert ccddata2.flag is None # even though we set ccddata.flag before saving From 8ff81757b2c73c4e4705f2ce824abb98d200ab2e Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Wed, 9 Apr 2014 14:23:08 -0500 Subject: [PATCH 19/24] User may not disable scaling when reading data from a FITS file. --- ccdproc_api.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 98d6e06..ae31b7d 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -109,10 +109,17 @@ image_unit=u.adu, hdu=2) -# any additional keywords are passed through to fits.open, e.g. +# any additional keywords are passed through to fits.open, with the exception +# of those related scaling (do_not_scale_image_data and scale_back) +ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', + image_unit=u.adu, + memmap=False) +# The call below raises a TypeErrror ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', image_unit=u.adu, do_not_scale_image_data=True) + + # This function should then be registered with astropy.io.registry, and the # FITS format auto-identified using the fits.connect.is_fits so # the standard way for reading in a fits image will be From f31e5dbe70ae976cd0b841f5e077377196caa4aa Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Wed, 9 Apr 2014 21:14:31 -0500 Subject: [PATCH 20/24] Describe method for converting CCDData object to an hdu. This was called toFITS in the initial draft of the CCDData code though it was not mentioned in earlier drafts of the API --- ccdproc_api.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/ccdproc_api.py b/ccdproc_api.py index ae31b7d..5a698c1 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -153,6 +153,10 @@ assert ccddata2.mask is None # even though we set ccddata.mask before saving assert ccddata2.flag is None # even though we set ccddata.flag before saving +# CCDData provides a convenience method to construct a FITS HDU from the data +# and metadata +hdu = ccddata.to_hdu() + ''' Keyword is an object that represents a key, value pair for use in passing data between functions in ``ccdproc``. The value is an astropy.units.Quantity, From 2ec4280455ac33bf7b6a96c7783b4d5eb2e2dacd Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Thu, 10 Apr 2014 15:41:16 -0500 Subject: [PATCH 21/24] Simplify handling of keywords that are strings --- ccdproc_api.py | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 5a698c1..27f9d1a 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -173,11 +173,17 @@ # the value of a Keyword can also be set directly: key.value = 20 * u.sec -# String values are accommodated by setting the unit to the python type str: +# String values are accommodated by not setting the unit and setting the value +# to a string -string_key = ccdproc.Keyword('filter', unit=str) +string_key = ccdproc.Keyword('filter') string_key.value = 'V' +# Setting a string value when a unit has been specified is an error + +bad_key = ccdproc.Keyword('exposure', unit=u.sec) +bad_key.value = '30' # raise a ValueError + ''' Combiner is an object to facilitate image combination. Its methods perform the steps that are bundled into the IRAF task combine From 5171b3a825f3e80b3b8d8da4b9bdb185da38d199 Mon Sep 17 00:00:00 2001 From: Matthew Craig Date: Thu, 10 Apr 2014 20:45:17 -0500 Subject: [PATCH 22/24] Fix typo --- ccdproc_api.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 27f9d1a..8d4361a 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -219,14 +219,14 @@ gain = 1.0 ccddata.electron / u.adu readnoise = 5.0 * u.electron # This succeeds because the units are consistent -ccddata = ccdproc.create_variance(ccddata, gain-gain, readnoise=readnoise) +ccddata = ccdproc.create_variance(ccddata, gain=gain, readnoise=readnoise) # with the gain units below there is a mismatch between the units of the # image after scaling by the gain and the units of the readnoise... gain = 1.0 u.photon / u.adu # ...so this should fail with an error. -ccddata = ccdproc.create_variance(ccddata, gain-gain, readnoise=readnoise) +ccddata = ccdproc.create_variance(ccddata, gain=gain, readnoise=readnoise) #Overscan subtract the data From 26e5493b15041b1f0ddde87ad5adeb160f05deb2 Mon Sep 17 00:00:00 2001 From: Steven Crawford Date: Tue, 22 Apr 2014 20:44:54 +0200 Subject: [PATCH 23/24] simplified and clarified combiner --- ccdproc_api.py | 20 ++++++-------------- 1 file changed, 6 insertions(+), 14 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 8d4361a..7d09ab9 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -414,23 +414,15 @@ # below under "Helper Function" # The combination process begins with the creation of a Combiner object, -# initialized with the images to be combined +# initialized with a list of the images to be combined -combine = ccdproc.Combiner(ccddata1, ccddata2, ccddata3) - -# 3. -# Masking: ccdpro.CCDData objects are already masked arrays, allowing -# automatic exclusion of masked pixels from all operations. -# -# Threshold rejection of all pixels with data value over 30000 or under -100: - -combine.threshold_reject(max=30000, min=-100) +combine = ccdproc.Combiner([ccddata1, ccddata2, ccddata3]) # automatic rejection by min/max, sigmaclip, ccdclip, etc. provided through # one method, with different helper functions # min/max -combine.clip(method=ccdproc.minmax) +combine.clip(method=ccdproc.minmax, max_data=30000, data_min=-100) # sigmaclip (relies on astropy.stats.funcs) combine.clip(method=sigma_clip, @@ -449,7 +441,6 @@ # 4. Image scaling/zero offset with scaling determined by the mode of each # image and the offset by the median. This method calculates what are # effectively weights for each image - combine.calc_weights(scale_by=np.mode, offset_by=np.median) # 5. The actual combination -- a couple of ways images can be combined @@ -458,8 +449,9 @@ # rejection, clipping, etc, are wrapped into a single mask for each image combined_image = combine(method=np.ma.median) -# average; in this case image weights can also be specified -combined_image = combine(method=np.ma.mean, weights=[0.5, 1, 2]) +# average; in this case image weights can also be specified if they +# have been +combined_image = combine(method=np.ma.mean) # ================ # Helper Functions From 9cf829c7c83731f19c564a2d4b98315d144904cc Mon Sep 17 00:00:00 2001 From: Steven Crawford Date: Tue, 6 May 2014 22:06:08 +0200 Subject: [PATCH 24/24] updates for current status of ccdproc --- ccdproc_api.py | 245 ++++++++++++++++++++++++++----------------------- 1 file changed, 130 insertions(+), 115 deletions(-) diff --git a/ccdproc_api.py b/ccdproc_api.py index 7d09ab9..8618220 100644 --- a/ccdproc_api.py +++ b/ccdproc_api.py @@ -10,7 +10,9 @@ from astropy import units as u from astropy.stats.funcs import sigma_clip -import ccdproc +from ccdproc import ccddata +from ccdproc import ccdproc + ''' The ccdproc package provides tools for the reduction and analysis of optical data captured with a CCD. The package @@ -78,44 +80,43 @@ ''' data = 100 + 10 * np.random.random((110, 100)) # initializing without a unit raises an error -ccddata = ccdproc.CCDData(data=data) # ValueError -ccddata = ccdproc.CCDData(data=data, unit=u.adu) -ccddata = ccdproc.CCDData(NDData.NDData(data), unit=u.photon) -ccddata = ccdproc.CCDData(fits.ImageHDU(data), unit=u.adu) +ccd = ccddata.CCDData(data=data) # ValueError + +ccd = ccddata.CCDData(data=data, unit=u.adu) +ccd = ccddata.CCDData(NDData(data), unit=u.photon) #Setting basic properties of the object # ---------------------- -ccddata.variance = data**0.5 -ccddata.mask = np.ones(110, 100) -ccddata.flags = np.zeros(110, 100) +ccddata.uncertainty = data**0.5 +ccddata.mask = np.ones((110, 100), dtype=bool) +ccddata.flags = np.zeros((110, 100)) ccddata.wcs = None ccddata.meta = {} - -# making the metadata a Header instance gets us a case-insensitive dictionary... -assert isinstance(ccddata.meta, fits.Header) +ccddata.header = {} #header and meta are interchangable ccddata.units = u.adu # is this valid? + #The ccddata class should have a functional form to create a CCDData #object directory from a fits file -ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits', image_unit=u.adu) +ccd = ccddata.fits_ccddata_reader('img.fits', unit=u.adu) # omitting a unit causes an error for now -- in the future an attempt should # be made to extract the unit from the FITS header. -ccddata = ccdproc.CCDData.fits_ccddata_reader('img.fits') # raises ValueError +ccd = ccddata.fits_ccddata_reader('img.fits') # raises ValueError # If a file has multiple extensions the desired hdu should be specified. It # defaults to zero. -ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', +ccd = ccddata.fits_ccddata_reader('multi_extension.fits', image_unit=u.adu, hdu=2) # any additional keywords are passed through to fits.open, with the exception # of those related scaling (do_not_scale_image_data and scale_back) -ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', +ccd = ccddata.fits_ccddata_reader('multi_extension.fits', image_unit=u.adu, memmap=False) # The call below raises a TypeErrror -ccddata = ccdproc.CCDData.fits_ccddata_reader('multi_extension.fits', +ccd = ccddata.fits_ccddata_reader('multi_extension.fits', image_unit=u.adu, do_not_scale_image_data=True) @@ -124,38 +125,38 @@ # FITS format auto-identified using the fits.connect.is_fits so # the standard way for reading in a fits image will be # the following: -ccddata = ccdproc.CCDData.read('img.fits', image_unit=u.adu) +ccd = ccdproc.CCDData.read('img.fits', image_unit=u.adu) # CCDData raises an error if no unit is provided; eventually an attempt to -# extract the unit from the FITS header should be made. +# extract the unit from the FITS header should be made. # Writing is handled in a similar fashion; the image ccddata is written to # the file img2.fits with: -ccdproc.CCDData.fits_ccddata_writer(ccddata, 'img2.fits') +ccdata.fits_ccddata_writer(ccd, 'img2.fits') # all additional keywords are passed on to the underlying FITS writer, e.g. -ccdproc.CCDData.fits_ccddata_writer(ccddata, 'img2.fits', clobber=True) +ccddata.fits_ccddata_writer(ccd, 'img2.fits', clobber=True) # The writer is registered with unified io so that in practice a user will do -ccddata.write('img2.fits') +ccd.write('img2.fits') # NOTE: for now any flag, mask and/or unit for ccddata is discarded when # writing. If you want all or some of that information preserved you must # create the FITS files manually. -# To be completely explicit about this: -ccddata.mask = np.ones(110, 100) -ccddata.flags = np.zeros(110, 100) -ccddata.write('img2.fits') +# To be completely explicit about not writing out addition information: +ccd.mask = np.ones(110, 100) +ccd.flags = np.zeros(110, 100) +ccd.write('img2.fits') -ccddata2 = ccdproc.CCDData.read('img2.fits', image_unit=u.adu) -assert ccddata2.mask is None # even though we set ccddata.mask before saving -assert ccddata2.flag is None # even though we set ccddata.flag before saving +ccd2 = ccddata.CCDData.read('img2.fits', image_unit=u.adu) +assert ccd2.mask is None # even though we set ccddata.mask before saving +assert ccd2.flag is None # even though we set ccddata.flag before saving # CCDData provides a convenience method to construct a FITS HDU from the data # and metadata -hdu = ccddata.to_hdu() +hdu = ccd.to_hdu() ''' Keyword is an object that represents a key, value pair for use in passing @@ -163,15 +164,15 @@ with the unit specified explicitly when the Keyword instance is created. The key is case-insensitive. ''' -key = ccdproc.Keyword('exposure', unit=u.sec) +key = ccdproc.Keyword('exposure', unit=u.s) header = fits.Header() header['exposure'] = 15.0 # value matched by keyword name exposure value = key.value_from(header) -assert value == 15 * u.sec +assert value == 15 * u.s # the value of a Keyword can also be set directly: -key.value = 20 * u.sec +key.value = 20 * u.s # String values are accommodated by not setting the unit and setting the value # to a string @@ -181,98 +182,92 @@ # Setting a string value when a unit has been specified is an error -bad_key = ccdproc.Keyword('exposure', unit=u.sec) +bad_key = ccdproc.Keyword('exposure', unit=u.s) bad_key.value = '30' # raise a ValueError +''' Functional Requirements + ---------------------- + A number of these different fucntions are convenient functions that + just outline the process that is needed. The important thing is that + the process is being logged and that a clear process is being handled + by each step to make building a pipeline easy. Then again, it might + not be easy to handle all possible steps which are needed, and the more + important steps will be things which aren't already handled by NDData. + + All functions should propogate throught to the variance frame and + bad pixel mask ''' -Combiner is an object to facilitate image combination. Its methods perform -the steps that are bundled into the IRAF task combine -It is discussed in detail below, in the section "Image combination" -''' -# Functional Requirements -# ---------------------- -# A number of these different fucntions are convenient functions that -# just outline the process that is needed. The important thing is that -# the process is being logged and that a clear process is being handled -# by each step to make building a pipeline easy. Then again, it might -# not be easy to handle all possible steps which are needed, and the more -# important steps will be things which aren't already handled by NDData. - -#All functions should propogate throught to the variance frame and -#bad pixel mask - -#convenience function based on a given value for -#the readnoise and gain. Units should be specified -#for these values though. -#Question: Do we have an object that is just a scalar -#and a unit? Or a scalar, unit and an error? ie, This -#could actually be handled by the gain and readnoise being -#specified as an NDData object - -ccddata.unit = u.adu + +ccd.unit = u.adu # The call below should raise an error because gain and readnoise are provided -# without units. +# without units. Gain and rdnoise should be Quantities ccddata = ccdproc.create_variance(ccddata, gain=1.0, readnoise=5.0) # The electron unit is provided by ccddata -gain = 1.0 ccddata.electron / u.adu -readnoise = 5.0 * u.electron +gain = 1.0 * ccddata.electron / u.adu +readnoise = 5.0 * ccddata.electron # This succeeds because the units are consistent -ccddata = ccdproc.create_variance(ccddata, gain=gain, readnoise=readnoise) +ccd = ccdproc.create_variance(ccd, gain=gain, readnoise=readnoise) # with the gain units below there is a mismatch between the units of the # image after scaling by the gain and the units of the readnoise... -gain = 1.0 u.photon / u.adu +gain = 1.0 * u.photon / u.adu # ...so this should fail with an error. -ccddata = ccdproc.create_variance(ccddata, gain=gain, readnoise=readnoise) +ccd = ccdproc.create_variance(ccd, gain=gain, readnoise=readnoise) + + +#Overscan subtract the data by providing the slice of ccd with +#the overscan section. This will subtract off the median +#of each row +ccd = ccdproc.subtract_overscan(ccd, overscan=ccd[:,100:110]) + +#Overscan subtract the data by providing the slice of ccd with +#the overscan section. +ccd = ccdproc.subtract_overscan(ccd, fits_section='[101:110,:]') + +#For the overscan region the astropy.model to fit to the data can +#be specified -#Overscan subtract the data -#Should be able to provide the meta data for -#the keyworkd or provide a section to define and -#possible an axis to specify the oritation of the -#Question: Best way to specify the section? Should it be given -#Error Checks: That the section is within the image -ccddata = ccdproc.subtract_overscan(ccddata, section='[:,100:110]', - function='polynomial', order=3) #trim the images--the section gives the part of the image to keep -#That the trim section is within the image -ccddata = ccdproc.trim_image(ccddata, section='[0:100,0:100]') +#That the trim section is within the image. +ccd = ccdproc.trim_image(ccd[0:100,0:100]) + +#Using a section defined by fits convention +ccd = ccdproc.trim_image(ccd, fits_section='[1:100,1:100]') #subtract the master bias. Although this is a convenience function as #subtracting the two arrays will do the same thing. This should be able -#to handle logging of subtracting it off (logging should be added to NDData -#and then this is really just a convenience function -#Error checks: the masterbias and image are the same shape -masterbias = NDData.NDData(np.zeros(100, 100)) -ccddata = ccdproc.subtract_bias(ccddata, masterbias) +#to handle logging of subtracting it off +#Error checks: the masterbias and image are the same shape and units +masterbias = ccddata.CCDData(np.zeros((100, 100)), unit=u.adu) +ccd = ccdproc.subtract_bias(ccd, masterbias) #correct for dark frames -#Error checks: the masterbias and image are the same shape #Options: Exposure time of the data image and the master dark image can be # specified as either an astropy.units.Quantity or as a ccdata.Keyword; # in the second case the exposure time will be extracted from the # metadata for each image. -masterdark = ccdproc.CCDData(np.zeros(100, 100)) +masterdark = ccddata.CCDData(np.zeros((100, 100))+10, unit=u.adu) masterdark.meta['exptime'] = 30.0 ccddata.meta['EXPOSURE'] = 15.0 exposure_time_key = ccdproc.Keyword('exposure', - unit=u.sec, + unit=u.s, synonyms=['exptime']) # explicitly specify exposure times -ccddata = ccdproc.subtract_dark(ccddata, masterdark, - data_exposure=15 * u.sec, - dark_exposure=30 * u.sec, +ccd = ccdproc.subtract_dark(ccd, masterdark, + data_exposure=15 * u.s, + dark_exposure=30 * u.s, scale=True ) # get exposure times from metadata -ccddata = ccdproc.subtract_dark(ccddata, masterdark, +ccd = ccdproc.subtract_dark(ccd, masterdark, exposure_time=exposure_time_key, scale=True) @@ -280,13 +275,15 @@ #associated with it. # gain can be specified as a Quantity... -ccddata = ccdproc.gain_correct(ccddata, gain=1.0 * u.ph / u.adu) +ccd = ccdproc.gain_correct(ccd, gain=1.0 * u.ph / u.adu) # ...or the gain can be specified as a ccdproc.Keyword: gain_key = ccdproc.Keyword('gain', unit=u.ph / u.adu) -ccddata = ccdproc.gain_correct(ccddata, gain=gain_key) +ccd = ccdproc.gain_correct(ccd, gain=gain_key) #Also the gain may be non-linear -ccddata = ccdproc.gain_correct(ccddata, gain=np.array([1.0, 0.5e-3])) +#TODO: Not impliement in v0.1 +ccd = ccdproc.gain_correct(ccd, gain=np.array([1.0, 0.5e-3])) + #although then this step should be apply before any other corrections #if it is non-linear, but that is more up to the person processing their #own data. @@ -296,8 +293,9 @@ #general because this will be dependent on the CCD and the set up of the CCD. #Not applicable for a single CCD situation #Error checks: the xtalkimage and image are the same shape -xtalkimage = NDData.NDData(np.zeros(100, 100)) -ccddata = ccdproc.xtalk_correct(ccddata, xtalkimage, coef=1e-3) +#TODO: Not impliement in v0.1 +xtalkimage = ccddata.CCDData(np.zeros((100, 100))+10, unit=u.adu) +ccd = ccdproc.xtalk_correct(ccddata, xtalkimage, coef=1e-3) #flat field correction--this can either be a dome flat, sky flat, or an #illumination corrected image. This step should normalize by the value of the @@ -305,14 +303,15 @@ #Error checks: the flatimage and image are the same shape #Error checks: check for divive by zero #Features: If the flat is less than minvalue, minvalue is used -flatimage = NDData.NDData(np.ones(100, 100)) -ccddata = ccdproc.flat_correct(ccddata, flatimage, minvalue=1) +flatimage = ccddata.CCDData(np.zeros((100, 100))+10, unit=u.adu) +ccd = ccdproc.flat_correct(ccd, flatimage, minvalue=1) #fringe correction or any correction that requires subtracting #off a potentially scaled image #Error checks: the flatimage and image are the same shape -fringeimage = NDData.NDData(np.ones(100, 100)) -ccddata = ccdproc.fringe_correct(ccddata, fringeimage, scale=1, +#TODO: Not impliement in v0.1 +fringeimage = ccddata.CCDData(np.zeros((100, 100))+10, unit=u.adu) +ccddata = ccdproc.fringe_correct(ccd, fringeimage, scale=1, operation='multiple') #cosmic ray cleaning step--this should have options for different @@ -321,11 +320,15 @@ #step should update the mask and flags. So the user could have options #to replace the cosmic rays, only flag the cosmic rays, or flag and #mask the cosmic rays, or all of the above. -ccddata = ccdproc.cosmicray_laplace(ccddata, method='laplace') -ccddata = ccdproc.cosmicray_median(ccddata, method='laplace') +#median correct the cosmic rays +ccd = ccdproc.cosmicray_clean(ccd, thresh = 5, cr_func=ccdproc.cosmicray_median) + +#remove cosmic rays following van dokkum method in LA COSMIC +ccddata = ccdproc.cosmicray_laplace(ccddata, thresh = 5, cr_func=ccdproc.cosmicray_lapace) #Apply distortion corrections #Either update the WCS or transform the frame +#TODO: Add after v0.1 when WCS is working ccddata = ccdproc.distortion_correct(ccddata, distortion) # ======= @@ -416,42 +419,50 @@ # The combination process begins with the creation of a Combiner object, # initialized with a list of the images to be combined -combine = ccdproc.Combiner([ccddata1, ccddata2, ccddata3]) +from ccdproc import combiner +combine = combiner.Combiner([ccddata1, ccddata2, ccddata3]) # automatic rejection by min/max, sigmaclip, ccdclip, etc. provided through # one method, with different helper functions # min/max -combine.clip(method=ccdproc.minmax, max_data=30000, data_min=-100) +combine.minmax_clipping(method=ccdproc.minmax, max_data=30000, data_min=-100) # sigmaclip (relies on astropy.stats.funcs) -combine.clip(method=sigma_clip, - sigma_high=3.0, sigma_low=2.0, - centerfunc=np.mean, - exclude_extrema=True) # min/max pixels can be excluded in the - # sigma_clip. IRAF's clip excludes them +combine.sigma_clipping(low_thresh = 3.0, + high_thresh=3.0, + func=np.mean, + dev_func=ma.std) + +# TODO: min/max pixels can be excluded in the sigma_clip. IRAF's clip excludes them +# Add exclude_extrema to min/max function # ccdclip -combine.clip(method=ccdproc.ccdclip, - sigma_high=3.0, sigma_low=2.0, +# TODO: Not implemented in v0.1. +#Can sigma_clipping be updated to use it? +combine.ccdclip_clipping( sigma_high=3.0, sigma_low=2.0, gain=1.0, read_noise=5.0, centerfunc=np.mean, exclude_extrema=True) -# 4. Image scaling/zero offset with scaling determined by the mode of each -# image and the offset by the median. This method calculates what are -# effectively weights for each image -combine.calc_weights(scale_by=np.mode, offset_by=np.median) +# 4. Image scaling/zero offset with scaling are set by the user prior +# to passing the arrays to the task +# TODO: Implement some scaling + +# Image weights for use in the average combination can also be +# set by the user +# TODO: Implement some method for calculating the weighting +combine.weights = weights # 5. The actual combination -- a couple of ways images can be combined # median; the excluded pixels based on the individual image masks, threshold # rejection, clipping, etc, are wrapped into a single mask for each image -combined_image = combine(method=np.ma.median) +combined_image = combine.median_combine() # average; in this case image weights can also be specified if they # have been -combined_image = combine(method=np.ma.mean) +combined_image = combine.average_combine() # ================ # Helper Functions @@ -461,19 +472,23 @@ #select different functions to fit. #other options are reject parameters, number of iteractions #and/or convergernce limit -coef = ccdproc.iterfit(x, y, function='polynomial', order=3) +g = models.Gaussian1D(amplitude=1.2, mean=0.9, stddev=0.5) +coef = ccdproc.iterfit(x, y, model=g) + #fit a 2-D function with iterative rejections and the ability to #select different functions to fit. #other options are reject parameters, number of iteractions #and/or convergernce limit -coef = ccdproc.iterfit(data, function='polynomial', order=3) +p = models.Polynomial2D(degree=2) +coef = ccdproc.iterfit(data, function=p) #in addition to these operations, basic addition, subtraction # multiplication, and division should work for CCDDATA objects ccddata = ccddata + ccddata ccddata2 = ccddata * 2 - +# TODO: This needs changes in base classes to work but +# .add works for CCDDATA objects #combine a set of NDData objects alldata = ccdproc.combine([ccddata, ccddata2], method='average',