diff options
Diffstat (limited to 'stwcs/distortion/models.py')
| -rw-r--r-- | stwcs/distortion/models.py | 185 |
1 files changed, 71 insertions, 114 deletions
diff --git a/stwcs/distortion/models.py b/stwcs/distortion/models.py index 231a9f1..136c6cf 100644 --- a/stwcs/distortion/models.py +++ b/stwcs/distortion/models.py @@ -1,14 +1,9 @@ -from __future__ import absolute_import, division, print_function # confidence high - +from __future__ import absolute_import, division, print_function import numpy as np - -# Import PyDrizzle utility modules from . import mutil from .mutil import combin -yes = True -no = False ################# # @@ -72,8 +67,8 @@ class GeometryModel: to the reference position of the chip. """ - _cxs = np.zeros(shape=self.cx.shape,dtype=self.cx.dtype) - _cys = np.zeros(shape=self.cy.shape,dtype=self.cy.dtype) + _cxs = np.zeros(shape=self.cx.shape, dtype=self.cx.dtype) + _cys = np.zeros(shape=self.cy.shape, dtype=self.cy.dtype) _k = self.norder + 1 # loop over each input coefficient for m in range(_k): @@ -83,15 +78,17 @@ class GeometryModel: _ilist = list(range(m, _k)) # sum from m to k for i in _ilist: - _jlist = list(range(n, i - (m-n)+1)) + _jlist = list(range(n, i - (m - n) + 1)) # sum from n to i-(m-n) for j in _jlist: - _cxs[m,n] += self.cx[i,j]*combin(j,n)*combin((i-j),(m-n))*pow(xs,(j-n))*pow(ys,((i-j)-(m-n))) - _cys[m,n] += self.cy[i,j]*combin(j,n)*combin((i-j),(m-n))*pow(xs,(j-n))*pow(ys,((i-j)-(m-n))) + _cxs[m, n] += self.cx[i, j] * combin(j, n) * combin((i - j), (m - n)) * \ + pow(xs, (j - n)) * pow(ys, ((i - j) - (m - n))) + _cys[m, n] += self.cy[i, j] * combin(j, n) * combin((i - j), (m - n)) * \ + pow(xs, (j - n)) * pow(ys, ((i - j) - (m - n))) self.cx = _cxs.copy() self.cy = _cys.copy() - def convert(self, tmpname, xref=None,yref=None,delta=yes): + def convert(self, tmpname, xref=None, yref=None, delta=True): """ Open up an ASCII file, output coefficients in drizzle format after converting them as necessary. @@ -104,89 +101,62 @@ class GeometryModel: Added ability to shift coefficients to new center for support of subarrays. """ - cx = self.cx/self.pscale - cy = self.cy/self.pscale - x0 = self.refpix['XDELTA'] + cx[0,0] - y0 = self.refpix['YDELTA'] + cy[0,0] - #xr = self.refpix['XREF'] - #yr = self.refpix['YREF'] + cx = self.cx / self.pscale + cy = self.cy / self.pscale + x0 = self.refpix['XDELTA'] + cx[0, 0] + y0 = self.refpix['YDELTA'] + cy[0, 0] xr = self.refpix['CHIP_XREF'] yr = self.refpix['CHIP_YREF'] - - - ''' - if xref != None: - # Shift coefficients for use with drizzle - _xs = xref - self.refpix['XREF'] + 1.0 - _ys = yref - self.refpix['YREF'] + 1.0 - - - if _xs != 0 or _ys != 0: - cxs,cys= self.shift(cx, cy, _xs, _ys) - cx = cxs - cy = cys - - # We only want to apply this shift to coeffs - # for subarray images. - if delta == no: - cxs[0,0] = cxs[0,0] - _xs - cys[0,0] = cys[0,0] - _ys - - # Now, apply only the difference introduced by the distortion.. - # i.e., (undistorted - original) shift. - x0 += cxs[0,0] - y0 += cys[0,0] - ''' - self.x0 = x0 #+ 1.0 - self.y0 = y0 #+ 1.0 + self.x0 = x0 + self.y0 = y0 # Now, write out the coefficients into an ASCII # file in 'drizzle' format. lines = [] - - lines.append('# Polynomial distortion coefficients\n') - lines.append('# Extracted from "%s" \n'%self.name) - lines.append('refpix %f %f \n'%(xr,yr)) - if self.norder==3: + lines.append('# Extracted from "%s" \n' % self.name) + lines.append('refpix %f %f \n' % (xr, yr)) + if self.norder == 3: lines.append('cubic\n') - elif self.norder==4: + elif self.norder == 4: lines.append('quartic\n') - elif self.norder==5: + elif self.norder == 5: lines.append('quintic\n') else: raise ValueError("Drizzle cannot handle poly distortions of order %d" % self.norder) - str = '%16.8f %16.8g %16.8g %16.8g %16.8g \n'% (x0,cx[1,1],cx[1,0],cx[2,2],cx[2,1]) + str = '%16.8f %16.8g %16.8g %16.8g %16.8g \n' % (x0, cx[1, 1], cx[1, 0], cx[2, 2], cx[2, 1]) lines.append(str) - str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cx[2,0],cx[3,3],cx[3,2],cx[3,1],cx[3,0]) + str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cx[2, 0], cx[3, 3], cx[3, 2], cx[3, 1], cx[3, 0]) lines.append(str) - if self.norder>3: - str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cx[4,4],cx[4,3],cx[4,2],cx[4,1],cx[4,0]) + if self.norder > 3: + str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cx[4, 4], cx[4, 3], cx[4, 2], cx[4, 1], cx[4, 0]) lines.append(str) - if self.norder>4: - str = '%16.8g %16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cx[5,5],cx[5,4],cx[5,3],cx[5,2],cx[5,1],cx[5,0]) + if self.norder > 4: + str = '%16.8g %16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cx[5, 5], cx[5, 4], cx[5, 3], cx[5, 2], cx[5, 1], cx[5, 0]) lines.append(str) lines.append("\n") - str = '%16.8f %16.8g %16.8g %16.8g %16.8g \n'% (y0,cy[1,1],cy[1,0],cy[2,2],cy[2,1]) + str = '%16.8f %16.8g %16.8g %16.8g %16.8g \n' % (y0, cy[1, 1], cy[1, 0], cy[2, 2], cy[2, 1]) lines.append(str) - str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cy[2,0],cy[3,3],cy[3,2],cy[3,1],cy[3,0]) + str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cy[2, 0], cy[3, 3], cy[3, 2], cy[3, 1], + cy[3, 0]) lines.append(str) - if self.norder>3: - str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cy[4,4],cy[4,3],cy[4,2],cy[4,1],cy[4,0]) + if self.norder > 3: + str = '%16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cy[4, 4], cy[4, 3], cy[4, 2], + cy[4, 1], cy[4, 0]) lines.append(str) - if self.norder>4: - str = '%16.8g %16.8g %16.8g %16.8g %16.8g %16.8g \n'% (cy[5,5],cy[5,4],cy[5,3],cy[5,2],cy[5,1],cy[5,0]) + if self.norder > 4: + str = '%16.8g %16.8g %16.8g %16.8g %16.8g %16.8g \n' % (cy[5, 5], cy[5, 4], cy[5, 3], + cy[5, 2], cy[5, 1], cy[5, 0]) lines.append(str) - output = open(tmpname,'w') + output = open(tmpname, 'w') output.writelines(lines) output.close() - - def apply(self, pixpos,scale=1.0,order=None): + def apply(self, pixpos, scale=1.0, order=None): """ Apply coefficients to a pixel position or a list of positions. This should be the same for all coefficients tables. @@ -195,13 +165,13 @@ class GeometryModel: Compute delta from reference position """ - + """ scale actually is a ratio of pscale/self.model.pscale what is pscale? """ - if self.cx == None: - return pixpos[:,0],pixpos[:,1] + if self.cx is None: + return pixpos[:, 0], pixpos[:, 1] if order is None: order = self.norder @@ -209,29 +179,29 @@ class GeometryModel: # Apply in the same way that 'drizzle' would... _cx = self.cx / (self.pscale * scale) _cy = self.cy / (self.pscale * scale) - _convert = no + _convert = False _p = pixpos # Do NOT include any zero-point terms in CX,CY here # as they should not be scaled by plate-scale like rest # of coeffs... This makes the computations consistent # with 'drizzle'. WJH 17-Feb-2004 - _cx[0,0] = 0. - _cy[0,0] = 0. + _cx[0, 0] = 0. + _cy[0, 0] = 0. if isinstance(_p, list) or isinstance(_p, tuple): - _p = np.array(_p,dtype=np.float64) - _convert = yes + _p = np.array(_p, dtype=np.float64) + _convert = True - dxy = _p - (self.refpix['XREF'],self.refpix['YREF']) + dxy = _p - (self.refpix['XREF'], self.refpix['YREF']) # Apply coefficients from distortion model here... c = _p * 0. - for i in range(order+1): - for j in range(i+1): - c[:,0] = c[:,0] + _cx[i][j] * pow(dxy[:,0],j) * pow(dxy[:,1],(i-j)) - c[:,1] = c[:,1] + _cy[i][j] * pow(dxy[:,0],j) * pow(dxy[:,1],(i-j)) - xc = c[:,0] - yc = c[:,1] + for i in range(order + 1): + for j in range(i + 1): + c[:, 0] = c[:, 0] + _cx[i][j] * pow(dxy[:, 0], j) * pow(dxy[:, 1], (i - j)) + c[:, 1] = c[:, 1] + _cy[i][j] * pow(dxy[:, 0], j) * pow(dxy[:, 1], (i - j)) + xc = c[:, 0] + yc = c[:, 1] # Convert results back to same form as original input if _convert: @@ -242,11 +212,11 @@ class GeometryModel: xc = xc[0] yc = yc[0] - return xc,yc + return xc, yc - def setPScaleCoeffs(self,pscale): - self.cx[1,1] = pscale - self.cy[1,0] = pscale + def setPScaleCoeffs(self, pscale): + self.cx[1, 1] = pscale + self.cy[1, 0] = pscale self.refpix['PSCALE'] = pscale self.pscale = pscale @@ -259,7 +229,7 @@ class IDCModel(GeometryModel): We also need to read in SCALE, XCOM,YCOM, XREF,YREF as well. """ def __init__(self, idcfile, date=None, chip=1, direction='forward', - filter1='CLEAR1',filter2='CLEAR2',offtab=None, binned=1): + filter1='CLEAR1', filter2='CLEAR2', offtab=None, binned=1): GeometryModel.__init__(self) # # Norder must be derived from the coeffs file itself, @@ -269,9 +239,9 @@ class IDCModel(GeometryModel): # Populate cx,cy,scale, and other variables here. # self.name = idcfile - self.cx,self.cy,self.refpix,self.norder = mutil.readIDCtab(idcfile, - chip=chip,direction=direction,filter1=filter1,filter2=filter2, - date=date, offtab=offtab) + self.cx, self.cy, self.refpix, self.norder = mutil.readIDCtab(idcfile, + chip=chip, direction=direction, filter1=filter1, filter2=filter2, + date=date, offtab=offtab) if 'empty_model' in self.refpix and self.refpix['empty_model']: pass @@ -285,61 +255,54 @@ class IDCModel(GeometryModel): self.refpix['YSIZE'] = self.refpix['YSIZE'] / binned self.pscale = self.refpix['PSCALE'] - + class WCSModel(GeometryModel): """ This class sets up a distortion model based on coefficients found in the image header. """ - def __init__(self,header,rootname): + def __init__(self, header, rootname): GeometryModel.__init__(self) - if 'rootname' in header: self.name = header['rootname'] else: self.name = rootname # Initialize all necessary distortion arrays with # default model... - #self.cx,self.cy,self.refpix,self.order = mutil.defaultModel() + # self.cx,self.cy,self.refpix,self.order = mutil.defaultModel() # Read in values from header, and update distortion arrays. - self.cx,self.cy,self.refpix,self.norder = mutil.readWCSCoeffs(header) + self.cx, self.cy, self.refpix, self.norder = mutil.readWCSCoeffs(header) self.pscale = self.refpix['PSCALE'] - class DrizzleModel(GeometryModel): """ This class will read in an ASCII Cubic drizzle coeffs file and populate the cx,cy arrays. """ - def __init__(self, idcfile, scale = None): + def __init__(self, idcfile, scale=None): GeometryModel.__init__(self) # # We now need to read in the file, populate cx,cy, and # other variables as necessary. # self.name = idcfile - self.cx,self.cy,self.refpix,self.norder = mutil.readCubicTable(idcfile) + self.cx, self.cy, self.refpix, self.norder = mutil.readCubicTable(idcfile) - # scale is the ratio wcs.pscale/model.pscale. + # scale is the ratio wcs.pscale/model.pscale. # model.pscale for WFPC2 is passed from REFDATA. # This is needed for WFPC2 binned data. - - if scale != None: + + if scale is not None: self.pscale = scale else: self.pscale = self.refpix['PSCALE'] - """ - The above definition looks wrong. - In one case it's a ratio in the other it's pscale. - - """ class TraugerModel(GeometryModel): """ @@ -349,14 +312,8 @@ class TraugerModel(GeometryModel): """ NORDER = 3 - def __init__(self, idcfile,lam): + def __init__(self, idcfile, lam): GeometryModel.__init__(self) self.name = idcfile - self.cx,self.cy,self.refpix,self.norder = mutil.readTraugerTable(idcfile,lam) + self.cx, self.cy, self.refpix, self.norder = mutil.readTraugerTable(idcfile, lam) self.pscale = self.refpix['PSCALE'] - # - # Read in file here. - # Populate cx,cy, and other variables. - # - - |