summaryrefslogtreecommitdiff
path: root/stwcs/updatewcs
diff options
context:
space:
mode:
Diffstat (limited to 'stwcs/updatewcs')
-rw-r--r--stwcs/updatewcs/__init__.py116
-rw-r--r--stwcs/updatewcs/apply_corrections.py64
-rw-r--r--stwcs/updatewcs/corrections.py193
-rw-r--r--stwcs/updatewcs/det2im.py95
-rw-r--r--stwcs/updatewcs/makewcs.py145
-rw-r--r--stwcs/updatewcs/npol.py117
-rw-r--r--stwcs/updatewcs/utils.py57
-rw-r--r--stwcs/updatewcs/wfpc2_dgeo.py19
8 files changed, 420 insertions, 386 deletions
diff --git a/stwcs/updatewcs/__init__.py b/stwcs/updatewcs/__init__.py
index eb5e850..a59f106 100644
--- a/stwcs/updatewcs/__init__.py
+++ b/stwcs/updatewcs/__init__.py
@@ -1,14 +1,16 @@
-from __future__ import absolute_import, division, print_function # confidence high
+from __future__ import absolute_import, division, print_function
+import atexit
from astropy.io import fits
-from stwcs import wcsutil
-from stwcs.wcsutil import HSTWCS
-import stwcs
+from .. import wcsutil
+#from ..wcsutil.hwstwcs import HSTWCS
+
+from .. import __version__
from astropy import wcs as pywcs
import astropy
-from . import utils, corrections, makewcs
+from . import utils, corrections
from . import npol, det2im
from stsci.tools import parseinput, fileutil
from . import apply_corrections
@@ -17,10 +19,10 @@ import time
import logging
logger = logging.getLogger('stwcs.updatewcs')
-import atexit
atexit.register(logging.shutdown)
-#Note: The order of corrections is important
+# Note: The order of corrections is important
+
def updatewcs(input, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True,
checkfiles=True, verbose=False):
@@ -62,7 +64,7 @@ def updatewcs(input, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True,
geis and waiver fits files will be converted to MEF format.
Default value is True for standalone mode.
"""
- if verbose == False:
+ if not verbose:
logger.setLevel(100)
else:
formatter = logging.Formatter("%(asctime)s - %(name)s - %(levelname)s - %(message)s")
@@ -74,12 +76,12 @@ def updatewcs(input, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True,
logger.setLevel(verbose)
args = "vacorr=%s, tddcorr=%s, npolcorr=%s, d2imcorr=%s, checkfiles=%s, \
" % (str(vacorr), str(tddcorr), str(npolcorr),
- str(d2imcorr), str(checkfiles))
+ str(d2imcorr), str(checkfiles))
logger.info('\n\tStarting UPDATEWCS: %s', time.asctime())
files = parseinput.parseinput(input)[0]
logger.info("\n\tInput files: %s, " % [i for i in files])
- logger.info("\n\tInput arguments: %s" %args)
+ logger.info("\n\tInput arguments: %s" % args)
if checkfiles:
files = checkFiles(files)
if not files:
@@ -87,8 +89,8 @@ def updatewcs(input, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True,
return
for f in files:
- acorr = apply_corrections.setCorrections(f, vacorr=vacorr, \
- tddcorr=tddcorr,npolcorr=npolcorr, d2imcorr=d2imcorr)
+ acorr = apply_corrections.setCorrections(f, vacorr=vacorr, tddcorr=tddcorr,
+ npolcorr=npolcorr, d2imcorr=d2imcorr)
if 'MakeWCS' in acorr and newIDCTAB(f):
logger.warning("\n\tNew IDCTAB file detected. All current WCSs will be deleted")
cleanWCS(f)
@@ -97,6 +99,7 @@ def updatewcs(input, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True,
return files
+
def makecorr(fname, allowed_corr):
"""
Purpose
@@ -111,11 +114,11 @@ def makecorr(fname, allowed_corr):
"""
logger.info("Allowed corrections: {0}".format(allowed_corr))
f = fits.open(fname, mode='update')
- #Determine the reference chip and create the reference HSTWCS object
+ # Determine the reference chip and create the reference HSTWCS object
nrefchip, nrefext = getNrefchip(f)
wcsutil.restoreWCS(f, nrefext, wcskey='O')
- rwcs = HSTWCS(fobj=f, ext=nrefext)
- rwcs.readModel(update=True,header=f[nrefext].header)
+ rwcs = wcsutil.HSTWCS(fobj=f, ext=nrefext)
+ rwcs.readModel(update=True, header=f[nrefext].header)
if 'DET2IMCorr' in allowed_corr:
kw2update = det2im.DET2IMCorr.updateWCS(f)
@@ -127,16 +130,16 @@ def makecorr(fname, allowed_corr):
if 'extname' in extn.header:
extname = extn.header['extname'].lower()
- if extname == 'sci':
+ if extname == 'sci':
wcsutil.restoreWCS(f, ext=i, wcskey='O')
sciextver = extn.header['extver']
ref_wcs = rwcs.deepcopy()
hdr = extn.header
- ext_wcs = HSTWCS(fobj=f, ext=i)
- ### check if it exists first!!!
+ ext_wcs = wcsutil.HSTWCS(fobj=f, ext=i)
+ # check if it exists first!!!
# 'O ' can be safely archived again because it has been restored first.
wcsutil.archiveWCS(f, ext=i, wcskey="O", wcsname="OPUS", reusekey=True)
- ext_wcs.readModel(update=True,header=hdr)
+ ext_wcs.readModel(update=True, header=hdr)
for c in allowed_corr:
if c != 'NPOLCorr' and c != 'DET2IMCorr':
corr_klass = corrections.__getattribute__(c)
@@ -147,14 +150,14 @@ def makecorr(fname, allowed_corr):
idcname = f[0].header.get('IDCTAB', " ")
if idcname.strip() and 'idc.fits' in idcname:
wname = ''.join(['IDC_',
- utils.extract_rootname(idcname,suffix='_idc')])
+ utils.extract_rootname(idcname, suffix='_idc')])
else: wname = " "
hdr['WCSNAME'] = wname
elif extname in ['err', 'dq', 'sdq', 'samp', 'time']:
cextver = extn.header['extver']
if cextver == sciextver:
- hdr = f[('SCI',sciextver)].header
+ hdr = f[('SCI', sciextver)].header
w = pywcs.WCS(hdr, f)
copyWCS(w, extn.header)
@@ -167,14 +170,14 @@ def makecorr(fname, allowed_corr):
f[1].header[kw] = kw2update[kw]
# Finally record the version of the software which updated the WCS
if 'HISTORY' in f[0].header:
- f[0].header.set('UPWCSVER', value=stwcs.__version__,
+ f[0].header.set('UPWCSVER', value=__version__,
comment="Version of STWCS used to updated the WCS",
before='HISTORY')
f[0].header.set('PYWCSVER', value=astropy.__version__,
comment="Version of PYWCS used to updated the WCS",
before='HISTORY')
elif 'ASN_MTYP' in f[0].header:
- f[0].header.set('UPWCSVER', value=stwcs.__version__,
+ f[0].header.set('UPWCSVER', value=__version__,
comment="Version of STWCS used to updated the WCS",
after='ASN_MTYP')
f[0].header.set('PYWCSVER', value=astropy.__version__,
@@ -185,20 +188,21 @@ def makecorr(fname, allowed_corr):
for i in range(len(f[0].header) - 1, 0, -1):
if f[0].header[i].strip() != '':
break
- f[0].header.set('UPWCSVER', stwcs.__version__,
+ f[0].header.set('UPWCSVER', __version__,
"Version of STWCS used to updated the WCS",
after=i)
f[0].header.set('PYWCSVER', astropy.__version__,
"Version of PYWCS used to updated the WCS",
after=i)
# add additional keywords to be used by headerlets
- distdict = utils.construct_distname(f,rwcs)
+ distdict = utils.construct_distname(f, rwcs)
f[0].header['DISTNAME'] = distdict['DISTNAME']
f[0].header['SIPNAME'] = distdict['SIPNAME']
# Make sure NEXTEND keyword remains accurate
- f[0].header['NEXTEND'] = len(f)-1
+ f[0].header['NEXTEND'] = len(f) - 1
f.close()
+
def copyWCS(w, ehdr):
"""
This is a convenience function to copy a WCS object
@@ -214,6 +218,7 @@ def copyWCS(w, ehdr):
key = k[:7]
ehdr[key] = hwcs[k]
+
def getNrefchip(fobj):
"""
@@ -235,15 +240,15 @@ def getNrefchip(fobj):
if instrument == 'WFPC2':
chipkw = 'DETECTOR'
- extvers = [("SCI",img.header['EXTVER']) for img in
- fobj[1:] if img.header['EXTNAME'].lower()=='sci']
+ extvers = [("SCI", img.header['EXTVER']) for img in
+ fobj[1:] if img.header['EXTNAME'].lower() == 'sci']
detectors = [img.header[chipkw] for img in fobj[1:] if
- img.header['EXTNAME'].lower()=='sci']
+ img.header['EXTNAME'].lower() == 'sci']
fitsext = [i for i in range(len(fobj))[1:] if
- fobj[i].header['EXTNAME'].lower()=='sci']
- det2ext=dict(list(zip(detectors, extvers)))
- ext2det=dict(list(zip(extvers, detectors)))
- ext2fitsext=dict(list(zip(extvers, fitsext)))
+ fobj[i].header['EXTNAME'].lower() == 'sci']
+ det2ext = dict(list(zip(detectors, extvers)))
+ ext2det = dict(list(zip(extvers, detectors)))
+ ext2fitsext = dict(list(zip(extvers, fitsext)))
if 3 not in detectors:
nrefchip = ext2det.pop(extvers[0])
@@ -256,15 +261,15 @@ def getNrefchip(fobj):
elif (instrument == 'ACS' and fobj[0].header['DETECTOR'] == 'WFC') or \
(instrument == 'WFC3' and fobj[0].header['DETECTOR'] == 'UVIS'):
chipkw = 'CCDCHIP'
- extvers = [("SCI",img.header['EXTVER']) for img in
- fobj[1:] if img.header['EXTNAME'].lower()=='sci']
+ extvers = [("SCI", img.header['EXTVER']) for img in
+ fobj[1:] if img.header['EXTNAME'].lower() == 'sci']
detectors = [img.header[chipkw] for img in fobj[1:] if
- img.header['EXTNAME'].lower()=='sci']
+ img.header['EXTNAME'].lower() == 'sci']
fitsext = [i for i in range(len(fobj))[1:] if
- fobj[i].header['EXTNAME'].lower()=='sci']
- det2ext=dict(list(zip(detectors, extvers)))
- ext2det=dict(list(zip(extvers, detectors)))
- ext2fitsext=dict(list(zip(extvers, fitsext)))
+ fobj[i].header['EXTNAME'].lower() == 'sci']
+ det2ext = dict(list(zip(detectors, extvers)))
+ ext2det = dict(list(zip(extvers, detectors)))
+ ext2fitsext = dict(list(zip(extvers, fitsext)))
if 2 not in detectors:
nrefchip = ext2det.pop(extvers[0])
@@ -276,12 +281,13 @@ def getNrefchip(fobj):
else:
for i in range(len(fobj)):
extname = fobj[i].header.get('EXTNAME', None)
- if extname != None and extname.lower == 'sci':
+ if extname is not None and extname.lower == 'sci':
nrefext = i
break
return nrefchip, nrefext
+
def checkFiles(input):
"""
Checks that input files are in the correct format.
@@ -292,13 +298,14 @@ def checkFiles(input):
removed_files = []
newfiles = []
if not isinstance(input, list):
- input=[input]
+ input = [input]
for file in input:
try:
- imgfits,imgtype = fileutil.isFits(file)
+ imgfits, imgtype = fileutil.isFits(file)
except IOError:
- logger.warning( "\n\tFile %s could not be found, removing it from the input list.\n" %file)
+ logger.warning("File {0} could not be found, removing it from the"
+ "input list.".format(file))
removed_files.append(file)
continue
# Check for existence of waiver FITS input, and quit if found.
@@ -306,8 +313,9 @@ def checkFiles(input):
if imgfits:
if imgtype == 'waiver':
newfilename = waiver2mef(file, convert_dq=True)
- if newfilename == None:
- logger.warning("\n\tRemoving file %s from input list - could not convert waiver to mef" %file)
+ if newfilename is None:
+ logger.warning("Could not convert waiver to mef."
+ "Removing file {0} from input list".format(file))
removed_files.append(file)
else:
newfiles.append(newfilename)
@@ -319,24 +327,26 @@ def checkFiles(input):
# Convert the corresponding data quality file if present
if not imgfits:
newfilename = geis2mef(file, convert_dq=True)
- if newfilename == None:
- logger.warning("\n\tRemoving file %s from input list - could not convert geis to mef" %file)
+ if newfilename is None:
+ logger.warning("Could not convert file {0} from geis to mef, removing it from input list".format(file))
removed_files.append(file)
else:
newfiles.append(newfilename)
if removed_files:
- logger.warning('\n\tThe following files will be removed from the list of files to be processed %s' % removed_files)
+ logger.warning('\n\tThe following files will be removed from the list of files',
+ 'to be processed %s' % removed_files)
newfiles = checkFiles(newfiles)[0]
logger.info("\n\tThese files passed the input check and will be processed: %s" % newfiles)
return newfiles
+
def newIDCTAB(fname):
- #When this is called we know there's a kw IDCTAB in the header
+ # When this is called we know there's a kw IDCTAB in the header
hdul = fits.open(fname)
idctab = fileutil.osfn(hdul[0].header['IDCTAB'])
try:
- #check for the presence of IDCTAB in the first extension
+ # check for the presence of IDCTAB in the first extension
oldidctab = fileutil.osfn(hdul[1].header['IDCTAB'])
except KeyError:
return False
@@ -345,6 +355,7 @@ def newIDCTAB(fname):
else:
return True
+
def cleanWCS(fname):
# A new IDCTAB means all previously computed WCS's are invalid
# We are deleting all of them except the original OPUS WCS.
@@ -363,6 +374,7 @@ def cleanWCS(fname):
# Some extensions don't have the alternate (or any) WCS keywords
continue
+
def getCorrections(instrument):
"""
Print corrections available for an instrument
@@ -373,4 +385,4 @@ def getCorrections(instrument):
acorr = apply_corrections.allowed_corrections[instrument]
print("The following corrections will be performed for instrument %s\n" % instrument)
- for c in acorr: print(c,': ' , apply_corrections.cnames[c])
+ for c in acorr: print(c, ': ', apply_corrections.cnames[c])
diff --git a/stwcs/updatewcs/apply_corrections.py b/stwcs/updatewcs/apply_corrections.py
index 86b623a..4a8e5ac 100644
--- a/stwcs/updatewcs/apply_corrections.py
+++ b/stwcs/updatewcs/apply_corrections.py
@@ -1,37 +1,35 @@
-from __future__ import division, print_function # confidence high
+from __future__ import absolute_import, division, print_function
-import os
+import os.path
from astropy.io import fits
-import time
from stsci.tools import fileutil
-import os.path
-from stwcs.wcsutil import altwcs
from . import utils
from . import wfpc2_dgeo
import logging
logger = logging.getLogger("stwcs.updatewcs.apply_corrections")
-#Note: The order of corrections is important
+# Note: The order of corrections is important
# A dictionary which lists the allowed corrections for each instrument.
# These are the default corrections applied also in the pipeline.
-allowed_corrections={'WFPC2': ['DET2IMCorr', 'MakeWCS','CompSIP', 'VACorr'],
- 'ACS': ['DET2IMCorr', 'TDDCorr', 'MakeWCS', 'CompSIP','VACorr', 'NPOLCorr'],
- 'STIS': ['MakeWCS', 'CompSIP','VACorr'],
- 'NICMOS': ['MakeWCS', 'CompSIP','VACorr'],
- 'WFC3': ['DET2IMCorr', 'MakeWCS', 'CompSIP','VACorr', 'NPOLCorr'],
- }
+allowed_corrections = {'WFPC2': ['DET2IMCorr', 'MakeWCS', 'CompSIP', 'VACorr'],
+ 'ACS': ['DET2IMCorr', 'TDDCorr', 'MakeWCS', 'CompSIP', 'VACorr', 'NPOLCorr'],
+ 'STIS': ['MakeWCS', 'CompSIP', 'VACorr'],
+ 'NICMOS': ['MakeWCS', 'CompSIP', 'VACorr'],
+ 'WFC3': ['DET2IMCorr', 'MakeWCS', 'CompSIP', 'VACorr', 'NPOLCorr'],
+ }
cnames = {'DET2IMCorr': 'Detector to Image Correction',
- 'TDDCorr': 'Time Dependent Distortion Correction',
- 'MakeWCS': 'Recalculate basic WCS keywords based on the distortion model',
- 'CompSIP': 'Given IDCTAB distortion model calculate the SIP coefficients',
- 'VACorr': 'Velocity Aberration Correction',
- 'NPOLCorr': 'Lookup Table Distortion'
- }
+ 'TDDCorr': 'Time Dependent Distortion Correction',
+ 'MakeWCS': 'Recalculate basic WCS keywords based on the distortion model',
+ 'CompSIP': 'Given IDCTAB distortion model calculate the SIP coefficients',
+ 'VACorr': 'Velocity Aberration Correction',
+ 'NPOLCorr': 'Lookup Table Distortion'
+ }
+
def setCorrections(fname, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=True):
"""
@@ -40,7 +38,7 @@ def setCorrections(fname, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=Tru
for the instrument.
"""
instrument = fits.getval(fname, 'INSTRUME')
- # make a copy of this list !
+ # make a copy of this list
acorr = allowed_corrections[instrument][:]
# For WFPC2 images, the old-style DGEOFILE needs to be
@@ -56,18 +54,22 @@ def setCorrections(fname, vacorr=True, tddcorr=True, npolcorr=True, d2imcorr=Tru
if 'MakeWCS' in acorr: acorr.remove('MakeWCS')
if 'CompSIP' in acorr: acorr.remove('CompSIP')
- if 'VACorr' in acorr and vacorr==False: acorr.remove('VACorr')
+ if 'VACorr' in acorr and not vacorr:
+ acorr.remove('VACorr')
if 'TDDCorr' in acorr:
tddcorr = applyTDDCorr(fname, tddcorr)
- if tddcorr == False: acorr.remove('TDDCorr')
+ if not tddcorr:
+ acorr.remove('TDDCorr')
if 'NPOLCorr' in acorr:
npolcorr = applyNpolCorr(fname, npolcorr)
- if npolcorr == False: acorr.remove('NPOLCorr')
+ if not npolcorr:
+ acorr.remove('NPOLCorr')
if 'DET2IMCorr' in acorr:
d2imcorr = applyD2ImCorr(fname, d2imcorr)
- if d2imcorr == False: acorr.remove('DET2IMCorr')
- logger.info("\n\tCorrections to be applied to %s: %s " % (fname, str(acorr)))
+ if not d2imcorr:
+ acorr.remove('DET2IMCorr')
+ logger.info("Corrections to be applied to {0} {1}".format(fname, acorr))
return acorr
@@ -119,22 +121,21 @@ def applyTDDCorr(fname, utddcorr):
except KeyError:
tddswitch = 'PERFORM'
- if instrument == 'ACS' and detector == 'WFC' and utddcorr == True and tddswitch == 'PERFORM':
+ if instrument == 'ACS' and detector == 'WFC' and utddcorr and tddswitch == 'PERFORM':
tddcorr = True
try:
idctab = phdr['IDCTAB']
except KeyError:
tddcorr = False
- #print "***IDCTAB keyword not found - not applying TDD correction***\n"
if os.path.exists(fileutil.osfn(idctab)):
tddcorr = True
else:
tddcorr = False
- #print "***IDCTAB file not found - not applying TDD correction***\n"
else:
tddcorr = False
return tddcorr
+
def applyNpolCorr(fname, unpolcorr):
"""
Determines whether non-polynomial distortion lookup tables should be added
@@ -156,9 +157,8 @@ def applyNpolCorr(fname, unpolcorr):
return False
fnpol0 = fileutil.osfn(fnpol0)
if not fileutil.findFile(fnpol0):
- msg = """\n\tKw "NPOLFILE" exists in primary header but file %s not found
- Non-polynomial distortion correction will not be applied\n
- """ % fnpol0
+ msg = '"NPOLFILE" exists in primary header but file {0} not found.'
+ 'Non-polynomial distortion correction will not be applied.'.format(file)
logger.critical(msg)
raise IOError("NPOLFILE {0} not found".format(fnpol0))
try:
@@ -176,7 +176,7 @@ def applyNpolCorr(fname, unpolcorr):
else:
# npl file defined in first extension may not be found
# but if a valid kw exists in the primary header, non-polynomial
- #distortion correction should be applied.
+ # distortion correction should be applied.
applyNPOLCorr = True
except KeyError:
# the case of "NPOLFILE" kw present in primary header but "NPOLEXT" missing
@@ -191,6 +191,7 @@ def applyNpolCorr(fname, unpolcorr):
applyNPOLCorr = False
return (applyNPOLCorr and unpolcorr)
+
def isOldStyleDGEO(fname, dgname):
# checks if the file defined in a NPOLFILE kw is a full size
# (old style) image
@@ -209,6 +210,7 @@ def isOldStyleDGEO(fname, dgname):
else:
return False
+
def applyD2ImCorr(fname, d2imcorr):
applyD2IMCorr = True
try:
diff --git a/stwcs/updatewcs/corrections.py b/stwcs/updatewcs/corrections.py
index d3641eb..975a806 100644
--- a/stwcs/updatewcs/corrections.py
+++ b/stwcs/updatewcs/corrections.py
@@ -1,8 +1,9 @@
-from __future__ import division, print_function # confidence high
+from __future__ import absolute_import, division, print_function
import copy
import datetime
-import logging, time
+import logging
+import time
import numpy as np
from numpy import linalg
from stsci.tools import fileutil
@@ -11,11 +12,12 @@ from . import npol
from . import makewcs
from .utils import diff_angles
-logger=logging.getLogger('stwcs.updatewcs.corrections')
+logger = logging.getLogger('stwcs.updatewcs.corrections')
MakeWCS = makewcs.MakeWCS
NPOLCorr = npol.NPOLCorr
+
class TDDCorr(object):
"""
Apply time dependent distortion correction to distortion coefficients and basic
@@ -70,39 +72,39 @@ class TDDCorr(object):
ext_wcs.idcmodel.ocx = copy.deepcopy(ext_wcs.idcmodel.cx)
ext_wcs.idcmodel.ocy = copy.deepcopy(ext_wcs.idcmodel.cy)
- newkw = {'TDDALPHA': None, 'TDDBETA':None,
- 'OCX10':ext_wcs.idcmodel.ocx[1,0],
- 'OCX11':ext_wcs.idcmodel.ocx[1,1],
- 'OCY10':ext_wcs.idcmodel.ocy[1,0],
- 'OCY11':ext_wcs.idcmodel.ocy[1,1],
- 'TDD_CTA':None, 'TDD_CTB':None,
- 'TDD_CYA':None, 'TDD_CYB':None,
- 'TDD_CXA':None, 'TDD_CXB':None}
-
- if ext_wcs.idcmodel.refpix['skew_coeffs']is not None and \
- ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTB'] is not None:
+ newkw = {'TDDALPHA': None, 'TDDBETA': None,
+ 'OCX10': ext_wcs.idcmodel.ocx[1, 0],
+ 'OCX11': ext_wcs.idcmodel.ocx[1, 1],
+ 'OCY10': ext_wcs.idcmodel.ocy[1, 0],
+ 'OCY11': ext_wcs.idcmodel.ocy[1, 1],
+ 'TDD_CTA': None, 'TDD_CTB': None,
+ 'TDD_CYA': None, 'TDD_CYB': None,
+ 'TDD_CXA': None, 'TDD_CXB': None
+ }
+
+ if ext_wcs.idcmodel.refpix['skew_coeffs'] is not None and \
+ ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTB'] is not None:
cls.apply_tdd2idc2015(ref_wcs)
cls.apply_tdd2idc2015(ext_wcs)
- newkw.update({
- 'TDD_CTA':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTA'],
- 'TDD_CYA':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CYA'],
- 'TDD_CXA':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CXA'],
- 'TDD_CTB':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTB'],
- 'TDD_CYB':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CYB'],
- 'TDD_CXB':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CXB']})
-
- elif ext_wcs.idcmodel.refpix['skew_coeffs']is not None and \
- ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_BETA'] is not None:
- logger.info("\n\t Applying 2014-calibrated TDD: %s" % time.asctime())
+ newkw.update({'TDD_CTA': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTA'],
+ 'TDD_CYA': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CYA'],
+ 'TDD_CXA': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CXA'],
+ 'TDD_CTB': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CTB'],
+ 'TDD_CYB': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CYB'],
+ 'TDD_CXB': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CXB']
+ })
+
+ elif ext_wcs.idcmodel.refpix['skew_coeffs'] is not None and \
+ ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_BETA'] is not None:
+ logger.info("Applying 2014-calibrated TDD: {0}".format(time.asctime()))
# We have 2014-calibrated TDD, not J.A.-style TDD
cls.apply_tdd2idc2(ref_wcs)
cls.apply_tdd2idc2(ext_wcs)
- newkw.update({'TDD_CYA':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_ALPHA'],
- 'TDD_CYB':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_BETA'],
- 'TDD_CXA':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CX_ALPHA'],
- 'TDD_CXB':ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CX_BETA']})
-
+ newkw.update({'TDD_CYA': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_ALPHA'],
+ 'TDD_CYB': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CY_BETA'],
+ 'TDD_CXA': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CX_ALPHA'],
+ 'TDD_CXB': ext_wcs.idcmodel.refpix['skew_coeffs']['TDD_CX_BETA']})
else:
alpha, beta = cls.compute_alpha_beta(ext_wcs)
cls.apply_tdd2idc(ref_wcs, alpha, beta)
@@ -111,8 +113,7 @@ class TDDCorr(object):
ext_wcs.idcmodel.refpix['TDDBETA'] = beta
ref_wcs.idcmodel.refpix['TDDALPHA'] = alpha
ref_wcs.idcmodel.refpix['TDDBETA'] = beta
-
- newkw.update( {'TDDALPHA': alpha, 'TDDBETA':beta} )
+ newkw.update({'TDDALPHA': alpha, 'TDDBETA': beta} )
return newkw
updateWCS = classmethod(updateWCS)
@@ -121,10 +122,10 @@ class TDDCorr(object):
""" Applies 2015-calibrated TDD correction to a couple of IDCTAB
coefficients for ACS/WFC observations.
"""
- if not isinstance(hwcs.date_obs,float):
- year,month,day = hwcs.date_obs.split('-')
- rdate = datetime.datetime(int(year),int(month),int(day))
- rday = float(rdate.strftime("%j"))/365.25 + rdate.year
+ if not isinstance(hwcs.date_obs, float):
+ year, month, day = hwcs.date_obs.split('-')
+ rdate = datetime.datetime(int(year), int(month), int(day))
+ rday = float(rdate.strftime("%j")) / 365.25 + rdate.year
else:
rday = hwcs.date_obs
@@ -132,12 +133,11 @@ class TDDCorr(object):
delta_date = rday - skew_coeffs['TDD_DATE']
if skew_coeffs['TDD_CXB'] is not None:
- hwcs.idcmodel.cx[1,1] += skew_coeffs['TDD_CXB']*delta_date
+ hwcs.idcmodel.cx[1, 1] += skew_coeffs['TDD_CXB'] * delta_date
if skew_coeffs['TDD_CTB'] is not None:
- hwcs.idcmodel.cy[1,1] += skew_coeffs['TDD_CTB']*delta_date
+ hwcs.idcmodel.cy[1, 1] += skew_coeffs['TDD_CTB'] * delta_date
if skew_coeffs['TDD_CYB'] is not None:
- hwcs.idcmodel.cy[1,0] += skew_coeffs['TDD_CYB']*delta_date
- #print("CX[1,1]_TDD={}, CY[1,1]_TDD={}, CY[1,0]_TDD={}".format(hwcs.idcmodel.cx[1,1],hwcs.idcmodel.cy[1,1],hwcs.idcmodel.cy[1,0]))
+ hwcs.idcmodel.cy[1, 0] += skew_coeffs['TDD_CYB'] * delta_date
apply_tdd2idc2015 = classmethod(apply_tdd2idc2015)
@@ -145,10 +145,10 @@ class TDDCorr(object):
""" Applies 2014-calibrated TDD correction to single IDCTAB coefficient
of an ACS/WFC observation.
"""
- if not isinstance(hwcs.date_obs,float):
- year,month,day = hwcs.date_obs.split('-')
- rdate = datetime.datetime(int(year),int(month),int(day))
- rday = float(rdate.strftime("%j"))/365.25 + rdate.year
+ if not isinstance(hwcs.date_obs, float):
+ year, month, day = hwcs.date_obs.split('-')
+ rdate = datetime.datetime(int(year), int(month), int(day))
+ rday = float(rdate.strftime("%j")) / 365.25 + rdate.year
else:
rday = hwcs.date_obs
@@ -156,21 +156,24 @@ class TDDCorr(object):
cy_beta = skew_coeffs['TDD_CY_BETA']
cy_alpha = skew_coeffs['TDD_CY_ALPHA']
delta_date = rday - skew_coeffs['TDD_DATE']
- print("DELTA_DATE: {0} based on rday: {1}, TDD_DATE: {2}".format(delta_date,rday,skew_coeffs['TDD_DATE']))
+ logger.info("DELTA_DATE: {0} based on rday: {1}, TDD_DATE: {2}".format(delta_date, rday,
+ skew_coeffs['TDD_DATE']))
if cy_alpha is None:
- hwcs.idcmodel.cy[1,1] += cy_beta*delta_date
+ hwcs.idcmodel.cy[1, 1] += cy_beta * delta_date
else:
- new_beta = cy_alpha + cy_beta*delta_date
- hwcs.idcmodel.cy[1,1] = new_beta
- print("CY11: {0} based on alpha: {1}, beta: {2}".format(hwcs.idcmodel.cy[1,1],cy_alpha,cy_beta))
-
+ new_beta = cy_alpha + cy_beta * delta_date
+ hwcs.idcmodel.cy[1, 1] = new_beta
+ logger.info("CY11: {0} based on alpha: {1}, beta: {2}".format(hwcs.idcmodel.cy[1, 1],
+ cy_alpha, cy_beta))
+
cx_beta = skew_coeffs['TDD_CX_BETA']
cx_alpha = skew_coeffs['TDD_CX_ALPHA']
if cx_alpha is not None:
- new_beta = cx_alpha + cx_beta*delta_date
- hwcs.idcmodel.cx[1,1] = new_beta
- print("CX11: {0} based on alpha: {1}, beta: {2}".format(new_beta,cx_alpha,cx_beta))
+ new_beta = cx_alpha + cx_beta * delta_date
+ hwcs.idcmodel.cx[1, 1] = new_beta
+ logger.info("CX11: {0} based on alpha: {1}, beta: {2}".format(new_beta,
+ cx_alpha, cx_beta))
apply_tdd2idc2 = classmethod(apply_tdd2idc2)
@@ -182,11 +185,12 @@ class TDDCorr(object):
theta_v2v3 = 2.234529
mrotp = fileutil.buildRotMatrix(theta_v2v3)
mrotn = fileutil.buildRotMatrix(-theta_v2v3)
- tdd_mat = np.array([[1+(beta/2048.), alpha/2048.],[alpha/2048.,1-(beta/2048.)]],np.float64)
+ tdd_mat = np.array([[1 + (beta / 2048.), alpha / 2048.],
+ [alpha / 2048., 1 - (beta / 2048.)]], np.float64)
abmat1 = np.dot(tdd_mat, mrotn)
- abmat2 = np.dot(mrotp,abmat1)
+ abmat2 = np.dot(mrotp, abmat1)
xshape, yshape = hwcs.idcmodel.cx.shape, hwcs.idcmodel.cy.shape
- icxy = np.dot(abmat2,[hwcs.idcmodel.cx.ravel(), hwcs.idcmodel.cy.ravel()])
+ icxy = np.dot(abmat2, [hwcs.idcmodel.cx.ravel(), hwcs.idcmodel.cy.ravel()])
hwcs.idcmodel.cx = icxy[0]
hwcs.idcmodel.cy = icxy[1]
hwcs.idcmodel.cx.shape = xshape
@@ -209,10 +213,10 @@ class TDDCorr(object):
alpha = 0.095 + 0.090*(rday-dday)/2.5
beta = -0.029 - 0.030*(rday-dday)/2.5
"""
- if not isinstance(ext_wcs.date_obs,float):
- year,month,day = ext_wcs.date_obs.split('-')
- rdate = datetime.datetime(int(year),int(month),int(day))
- rday = float(rdate.strftime("%j"))/365.25 + rdate.year
+ if not isinstance(ext_wcs.date_obs, float):
+ year, month, day = ext_wcs.date_obs.split('-')
+ rdate = datetime.datetime(int(year), int(month), int(day))
+ rday = float(rdate.strftime("%j")) / 365.25 + rdate.year
else:
rday = ext_wcs.date_obs
@@ -220,26 +224,27 @@ class TDDCorr(object):
if skew_coeffs is None:
# Only print out warning for post-SM4 data where this may matter
if rday > 2009.0:
- err_str = "------------------------------------------------------------------------ \n"
+ err_str = "------------------------------------------------------------------------ \n"
err_str += "WARNING: the IDCTAB geometric distortion file specified in the image \n"
err_str += " header did not have the time-dependent distortion coefficients. \n"
err_str += " The pre-SM4 time-dependent skew solution will be used by default.\n"
err_str += " Please update IDCTAB with new reference file from HST archive. \n"
- err_str += "------------------------------------------------------------------------ \n"
+ err_str += "------------------------------------------------------------------------ \n"
print(err_str)
# Using default pre-SM4 coefficients
- skew_coeffs = {'TDD_A':[0.095,0.090/2.5],
- 'TDD_B':[-0.029,-0.030/2.5],
- 'TDD_DATE':2004.5,'TDDORDER':1}
+ skew_coeffs = {'TDD_A': [0.095, 0.090 / 2.5],
+ 'TDD_B': [-0.029, -0.030 / 2.5],
+ 'TDD_DATE': 2004.5,
+ 'TDDORDER': 1}
alpha = 0
beta = 0
# Compute skew terms, allowing for non-linear coefficients as well
- for c in range(skew_coeffs['TDDORDER']+1):
- alpha += skew_coeffs['TDD_A'][c]* np.power((rday-skew_coeffs['TDD_DATE']),c)
- beta += skew_coeffs['TDD_B'][c]*np.power((rday-skew_coeffs['TDD_DATE']),c)
+ for c in range(skew_coeffs['TDDORDER'] + 1):
+ alpha += skew_coeffs['TDD_A'][c] * np.power((rday - skew_coeffs['TDD_DATE']), c)
+ beta += skew_coeffs['TDD_B'][c] * np.power((rday - skew_coeffs['TDD_DATE']), c)
- return alpha,beta
+ return alpha, beta
compute_alpha_beta = classmethod(compute_alpha_beta)
@@ -254,21 +259,21 @@ class VACorr(object):
"""
def updateWCS(cls, ext_wcs, ref_wcs):
- logger.info("\n\tStarting VACorr: %s" % time.asctime())
+ logger.info("Starting VACorr: %s" % time.asctime())
if ext_wcs.vafactor != 1:
ext_wcs.wcs.cd = ext_wcs.wcs.cd * ext_wcs.vafactor
- crval0 = ref_wcs.wcs.crval[0] + ext_wcs.vafactor*diff_angles(ext_wcs.wcs.crval[0],
- ref_wcs.wcs.crval[0])
- crval1 = ref_wcs.wcs.crval[1] + ext_wcs.vafactor*diff_angles(ext_wcs.wcs.crval[1],
- ref_wcs.wcs.crval[1])
- crval = np.array([crval0,crval1])
+ crval0 = ref_wcs.wcs.crval[0] + ext_wcs.vafactor * diff_angles(ext_wcs.wcs.crval[0],
+ ref_wcs.wcs.crval[0])
+ crval1 = ref_wcs.wcs.crval[1] + ext_wcs.vafactor * diff_angles(ext_wcs.wcs.crval[1],
+ ref_wcs.wcs.crval[1])
+ crval = np.array([crval0, crval1])
ext_wcs.wcs.crval = crval
ext_wcs.wcs.set()
else:
pass
- kw2update={'CD1_1': ext_wcs.wcs.cd[0,0], 'CD1_2':ext_wcs.wcs.cd[0,1],
- 'CD2_1':ext_wcs.wcs.cd[1,0], 'CD2_2':ext_wcs.wcs.cd[1,1],
- 'CRVAL1':ext_wcs.wcs.crval[0], 'CRVAL2':ext_wcs.wcs.crval[1]}
+ kw2update = {'CD1_1': ext_wcs.wcs.cd[0, 0], 'CD1_2': ext_wcs.wcs.cd[0, 1],
+ 'CD2_1': ext_wcs.wcs.cd[1, 0], 'CD2_2': ext_wcs.wcs.cd[1, 1],
+ 'CRVAL1': ext_wcs.wcs.crval[0], 'CRVAL2': ext_wcs.wcs.crval[1]}
return kw2update
updateWCS = classmethod(updateWCS)
@@ -291,34 +296,34 @@ class CompSIP(object):
"""
def updateWCS(cls, ext_wcs, ref_wcs):
- logger.info("\n\tStarting CompSIP: %s" %time.asctime())
+ logger.info("Starting CompSIP: {0}".format(time.asctime()))
kw2update = {}
if not ext_wcs.idcmodel:
- logger.info("IDC model not found, SIP coefficient will not be computed")
+ logger.info("IDC model not found, SIP coefficient will not be computed.")
return kw2update
order = ext_wcs.idcmodel.norder
kw2update['A_ORDER'] = order
kw2update['B_ORDER'] = order
- pscale = ext_wcs.idcmodel.refpix['PSCALE']
+ # pscale = ext_wcs.idcmodel.refpix['PSCALE']
cx = ext_wcs.idcmodel.cx
cy = ext_wcs.idcmodel.cy
- matr = np.array([[cx[1,1],cx[1,0]], [cy[1,1],cy[1,0]]], dtype=np.float64)
+ matr = np.array([[cx[1, 1], cx[1, 0]], [cy[1, 1], cy[1, 0]]], dtype=np.float64)
imatr = linalg.inv(matr)
- akeys1 = np.zeros((order+1,order+1), dtype=np.float64)
- bkeys1 = np.zeros((order+1,order+1), dtype=np.float64)
- for n in range(order+1):
- for m in range(order+1):
- if n >= m and n>=2:
- idcval = np.array([[cx[n,m]],[cy[n,m]]])
+ akeys1 = np.zeros((order + 1, order + 1), dtype=np.float64)
+ bkeys1 = np.zeros((order + 1, order + 1), dtype=np.float64)
+ for n in range(order + 1):
+ for m in range(order + 1):
+ if n >= m and n >= 2:
+ idcval = np.array([[cx[n, m]], [cy[n, m]]])
sipval = np.dot(imatr, idcval)
- akeys1[m,n-m] = sipval[0]
- bkeys1[m,n-m] = sipval[1]
- Akey="A_%d_%d" % (m,n-m)
- Bkey="B_%d_%d" % (m,n-m)
- kw2update[Akey] = sipval[0,0] * ext_wcs.binned
- kw2update[Bkey] = sipval[1,0] * ext_wcs.binned
+ akeys1[m, n - m] = sipval[0]
+ bkeys1[m, n - m] = sipval[1]
+ Akey = "A_%d_%d" % (m, n - m)
+ Bkey = "B_%d_%d" % (m, n - m)
+ kw2update[Akey] = sipval[0, 0] * ext_wcs.binned
+ kw2update[Bkey] = sipval[1, 0] * ext_wcs.binned
kw2update['CTYPE1'] = 'RA---TAN-SIP'
kw2update['CTYPE2'] = 'DEC--TAN-SIP'
return kw2update
diff --git a/stwcs/updatewcs/det2im.py b/stwcs/updatewcs/det2im.py
index bddb37b..bb4cd9c 100644
--- a/stwcs/updatewcs/det2im.py
+++ b/stwcs/updatewcs/det2im.py
@@ -1,14 +1,13 @@
-from __future__ import absolute_import, division # confidence high
+from __future__ import absolute_import, division, print_function
-import numpy as np
from astropy.io import fits
from stsci.tools import fileutil
-from . import utils
-
-import logging, time
+import logging
+import time
logger = logging.getLogger('stwcs.updatewcs.d2im')
+
class DET2IMCorr(object):
"""
Defines a Lookup table prior distortion correction as per WCS paper IV.
@@ -38,11 +37,11 @@ class DET2IMCorr(object):
Science file, for which a distortion correction in a NPOLFILE is available
"""
- logger.info("\n\tStarting DET2IM: %s" %time.asctime())
+ logger.info("Starting DET2IM: {0}".format(time.asctime()))
try:
assert isinstance(fobj, fits.HDUList)
except AssertionError:
- logger.exception('\n\tInput must be a fits.HDUList object')
+ logger.exception('Input must be a fits.HDUList object')
raise
cls.applyDet2ImCorr(fobj)
@@ -84,7 +83,7 @@ class DET2IMCorr(object):
hdu = cls.createD2ImHDU(header, d2imfile=d2imfile,
wdvarr_ver=d2im_num_ext,
d2im_extname=ename[0],
- data=ename[1],ccdchip=ccdchip)
+ data=ename[1], ccdchip=ccdchip)
if wcsdvarr_ind and d2im_num_ext in wcsdvarr_ind:
fobj[wcsdvarr_ind[d2im_num_ext]] = hdu
else:
@@ -108,7 +107,7 @@ class DET2IMCorr(object):
continue
if ename == 'D2IMARR':
wcsd[fobj[e].header['EXTVER']] = e
- logger.debug("A map of D2IMARR extensions %s" % wcsd)
+ logger.debug("A map of D2IMARR extensions {0}".format(wcsd))
return wcsd
getWCSIndex = classmethod(getWCSIndex)
@@ -118,17 +117,17 @@ class DET2IMCorr(object):
Adds kw to sci extension to define WCSDVARR lookup table extensions
"""
- if d2im_extname =='DX':
- j=1
+ if d2im_extname == 'DX':
+ j = 1
else:
- j=2
-
- d2imerror = 'D2IMERR%s' %j
- d2imdis = 'D2IMDIS%s' %j
- d2imext = 'D2IM%s.' %j + 'EXTVER'
- d2imnaxes = 'D2IM%s.' %j +'NAXES'
- d2imaxis1 = 'D2IM%s.' %j+'AXIS.1'
- d2imaxis2 = 'D2IM%s.' %j+'AXIS.2'
+ j = 2
+
+ d2imerror = 'D2IMERR%s' % j
+ d2imdis = 'D2IMDIS%s' % j
+ d2imext = 'D2IM%s.' % j + 'EXTVER'
+ d2imnaxes = 'D2IM%s.' % j + 'NAXES'
+ d2imaxis1 = 'D2IM%s.' % j + 'AXIS.1'
+ d2imaxis2 = 'D2IM%s.' % j + 'AXIS.2'
keys = [d2imerror, d2imdis, d2imext, d2imnaxes, d2imaxis1, d2imaxis2]
values = {d2imerror: error_val,
d2imdis: 'Lookup',
@@ -154,7 +153,7 @@ class DET2IMCorr(object):
addSciExtKw = classmethod(addSciExtKw)
- def getData(cls,d2imfile, ccdchip):
+ def getData(cls, d2imfile, ccdchip):
"""
Get the data arrays from the reference D2I files
Make sure 'CCDCHIP' in the npolfile matches "CCDCHIP' in the science file.
@@ -162,8 +161,8 @@ class DET2IMCorr(object):
xdata, ydata = (None, None)
d2im = fits.open(d2imfile)
for ext in d2im:
- d2imextname = ext.header.get('EXTNAME',"")
- d2imccdchip = ext.header.get('CCDCHIP',1)
+ d2imextname = ext.header.get('EXTNAME', "")
+ d2imccdchip = ext.header.get('CCDCHIP', 1)
if d2imextname == 'DX' and d2imccdchip == ccdchip:
xdata = ext.data.copy()
continue
@@ -177,7 +176,7 @@ class DET2IMCorr(object):
getData = classmethod(getData)
def createD2ImHDU(cls, sciheader, d2imfile=None, wdvarr_ver=1,
- d2im_extname=None,data = None, ccdchip=1):
+ d2im_extname=None, data=None, ccdchip=1):
"""
Creates an HDU to be added to the file object.
"""
@@ -216,26 +215,26 @@ class DET2IMCorr(object):
naxis = d2im[1].header['NAXIS']
ccdchip = d2imextname
- kw = { 'NAXIS': 'Size of the axis',
- 'CDELT': 'Coordinate increment along axis',
- 'CRPIX': 'Coordinate system reference pixel',
- 'CRVAL': 'Coordinate system value at reference pixel',
- }
+ kw = {'NAXIS': 'Size of the axis',
+ 'CDELT': 'Coordinate increment along axis',
+ 'CRPIX': 'Coordinate system reference pixel',
+ 'CRVAL': 'Coordinate system value at reference pixel',
+ }
kw_comm1 = {}
kw_val1 = {}
for key in kw.keys():
- for i in range(1, naxis+1):
+ for i in range(1, naxis + 1):
si = str(i)
- kw_comm1[key+si] = kw[key]
+ kw_comm1[key + si] = kw[key]
- for i in range(1, naxis+1):
+ for i in range(1, naxis + 1):
si = str(i)
- kw_val1['NAXIS'+si] = d2im_header.get('NAXIS'+si)
- kw_val1['CDELT'+si] = d2im_header.get('CDELT'+si, 1.0) * sciheader.get('LTM'+si+'_'+si, 1)
- kw_val1['CRPIX'+si] = d2im_header.get('CRPIX'+si, 0.0)
- kw_val1['CRVAL'+si] = (d2im_header.get('CRVAL'+si, 0.0) - \
- sciheader.get('LTV'+str(i), 0))
+ kw_val1['NAXIS' + si] = d2im_header.get('NAXIS' + si)
+ kw_val1['CDELT' + si] = d2im_header.get('CDELT' + si, 1.0) * sciheader.get('LTM' + si + '_' + si, 1)
+ kw_val1['CRPIX' + si] = d2im_header.get('CRPIX' + si, 0.0)
+ kw_val1['CRVAL' + si] = (d2im_header.get('CRVAL' + si, 0.0) -
+ sciheader.get('LTV' + str(i), 0))
kw_comm0 = {'XTENSION': 'Image extension',
'BITPIX': 'IEEE floating point',
'NAXIS': 'Number of axes',
@@ -245,15 +244,15 @@ class DET2IMCorr(object):
'GCOUNT': 'One data group',
}
- kw_val0 = { 'XTENSION': 'IMAGE',
- 'BITPIX': -32,
- 'NAXIS': naxis,
- 'EXTNAME': 'D2IMARR',
- 'EXTVER': wdvarr_ver,
- 'PCOUNT': 0,
- 'GCOUNT': 1,
- 'CCDCHIP': ccdchip,
- }
+ kw_val0 = {'XTENSION': 'IMAGE',
+ 'BITPIX': -32,
+ 'NAXIS': naxis,
+ 'EXTNAME': 'D2IMARR',
+ 'EXTVER': wdvarr_ver,
+ 'PCOUNT': 0,
+ 'GCOUNT': 1,
+ 'CCDCHIP': ccdchip,
+ }
cdl = []
for key in kw_comm0.keys():
cdl.append((key, kw_val0[key], kw_comm0[key]))
@@ -266,8 +265,8 @@ class DET2IMCorr(object):
for i, c in enumerate(d2im_phdr):
if c == 'FILENAME':
start_indx = i
- if c == '': # remove blanks from end of header
- end_indx = i+1
+ if c == '': # remove blanks from end of header
+ end_indx = i + 1
break
if start_indx >= 0:
for card in d2im_phdr.cards[start_indx:end_indx]:
@@ -287,7 +286,7 @@ class DET2IMCorr(object):
if fobj[0].header['INSTRUME'] == 'ACS' and fobj[0].header['DETECTOR'] == 'WFC':
ccdchip = fobj[extname, extver].header['CCDCHIP']
elif fobj[0].header['INSTRUME'] == 'WFC3' and fobj[0].header['DETECTOR'] == 'UVIS':
- ccdchip = fobj[extname, extver].header['CCDCHIP']
+ ccdchip = fobj[extname, extver].header['CCDCHIP']
elif fobj[0].header['INSTRUME'] == 'WFPC2':
ccdchip = fobj[extname, extver].header['DETECTOR']
elif fobj[0].header['INSTRUME'] == 'STIS':
diff --git a/stwcs/updatewcs/makewcs.py b/stwcs/updatewcs/makewcs.py
index 06c6f9c..4c6b3df 100644
--- a/stwcs/updatewcs/makewcs.py
+++ b/stwcs/updatewcs/makewcs.py
@@ -1,16 +1,16 @@
-from __future__ import absolute_import, division # confidence high
-
-import datetime
+from __future__ import absolute_import, division, print_function
import numpy as np
-import logging, time
-from math import sin, sqrt, pow, cos, asin, atan2,pi
+import logging
+import time
+from math import sin, sqrt, pow, cos, asin, atan2, pi
from stsci.tools import fileutil
from . import utils
logger = logging.getLogger(__name__)
+
class MakeWCS(object):
"""
Recompute basic WCS keywords based on PA_V3 and distortion model.
@@ -36,7 +36,8 @@ class MakeWCS(object):
- Time dependent distortion correction is applied to both chips' V2/V3 values.
"""
- tdd_xyref = {1: [2048, 3072], 2:[2048, 1024]}
+ tdd_xyref = {1: [2048, 3072], 2: [2048, 1024]}
+
def updateWCS(cls, ext_wcs, ref_wcs):
"""
recomputes the basic WCS kw
@@ -53,20 +54,20 @@ class MakeWCS(object):
cls.uprefwcs(ext_wcs, ref_wcs, rv23_corr, ltvoff, offshift)
cls.upextwcs(ext_wcs, ref_wcs, v23_corr, rv23_corr, ltvoff, offshift)
- kw2update = {'CD1_1': ext_wcs.wcs.cd[0,0],
- 'CD1_2': ext_wcs.wcs.cd[0,1],
- 'CD2_1': ext_wcs.wcs.cd[1,0],
- 'CD2_2': ext_wcs.wcs.cd[1,1],
- 'CRVAL1': ext_wcs.wcs.crval[0],
- 'CRVAL2': ext_wcs.wcs.crval[1],
- 'CRPIX1': ext_wcs.wcs.crpix[0],
- 'CRPIX2': ext_wcs.wcs.crpix[1],
- 'IDCTAB': ext_wcs.idctab,
- 'OCX10' : ext_wcs.idcmodel.cx[1,0],
- 'OCX11' : ext_wcs.idcmodel.cx[1,1],
- 'OCY10' : ext_wcs.idcmodel.cy[1,0],
- 'OCY11' : ext_wcs.idcmodel.cy[1,1]
- }
+ kw2update = {'CD1_1': ext_wcs.wcs.cd[0, 0],
+ 'CD1_2': ext_wcs.wcs.cd[0, 1],
+ 'CD2_1': ext_wcs.wcs.cd[1, 0],
+ 'CD2_2': ext_wcs.wcs.cd[1, 1],
+ 'CRVAL1': ext_wcs.wcs.crval[0],
+ 'CRVAL2': ext_wcs.wcs.crval[1],
+ 'CRPIX1': ext_wcs.wcs.crpix[0],
+ 'CRPIX2': ext_wcs.wcs.crpix[1],
+ 'IDCTAB': ext_wcs.idctab,
+ 'OCX10': ext_wcs.idcmodel.cx[1, 0],
+ 'OCX11': ext_wcs.idcmodel.cx[1, 1],
+ 'OCY10': ext_wcs.idcmodel.cy[1, 0],
+ 'OCY11': ext_wcs.idcmodel.cy[1, 1]
+ }
return kw2update
updateWCS = classmethod(updateWCS)
@@ -82,41 +83,42 @@ class MakeWCS(object):
if ltv1 != 0. or ltv2 != 0.:
offsetx = ext_wcs.wcs.crpix[0] - ltv1 - ext_wcs.idcmodel.refpix['XREF']
offsety = ext_wcs.wcs.crpix[1] - ltv2 - ext_wcs.idcmodel.refpix['YREF']
- ext_wcs.idcmodel.shift(offsetx,offsety)
+ ext_wcs.idcmodel.shift(offsetx, offsety)
fx, fy = ext_wcs.idcmodel.cx, ext_wcs.idcmodel.cy
ext_wcs.setPscale()
- tddscale = (ref_wcs.pscale/fx[1,1])
- v2 = ext_wcs.idcmodel.refpix['V2REF'] + v23_corr[0,0] * tddscale
- v3 = ext_wcs.idcmodel.refpix['V3REF'] - v23_corr[1,0] * tddscale
- v2ref = ref_wcs.idcmodel.refpix['V2REF'] + rv23_corr[0,0] * tddscale
- v3ref = ref_wcs.idcmodel.refpix['V3REF'] - rv23_corr[1,0] * tddscale
+ tddscale = (ref_wcs.pscale / fx[1, 1])
+ v2 = ext_wcs.idcmodel.refpix['V2REF'] + v23_corr[0, 0] * tddscale
+ v3 = ext_wcs.idcmodel.refpix['V3REF'] - v23_corr[1, 0] * tddscale
+ v2ref = ref_wcs.idcmodel.refpix['V2REF'] + rv23_corr[0, 0] * tddscale
+ v3ref = ref_wcs.idcmodel.refpix['V3REF'] - rv23_corr[1, 0] * tddscale
- R_scale = ref_wcs.idcmodel.refpix['PSCALE']/3600.0
- off = sqrt((v2-v2ref)**2 + (v3-v3ref)**2)/(R_scale*3600.0)
+ R_scale = ref_wcs.idcmodel.refpix['PSCALE'] / 3600.0
+ off = sqrt((v2 - v2ref) ** 2 + (v3 - v3ref) ** 2) / (R_scale * 3600.0)
if v3 == v3ref:
- theta=0.0
+ theta = 0.0
else:
- theta = atan2(ext_wcs.parity[0][0]*(v2-v2ref), ext_wcs.parity[1][1]*(v3-v3ref))
+ theta = atan2(ext_wcs.parity[0][0] * (v2 - v2ref),
+ ext_wcs.parity[1][1] * (v3 - v3ref))
- if ref_wcs.idcmodel.refpix['THETA']: theta += ref_wcs.idcmodel.refpix['THETA']*pi/180.0
+ if ref_wcs.idcmodel.refpix['THETA']: theta += ref_wcs.idcmodel.refpix['THETA'] * pi / 180.0
- dX=(off*sin(theta)) + offshiftx
- dY=(off*cos(theta)) + offshifty
+ dX = (off * sin(theta)) + offshiftx
+ dY = (off * cos(theta)) + offshifty
- px = np.array([[dX,dY]])
+ px = np.array([[dX, dY]])
newcrval = ref_wcs.wcs.p2s(px, 1)['world'][0]
newcrpix = np.array([ext_wcs.idcmodel.refpix['XREF'] + ltvoffx,
- ext_wcs.idcmodel.refpix['YREF'] + ltvoffy])
+ ext_wcs.idcmodel.refpix['YREF'] + ltvoffy])
ext_wcs.wcs.crval = newcrval
ext_wcs.wcs.crpix = newcrpix
ext_wcs.wcs.set()
# Create a small vector, in reference image pixel scale
- delmat = np.array([[fx[1,1], fy[1,1]], \
- [fx[1,0], fy[1,0]]]) / R_scale/3600.
+ delmat = np.array([[fx[1, 1], fy[1, 1]],
+ [fx[1, 0], fy[1, 0]]]) / R_scale / 3600.
# Account for subarray offset
# Angle of chip relative to chip
@@ -131,10 +133,10 @@ class MakeWCS(object):
wc = ref_wcs.wcs.p2s((px + dxy), 1)['world']
# Calculate the new CDs and convert to degrees
- cd11 = utils.diff_angles(wc[0,0],newcrval[0])*cos(newcrval[1]*pi/180.0)
- cd12 = utils.diff_angles(wc[1,0],newcrval[0])*cos(newcrval[1]*pi/180.0)
- cd21 = utils.diff_angles(wc[0,1],newcrval[1])
- cd22 = utils.diff_angles(wc[1,1],newcrval[1])
+ cd11 = utils.diff_angles(wc[0, 0], newcrval[0]) * cos(newcrval[1] * pi / 180.0)
+ cd12 = utils.diff_angles(wc[1, 0], newcrval[0]) * cos(newcrval[1] * pi / 180.0)
+ cd21 = utils.diff_angles(wc[0, 1], newcrval[1])
+ cd22 = utils.diff_angles(wc[1, 1], newcrval[1])
cd = np.array([[cd11, cd12], [cd21, cd22]])
ext_wcs.wcs.cd = cd
ext_wcs.wcs.set()
@@ -151,39 +153,38 @@ class MakeWCS(object):
if ref_wcs.ltv1 != 0. or ref_wcs.ltv2 != 0.:
offsetx = ref_wcs.wcs.crpix[0] - ltv1 - ref_wcs.idcmodel.refpix['XREF']
offsety = ref_wcs.wcs.crpix[1] - ltv2 - ref_wcs.idcmodel.refpix['YREF']
- ref_wcs.idcmodel.shift(offsetx,offsety)
+ ref_wcs.idcmodel.shift(offsetx, offsety)
- rfx, rfy = ref_wcs.idcmodel.cx, ref_wcs.idcmodel.cy
+ # rfx, rfy = ref_wcs.idcmodel.cx, ref_wcs.idcmodel.cy
- offshift = offsh
+ # offshift = offsh
dec = ref_wcs.wcs.crval[1]
- tddscale = (ref_wcs.pscale/ref_wcs.idcmodel.cx[1,1])
- rv23 = [ref_wcs.idcmodel.refpix['V2REF'] + (rv23_corr_tdd[0,0] *tddscale),
- ref_wcs.idcmodel.refpix['V3REF'] - (rv23_corr_tdd[1,0] * tddscale)]
+ tddscale = (ref_wcs.pscale / ref_wcs.idcmodel.cx[1, 1])
+ rv23 = [ref_wcs.idcmodel.refpix['V2REF'] + (rv23_corr_tdd[0, 0] * tddscale),
+ ref_wcs.idcmodel.refpix['V3REF'] - (rv23_corr_tdd[1, 0] * tddscale)]
# Get an approximate reference position on the sky
- rref = np.array([[ref_wcs.idcmodel.refpix['XREF']+ltvoffx ,
- ref_wcs.idcmodel.refpix['YREF']+ltvoffy]])
+ rref = np.array([[ref_wcs.idcmodel.refpix['XREF'] + ltvoffx,
+ ref_wcs.idcmodel.refpix['YREF'] + ltvoffy]])
crval = ref_wcs.wcs.p2s(rref, 1)['world'][0]
# Convert the PA_V3 orientation to the orientation at the aperture
# This is for the reference chip only - we use this for the
# reference tangent plane definition
# It has the same orientation as the reference chip
- pv = troll(ext_wcs.pav3,dec,rv23[0], rv23[1])
+ pv = troll(ext_wcs.pav3, dec, rv23[0], rv23[1])
# Add the chip rotation angle
if ref_wcs.idcmodel.refpix['THETA']:
pv += ref_wcs.idcmodel.refpix['THETA']
-
# Set values for the rest of the reference WCS
ref_wcs.wcs.crval = crval
- ref_wcs.wcs.crpix = np.array([0.0,0.0])+offsh
+ ref_wcs.wcs.crpix = np.array([0.0, 0.0]) + offsh
parity = ref_wcs.parity
- R_scale = ref_wcs.idcmodel.refpix['PSCALE']/3600.0
- cd11 = parity[0][0] * cos(pv*pi/180.0)*R_scale
- cd12 = parity[0][0] * -sin(pv*pi/180.0)*R_scale
- cd21 = parity[1][1] * sin(pv*pi/180.0)*R_scale
- cd22 = parity[1][1] * cos(pv*pi/180.0)*R_scale
+ R_scale = ref_wcs.idcmodel.refpix['PSCALE'] / 3600.0
+ cd11 = parity[0][0] * cos(pv * pi / 180.0) * R_scale
+ cd12 = parity[0][0] * -sin(pv * pi / 180.0) * R_scale
+ cd21 = parity[1][1] * sin(pv * pi / 180.0) * R_scale
+ cd22 = parity[1][1] * cos(pv * pi / 180.0) * R_scale
rcd = np.array([[cd11, cd12], [cd21, cd22]])
ref_wcs.wcs.cd = rcd
@@ -191,9 +192,9 @@ class MakeWCS(object):
uprefwcs = classmethod(uprefwcs)
- def zero_point_corr(cls,hwcs):
+ def zero_point_corr(cls, hwcs):
if hwcs.idcmodel.refpix['skew_coeffs'] is not None and 'TDD_CY_BETA' in hwcs.idcmodel.refpix['skew_coeffs']:
- v23_corr = np.array([[0.],[0.]])
+ v23_corr = np.array([[0.], [0.]])
return v23_corr
else:
try:
@@ -202,15 +203,17 @@ class MakeWCS(object):
except KeyError:
alpha = 0.0
beta = 0.0
- v23_corr = np.array([[0.],[0.]])
- logger.debug("\n\tTDD Zero point correction for chip %s defaulted to: %s" % (hwcs.chip, v23_corr))
+ v23_corr = np.array([[0.], [0.]])
+ logger.debug("TDD Zero point correction for chip"
+ "{0} defaulted to: {1}".format(hwcs.chip, v23_corr))
return v23_corr
tdd = np.array([[beta, alpha], [alpha, -beta]])
- mrotp = fileutil.buildRotMatrix(2.234529)/2048.
- xy0 = np.array([[cls.tdd_xyref[hwcs.chip][0]-2048.], [cls.tdd_xyref[hwcs.chip][1]-2048.]])
- v23_corr = np.dot(mrotp,np.dot(tdd,xy0)) * 0.05
- logger.debug("\n\tTDD Zero point correction for chip %s: %s" % (hwcs.chip, v23_corr))
+ mrotp = fileutil.buildRotMatrix(2.234529) / 2048.
+ xy0 = np.array([[cls.tdd_xyref[hwcs.chip][0] - 2048.],
+ [cls.tdd_xyref[hwcs.chip][1] - 2048.]])
+ v23_corr = np.dot(mrotp, np.dot(tdd, xy0)) * 0.05
+ logger.debug("TDD Zero point correction for chip {0}: {1}".format(hwcs.chip, v23_corr))
return v23_corr
zero_point_corr = classmethod(zero_point_corr)
@@ -258,16 +261,16 @@ def troll(roll, dec, v2, v3):
_v3 = np.deg2rad(v3 / 3600.)
# compute components
- sin_rho = sqrt((pow(sin(_v2),2)+pow(sin(_v3),2)) - (pow(sin(_v2),2)*pow(sin(_v3),2)))
+ sin_rho = sqrt((pow(sin(_v2), 2) + pow(sin(_v3), 2)) - (pow(sin(_v2), 2) * pow(sin(_v3), 2)))
rho = asin(sin_rho)
- beta = asin(sin(_v3)/sin_rho)
+ beta = asin(sin(_v3) / sin_rho)
if _v2 < 0: beta = pi - beta
- gamma = asin(sin(_v2)/sin_rho)
+ gamma = asin(sin(_v2) / sin_rho)
if _v3 < 0: gamma = pi - gamma
- A = pi/2. + _roll - beta
- B = atan2( sin(A)*cos(_dec), (sin(_dec)*sin_rho - cos(_dec)*cos(rho)*cos(A)))
+ A = pi / 2. + _roll - beta
+ B = atan2(sin(A) * cos(_dec), (sin(_dec) * sin_rho - cos(_dec) * cos(rho) * cos(A)))
# compute final value
- troll = np.rad2deg(pi - (gamma+B))
+ troll = np.rad2deg(pi - (gamma + B))
return troll
diff --git a/stwcs/updatewcs/npol.py b/stwcs/updatewcs/npol.py
index 33579f0..9c88150 100644
--- a/stwcs/updatewcs/npol.py
+++ b/stwcs/updatewcs/npol.py
@@ -1,15 +1,16 @@
-from __future__ import absolute_import, division # confidence high
+from __future__ import absolute_import, division, print_function
-import logging, time
+import logging
+import time
import numpy as np
from astropy.io import fits
from stsci.tools import fileutil
-from . import utils
logger = logging.getLogger('stwcs.updatewcs.npol')
+
class NPOLCorr(object):
"""
Defines a Lookup table prior distortion correction as per WCS paper IV.
@@ -43,7 +44,7 @@ class NPOLCorr(object):
Science file, for which a distortion correction in a NPOLFILE is available
"""
- logger.info("\n\tStarting NPOL: %s" %time.asctime())
+ logger.info("\n\tStarting NPOL: %s" % time.asctime())
try:
assert isinstance(fobj, fits.HDUList)
except AssertionError:
@@ -79,28 +80,28 @@ class NPOLCorr(object):
header = ext.header
# get the data arrays from the reference file and transform
# them for use with SIP
- dx,dy = cls.getData(nplfile, ccdchip)
+ dx, dy = cls.getData(nplfile, ccdchip)
idccoeffs = cls.getIDCCoeffs(header)
if idccoeffs is not None:
- dx, dy = cls.transformData(dx,dy, idccoeffs)
+ dx, dy = cls.transformData(dx, dy, idccoeffs)
# Determine EXTVER for the WCSDVARR extension from the
# NPL file (EXTNAME, EXTVER) kw.
# This is used to populate DPj.EXTVER kw
- wcsdvarr_x_version = 2 * extversion -1
+ wcsdvarr_x_version = 2 * extversion - 1
wcsdvarr_y_version = 2 * extversion
- for ename in zip(['DX', 'DY'], [wcsdvarr_x_version,wcsdvarr_y_version],[dx, dy]):
+ for ename in zip(['DX', 'DY'], [wcsdvarr_x_version, wcsdvarr_y_version], [dx, dy]):
error_val = ename[2].max()
cls.addSciExtKw(header, wdvarr_ver=ename[1], npol_extname=ename[0], error_val=error_val)
- hdu = cls.createNpolHDU(header, npolfile=nplfile, \
- wdvarr_ver=ename[1], npl_extname=ename[0], data=ename[2],ccdchip=ccdchip)
+ hdu = cls.createNpolHDU(header, npolfile=nplfile,
+ wdvarr_ver=ename[1], npl_extname=ename[0],
+ data=ename[2], ccdchip=ccdchip)
if wcsdvarr_ind:
fobj[wcsdvarr_ind[ename[1]]] = hdu
else:
fobj.append(hdu)
-
applyNPOLCorr = classmethod(applyNPOLCorr)
def getWCSIndex(cls, fobj):
@@ -129,17 +130,17 @@ class NPOLCorr(object):
Adds kw to sci extension to define WCSDVARR lookup table extensions
"""
- if npol_extname =='DX':
- j=1
+ if npol_extname == 'DX':
+ j = 1
else:
- j=2
-
- cperror = 'CPERR%s' %j
- cpdis = 'CPDIS%s' %j
- dpext = 'DP%s.' %j + 'EXTVER'
- dpnaxes = 'DP%s.' %j +'NAXES'
- dpaxis1 = 'DP%s.' %j+'AXIS.1'
- dpaxis2 = 'DP%s.' %j+'AXIS.2'
+ j = 2
+
+ cperror = 'CPERR%s' % j
+ cpdis = 'CPDIS%s' % j
+ dpext = 'DP%s.' % j + 'EXTVER'
+ dpnaxes = 'DP%s.' % j + 'NAXES'
+ dpaxis1 = 'DP%s.' % j + 'AXIS.1'
+ dpaxis2 = 'DP%s.' % j + 'AXIS.2'
keys = [cperror, cpdis, dpext, dpnaxes, dpaxis1, dpaxis2]
values = {cperror: error_val,
cpdis: 'Lookup',
@@ -165,15 +166,15 @@ class NPOLCorr(object):
addSciExtKw = classmethod(addSciExtKw)
- def getData(cls,nplfile, ccdchip):
+ def getData(cls, nplfile, ccdchip):
"""
Get the data arrays from the reference NPOL files
Make sure 'CCDCHIP' in the npolfile matches "CCDCHIP' in the science file.
"""
npl = fits.open(nplfile)
for ext in npl:
- nplextname = ext.header.get('EXTNAME',"")
- nplccdchip = ext.header.get('CCDCHIP',1)
+ nplextname = ext.header.get('EXTNAME', "")
+ nplccdchip = ext.header.get('CCDCHIP', 1)
if nplextname == 'DX' and nplccdchip == ccdchip:
xdata = ext.data.copy()
continue
@@ -192,7 +193,7 @@ class NPOLCorr(object):
"""
ndx, ndy = np.dot(coeffs, [dx.ravel(), dy.ravel()]).astype(np.float32)
ndx.shape = dx.shape
- ndy.shape=dy.shape
+ ndy.shape = dy.shape
return ndx, ndy
transformData = classmethod(transformData)
@@ -206,7 +207,7 @@ class NPOLCorr(object):
ocx11 = header['OCX11']
ocy10 = header['OCY10']
ocy11 = header['OCY11']
- coeffs = np.array([[ocx11, ocx10], [ocy11,ocy10]], dtype=np.float32)
+ coeffs = np.array([[ocx11, ocx10], [ocy11, ocy10]], dtype=np.float32)
except KeyError:
logger.exception('\n\tFirst order IDCTAB coefficients are not available. \n\
Cannot convert SIP to IDC coefficients.')
@@ -217,15 +218,17 @@ class NPOLCorr(object):
logger.exception("IDCSCALE not found in header - setting it to 1.")
idcscale = 1
- return np.linalg.inv(coeffs/idcscale)
+ return np.linalg.inv(coeffs / idcscale)
getIDCCoeffs = classmethod(getIDCCoeffs)
- def createNpolHDU(cls, sciheader, npolfile=None, wdvarr_ver=1, npl_extname=None,data = None, ccdchip=1):
+ def createNpolHDU(cls, sciheader, npolfile=None, wdvarr_ver=1, npl_extname=None,
+ data=None, ccdchip=1):
"""
Creates an HDU to be added to the file object.
"""
- hdr = cls.createNpolHdr(sciheader, npolfile=npolfile, wdvarr_ver=wdvarr_ver, npl_extname=npl_extname, ccdchip=ccdchip)
+ hdr = cls.createNpolHdr(sciheader, npolfile=npolfile, wdvarr_ver=wdvarr_ver,
+ npl_extname=npl_extname, ccdchip=ccdchip)
hdu = fits.ImageHDU(header=hdr, data=data)
return hdu
@@ -256,29 +259,29 @@ class NPOLCorr(object):
npl.close()
naxis = npl[1].header['NAXIS']
- ccdchip = nplextname #npol_header['CCDCHIP']
+ ccdchip = nplextname # npol_header['CCDCHIP']
- kw = { 'NAXIS': 'Size of the axis',
- 'CDELT': 'Coordinate increment along axis',
- 'CRPIX': 'Coordinate system reference pixel',
- 'CRVAL': 'Coordinate system value at reference pixel',
- }
+ kw = {'NAXIS': 'Size of the axis',
+ 'CDELT': 'Coordinate increment along axis',
+ 'CRPIX': 'Coordinate system reference pixel',
+ 'CRVAL': 'Coordinate system value at reference pixel',
+ }
kw_comm1 = {}
kw_val1 = {}
for key in kw.keys():
- for i in range(1, naxis+1):
+ for i in range(1, naxis + 1):
si = str(i)
- kw_comm1[key+si] = kw[key]
+ kw_comm1[key + si] = kw[key]
- for i in range(1, naxis+1):
+ for i in range(1, naxis + 1):
si = str(i)
- kw_val1['NAXIS'+si] = npol_header.get('NAXIS'+si)
- kw_val1['CDELT'+si] = npol_header.get('CDELT'+si, 1.0) * \
- sciheader.get('LTM'+si+'_'+si, 1)
- kw_val1['CRPIX'+si] = npol_header.get('CRPIX'+si, 0.0)
- kw_val1['CRVAL'+si] = (npol_header.get('CRVAL'+si, 0.0) - \
- sciheader.get('LTV'+str(i), 0))
+ kw_val1['NAXIS' + si] = npol_header.get('NAXIS' + si)
+ kw_val1['CDELT' + si] = npol_header.get('CDELT' + si, 1.0) * \
+ sciheader.get('LTM' + si + '_' + si, 1)
+ kw_val1['CRPIX' + si] = npol_header.get('CRPIX' + si, 0.0)
+ kw_val1['CRVAL' + si] = (npol_header.get('CRVAL' + si, 0.0) -
+ sciheader.get('LTV' + str(i), 0))
kw_comm0 = {'XTENSION': 'Image extension',
'BITPIX': 'IEEE floating point',
@@ -289,15 +292,15 @@ class NPOLCorr(object):
'GCOUNT': 'One data group',
}
- kw_val0 = { 'XTENSION': 'IMAGE',
- 'BITPIX': -32,
- 'NAXIS': naxis,
- 'EXTNAME': 'WCSDVARR',
- 'EXTVER': wdvarr_ver,
- 'PCOUNT': 0,
- 'GCOUNT': 1,
- 'CCDCHIP': ccdchip,
- }
+ kw_val0 = {'XTENSION': 'IMAGE',
+ 'BITPIX': -32,
+ 'NAXIS': naxis,
+ 'EXTNAME': 'WCSDVARR',
+ 'EXTVER': wdvarr_ver,
+ 'PCOUNT': 0,
+ 'GCOUNT': 1,
+ 'CCDCHIP': ccdchip,
+ }
cdl = []
for key in kw_comm0.keys():
cdl.append((key, kw_val0[key], kw_comm0[key]))
@@ -310,11 +313,11 @@ class NPOLCorr(object):
for i, c in enumerate(npol_phdr):
if c == 'FILENAME':
start_indx = i
- if c == '': # remove blanks from end of header
- end_indx = i+1
+ if c == '': # remove blanks from end of header
+ end_indx = i + 1
break
if start_indx >= 0:
- for card in npol_phdr.cards[start_indx:end_indx]:
+ for card in npol_phdr.cards[start_indx: end_indx]:
cdl.append(card)
hdr = fits.Header(cards=cdl)
@@ -331,7 +334,7 @@ class NPOLCorr(object):
if fobj[0].header['INSTRUME'] == 'ACS' and fobj[0].header['DETECTOR'] == 'WFC':
ccdchip = fobj[extname, extver].header['CCDCHIP']
elif fobj[0].header['INSTRUME'] == 'WFC3' and fobj[0].header['DETECTOR'] == 'UVIS':
- ccdchip = fobj[extname, extver].header['CCDCHIP']
+ ccdchip = fobj[extname, extver].header['CCDCHIP']
elif fobj[0].header['INSTRUME'] == 'WFPC2':
ccdchip = fobj[extname, extver].header['DETECTOR']
elif fobj[0].header['INSTRUME'] == 'STIS':
diff --git a/stwcs/updatewcs/utils.py b/stwcs/updatewcs/utils.py
index ebfc701..8c8d416 100644
--- a/stwcs/updatewcs/utils.py
+++ b/stwcs/updatewcs/utils.py
@@ -1,4 +1,4 @@
-from __future__ import division # confidence high
+from __future__ import absolute_import, division, print_function
import os
from astropy.io import fits
from stsci.tools import fileutil
@@ -7,7 +7,7 @@ import logging
logger = logging.getLogger("stwcs.updatewcs.utils")
-def diff_angles(a,b):
+def diff_angles(a, b):
"""
Perform angle subtraction a-b taking into account
small-angle differences across 360degree line.
@@ -23,6 +23,7 @@ def diff_angles(a,b):
return diff
+
def getBinning(fobj, extver=1):
# Return the binning factor
binned = 1
@@ -30,9 +31,10 @@ def getBinning(fobj, extver=1):
mode = fobj[0].header.get('MODE', "")
if mode == 'AREA': binned = 2
else:
- binned = fobj['SCI', extver].header.get('BINAXIS',1)
+ binned = fobj['SCI', extver].header.get('BINAXIS', 1)
return binned
+
def updateNEXTENDKw(fobj):
"""
Updates PRIMARY header with correct value for NEXTEND, if present.
@@ -44,9 +46,10 @@ def updateNEXTENDKw(fobj):
"""
if 'nextend' in fobj[0].header:
- fobj[0].header['nextend'] = len(fobj) -1
+ fobj[0].header['nextend'] = len(fobj) - 1
+
-def extract_rootname(kwvalue,suffix=""):
+def extract_rootname(kwvalue, suffix=""):
""" Returns the rootname from a full reference filename
If a non-valid value (any of ['','N/A','NONE','INDEF',None]) is input,
@@ -55,13 +58,13 @@ def extract_rootname(kwvalue,suffix=""):
This function will also replace any 'suffix' specified with a blank.
"""
# check to see whether a valid kwvalue has been provided as input
- if kwvalue.strip() in ['','N/A','NONE','INDEF',None]:
- return 'NONE' # no valid value, so return 'NONE'
+ if kwvalue.strip() in ['', 'N/A', 'NONE', 'INDEF', None]:
+ return 'NONE' # no valid value, so return 'NONE'
# for a valid kwvalue, parse out the rootname
# strip off any environment variable from input filename, if any are given
if '$' in kwvalue:
- fullval = kwvalue[kwvalue.find('$')+1:]
+ fullval = kwvalue[kwvalue.find('$') + 1: ]
else:
fullval = kwvalue
# Extract filename without path from kwvalue
@@ -71,10 +74,11 @@ def extract_rootname(kwvalue,suffix=""):
rootname = fileutil.buildNewRootname(fname)
# Now, remove any known suffix from rootname
- rootname = rootname.replace(suffix,'')
+ rootname = rootname.replace(suffix, '')
return rootname.strip()
-def construct_distname(fobj,wcsobj):
+
+def construct_distname(fobj, wcsobj):
"""
This function constructs the value for the keyword 'DISTNAME'.
It relies on the reference files specified by the keywords 'IDCTAB',
@@ -85,8 +89,8 @@ def construct_distname(fobj,wcsobj):
and have a value of 'NONE' if no reference files are specified.
"""
idcname = extract_rootname(fobj[0].header.get('IDCTAB', "NONE"),
- suffix='_idc')
- if (idcname is None or idcname=='NONE') and wcsobj.sip is not None:
+ suffix='_idc')
+ if (idcname is None or idcname == 'NONE') and wcsobj.sip is not None:
idcname = 'UNKNOWN'
npolname, npolfile = build_npolname(fobj)
@@ -99,29 +103,30 @@ def construct_distname(fobj,wcsobj):
sipname, idctab = build_sipname(fobj)
- distname = build_distname(sipname,npolname,d2imname)
- return {'DISTNAME':distname,'SIPNAME':sipname}
+ distname = build_distname(sipname, npolname, d2imname)
+ return {'DISTNAME': distname, 'SIPNAME': sipname}
+
-def build_distname(sipname,npolname,d2imname):
+def build_distname(sipname, npolname, d2imname):
"""
Core function to build DISTNAME keyword value without the HSTWCS input.
"""
-
distname = sipname.strip()
if npolname != 'NONE' or d2imname != 'NONE':
if d2imname != 'NONE':
- distname+= '-'+npolname.strip() + '-'+d2imname.strip()
+ distname += '-' + npolname.strip() + '-' + d2imname.strip()
else:
- distname+='-'+npolname.strip()
+ distname += '-' + npolname.strip()
return distname
-def build_default_wcsname(idctab):
- idcname = extract_rootname(idctab,suffix='_idc')
+def build_default_wcsname(idctab):
+ idcname = extract_rootname(idctab, suffix='_idc')
wcsname = 'IDC_' + idcname
return wcsname
+
def build_sipname(fobj, fname=None, sipname=None):
"""
Build a SIPNAME from IDCTAB
@@ -141,10 +146,10 @@ def build_sipname(fobj, fname=None, sipname=None):
"""
try:
idctab = fobj[0].header['IDCTAB']
- idcname = extract_rootname(idctab,suffix='_idc')
+ idcname = extract_rootname(idctab, suffix='_idc')
except KeyError:
idctab = 'N/A'
- idcname= 'N/A'
+ idcname = 'N/A'
if not fname:
try:
fname = fobj.filename()
@@ -161,7 +166,7 @@ def build_sipname(fobj, fname=None, sipname=None):
rootname = fobj[0].header['rootname']
except KeyError:
rootname = fname
- sipname = rootname +'_'+ idcname
+ sipname = rootname + '_' + idcname
else:
if 'A_ORDER' in fobj[1].header or 'B_ORDER' in fobj[1].header:
sipname = 'UNKNOWN'
@@ -170,6 +175,7 @@ def build_sipname(fobj, fname=None, sipname=None):
return sipname, idctab
+
def build_npolname(fobj, npolfile=None):
"""
Build a NPOLNAME from NPOLFILE
@@ -203,6 +209,7 @@ def build_npolname(fobj, npolfile=None):
npolname = 'NOMODEL'
return npolname, npolfile
+
def build_d2imname(fobj, d2imfile=None):
"""
Build a D2IMNAME from D2IMFILE
@@ -227,11 +234,11 @@ def build_d2imname(fobj, d2imfile=None):
d2imname = 'UNKNOWN'
else:
d2imname = 'NOMODEL'
- d2imname = extract_rootname(d2imfile,suffix='_d2i')
+ d2imname = extract_rootname(d2imfile, suffix='_d2i')
if d2imname == 'NONE':
d2imname = 'NOMODEL'
else:
- d2imname = extract_rootname(d2imfile,suffix='_d2i')
+ d2imname = extract_rootname(d2imfile, suffix='_d2i')
if d2imname == 'NONE':
d2imname = 'NOMODEL'
return d2imname, d2imfile
diff --git a/stwcs/updatewcs/wfpc2_dgeo.py b/stwcs/updatewcs/wfpc2_dgeo.py
index e57bb5c..b76198a 100644
--- a/stwcs/updatewcs/wfpc2_dgeo.py
+++ b/stwcs/updatewcs/wfpc2_dgeo.py
@@ -1,6 +1,7 @@
""" wfpc2_dgeo - Functions to convert WFPC2 DGEOFILE into D2IMFILE
"""
+from __future__ import absolute_import, division, print_function
import os
import datetime
@@ -13,6 +14,7 @@ from stsci.tools import fileutil
import logging
logger = logging.getLogger("stwcs.updatewcs.apply_corrections")
+
def update_wfpc2_d2geofile(filename, fhdu=None):
"""
Creates a D2IMFILE from the DGEOFILE for a WFPC2 image (input), and
@@ -35,7 +37,7 @@ def update_wfpc2_d2geofile(filename, fhdu=None):
image header will be updated/added to point to this newly created file.
"""
-
+
close_fhdu = False
if fhdu is None:
fhdu = fileutil.openImage(filename, mode='update')
@@ -48,7 +50,7 @@ def update_wfpc2_d2geofile(filename, fhdu=None):
dgeofile = fhdu['PRIMARY'].header.get('ODGEOFIL', None)
logger.info('Converting DGEOFILE %s into D2IMFILE...' % dgeofile)
rootname = filename[:filename.find('.fits')]
- d2imfile = convert_dgeo_to_d2im(dgeofile,rootname)
+ d2imfile = convert_dgeo_to_d2im(dgeofile, rootname)
fhdu['PRIMARY'].header['ODGEOFIL'] = dgeofile
fhdu['PRIMARY'].header['DGEOFILE'] = 'N/A'
fhdu['PRIMARY'].header['D2IMFILE'] = d2imfile
@@ -67,25 +69,26 @@ def update_wfpc2_d2geofile(filename, fhdu=None):
# (multidrizzle clean=True mode of operation)
return d2imfile
-def convert_dgeo_to_d2im(dgeofile,output,clobber=True):
+
+def convert_dgeo_to_d2im(dgeofile, output, clobber=True):
""" Routine that converts the WFPC2 DGEOFILE into a D2IMFILE.
"""
dgeo = fileutil.openImage(dgeofile)
- outname = output+'_d2im.fits'
+ outname = output + '_d2im.fits'
removeFileSafely(outname)
- data = np.array([dgeo['dy',1].data[:,0]])
+ data = np.array([dgeo['dy', 1].data[:, 0]])
scihdu = fits.ImageHDU(data=data)
dgeo.close()
# add required keywords for D2IM header
scihdu.header['EXTNAME'] = ('DY', 'Extension name')
scihdu.header['EXTVER'] = (1, 'Extension version')
- fits_str = 'PYFITS Version '+str(astropy.__version__)
+ fits_str = 'PYFITS Version ' + str(astropy.__version__)
scihdu.header['ORIGIN'] = (fits_str, 'FITS file originator')
scihdu.header['INHERIT'] = (False, 'Inherits global header')
dnow = datetime.datetime.now()
- scihdu.header['DATE'] = (str(dnow).replace(' ','T'),
+ scihdu.header['DATE'] = (str(dnow).replace(' ', 'T'),
'Date FITS file was generated')
scihdu.header['CRPIX1'] = (0, 'Distortion array reference pixel')
@@ -116,7 +119,7 @@ def convert_dgeo_to_d2im(dgeofile,output,clobber=True):
return outname
-def removeFileSafely(filename,clobber=True):
+def removeFileSafely(filename, clobber=True):
""" Delete the file specified, but only if it exists and clobber is True.
"""
if filename is not None and filename.strip() != '':
generated by cgit (git 2.34.1) at 2025-05-02 11:57:29 -0400