diff options
-rw-r--r-- | hstwcs/__init__.py | 6 | ||||
-rw-r--r-- | hstwcs/det2im.py | 140 | ||||
-rw-r--r-- | wcsutil/__init__.py | 14 |
3 files changed, 153 insertions, 7 deletions
diff --git a/hstwcs/__init__.py b/hstwcs/__init__.py index 8b83e1e..a843888 100644 --- a/hstwcs/__init__.py +++ b/hstwcs/__init__.py @@ -84,7 +84,8 @@ def makecorr(fname, allowed_corr): if 'DET2IMCorr' in allowed_corr: kw2update = det2im.DET2IMCorr.updateWCS(f) for kw in kw2update: - f[1].header.update(kw, kw2update[kw]) + #these kw are written to the primary header + f[0].header.update(kw, kw2update[kw]) for i in range(len(f))[1:]: # Perhaps all ext headers should be corrected (to be consistent) @@ -106,6 +107,9 @@ def makecorr(fname, allowed_corr): kw2update = dgeo.DGEOCorr.updateWCS(f) for kw in kw2update: f[1].header.update(kw, kw2update[kw]) + + + f.close() def getNrefchip(fobj): diff --git a/hstwcs/det2im.py b/hstwcs/det2im.py new file mode 100644 index 0000000..ec8e306 --- /dev/null +++ b/hstwcs/det2im.py @@ -0,0 +1,140 @@ +import pyfits +from pytools import fileutil +#from hstwcs.mappings import dgeo_vals +import numpy + + +class DET2IMCorr(object): + def updateWCS(cls, fobj): + """ + :Parameters: + `fobj`: pyfits object + Science file, for which a detector to image correction + is available + + """ + assert isinstance(fobj, pyfits.NP_pyfits.HDUList) + + d2imfile = fobj[0].header['D2IMFILE'] + axiscorr = cls.getAxisCorr(d2imfile) + if axiscorr == None: + new_kw = {} + else: + new_kw = {'D2IMEXT': d2imfile, 'AXISCORR': axiscorr} + cls.applyDet2ImCorr(fobj, axiscorr) + return new_kw + + updateWCS = classmethod(updateWCS) + + def getAxisCorr(cls, refname): + try: + direction = pyfits.getval(refname, ext=1, key='EXTNAME') + if direction == 'DX': return 1 + elif direction == 'DY': return 2 + else: + print '\tDET2IM correction expects the reference file to have' + print '\tan EXTNAME keyword of value "DX" or "DY"' + return None + except AttributeError: + print "\tAxis to which to apply the detector to image " + print "\tcorrection cannot be determined because the reference " + print "\tfile %s is missing a keyword EXTNAME" % refname + direction = None + return direction + getAxisCorr = classmethod(getAxisCorr) + + def applyDet2ImCorr(cls,fobj, axiscorr): + d2imfile = fileutil.osfn(fobj[0].header['D2IMFILE']) + hdu = cls.createDgeoHDU(d2imfile, axiscorr) + d2imarr_ind = cls.getD2imIndex(fobj) + if d2imarr_ind: + fobj[d2imarr_ind] = hdu + else: + fobj.append(hdu) + applyDet2ImCorr = classmethod(applyDet2ImCorr) + + def getD2imIndex(cls,fobj): + index = None + for e in range(len(fobj)): + try: + ename = fobj[e].header['EXTNAME'] + except KeyError: + continue + if ename == 'D2IMARR': + index = e + return index + getD2imIndex = classmethod(getD2imIndex) + + def createDgeoHDU(cls, d2imfile, axiscorr): + + d2im_data = pyfits.getdata(d2imfile, ext=1) + d2im_hdr = cls.createDet2ImHdr(d2im_data.shape, axiscorr) + hdu = pyfits.ImageHDU(header=d2im_hdr, data=d2im_data) + + return hdu + + createDgeoHDU = classmethod(createDgeoHDU) + + def createDet2ImHdr(cls, data_shape, axiscorr): + """ + Creates a header for the D2IMARR extension based on the + reference file recorded in D2IMFILE keyword in the primary header. + """ + + naxis1 = data_shape[0] + naxis2 = 0 + crpix1 = 0.0 + crpix2 = 0.0 + cdelt1 = 1.0 + cdelt2 = 1.0 + crval1 = 0.0 + crval2 = 0.0 + keys = ['XTENSION','BITPIX','NAXIS','NAXIS1','NAXIS2', + 'EXTNAME','EXTVER','PCOUNT','GCOUNT','CRPIX1', + 'CDELT1','CRVAL1','CRPIX2','CDELT2','CRVAL2', 'AXISCORR'] + + comments = {'XTENSION': 'Image extension', + 'BITPIX': 'IEEE floating point', + 'NAXIS': 'Number of axes', + 'NAXIS1': 'Number of image columns', + 'NAXIS2': 'Number of image rows', + 'EXTNAME': 'WCS distortion array', + 'EXTVER': 'Distortion array version number', + 'PCOUNT': 'Special data area of size 0', + 'GCOUNT': 'One data group', + 'CRPIX1': 'Distortion array reference pixel', + 'CDELT1': 'Grid step size in first coordinate', + 'CRVAL1': 'Image array pixel coordinate', + 'CRPIX2': 'Distortion array reference pixel', + 'CDELT2': 'Grid step size in second coordinate', + 'CRVAL2': 'Image array pixel coordinate', + 'AXISCORR': 'Direction in which to apply the det2im correction'} + + values = {'XTENSION': 'IMAGE', + 'BITPIX': -32, + 'NAXIS': 1, + 'NAXIS1': naxis1, + 'NAXIS2': naxis2, + 'EXTNAME': 'D2IMARR', + 'EXTVER': 1, + 'PCOUNT': 0, + 'GCOUNT': 1, + 'CRPIX1': crpix1, + 'CDELT1': cdelt1, + 'CRVAL1': crval1, + 'CRPIX2': crpix2, + 'CDELT2': cdelt2, + 'CRVAL2': crval2, + 'AXISCORR': axiscorr + } + + + cdl = pyfits.CardList() + for c in keys: + cdl.append(pyfits.Card(key=c, value=values[c], comment=comments[c])) + + hdr = pyfits.Header(cards=cdl) + return hdr + + createDet2ImHdr = classmethod(createDet2ImHdr) +
\ No newline at end of file diff --git a/wcsutil/__init__.py b/wcsutil/__init__.py index 05790f8..57e2460 100644 --- a/wcsutil/__init__.py +++ b/wcsutil/__init__.py @@ -116,6 +116,7 @@ class HSTWCS(WCS): self.dec_targ = primhdr.get('DEC_TARG', None) self.det2imfile = primhdr.get('D2IMFILE', None) self.det2imext = primhdr.get('D2IMEXT', None) + self.axiscorr = primhdr.get('AXISCORR', None) try: self.pav3 = primhdr['PA_V3'] @@ -260,10 +261,12 @@ class HSTWCS(WCS): def det2im(self, *args): - cpdis2 = self.get_d2im_lookup() + cpdis = self.get_d2im_lookup() d2im_wcs = WCS() - d2im_wcs.cpdis2 = cpdis2 - + if self.axiscorr == 1: + d2im_wcs.cpdis1 = cpdis + else: + d2im_wcs.cpdis2 = cpdis if len(args) == 2: xy, origin = args try: @@ -293,10 +296,9 @@ class HSTWCS(WCS): return [img[:, i] for i in range(img.shape[1])] def get_d2im_lookup(self): - d2im_data = pyfits.getdata(self.filename, ext=self.det2imext) - #d2im_data = d2im_data[:,N.newaxis] + d2im_data = pyfits.getdata(self.filename, ext=('D2IMARR', 1)) d2im_data = N.array([d2im_data]) - d2im_hdr = pyfits.getheader(self.filename, ext=self.det2imext) + d2im_hdr = pyfits.getheader(self.filename, ext=('D2IMARR', 1)) crpix = (d2im_hdr['CRPIX1'],d2im_hdr['CRPIX2']) crval = (d2im_hdr['CRVAL1'],d2im_hdr['CRVAL2']) |