HSTWCS API
- class stwcs.wcsutil.hstwcs.HSTWCS(fobj=None, ext=None, minerr=0.0, wcskey=' ')
- all_world2pix(*args, **kwargs)
all_world2pix(*arg, tolerance=1.0e-4, maxiter=20, adaptive=False, detect_divergence=True, quiet=False)
Performs full inverse transformation using iterative solution on full forward transformation with complete distortion model.
- Parameters:
tolerance (float, optional (Default = 1.0e-4)) – Absolute tolerance required the solution. Iteration terminates when the correction to the solution found during the previous iteration is smaller (in the sence of the L2 norm) than
accuracy
.maxiter (int, optional (Default = 20)) – Maximum number of iterations allowed to reach the solution.
adaptive (bool, optional (Default = False)) –
Specifies whether to adaptively select only points that did not converge to a solution whithin the required accuracy for the next iteration. Default is recommended for HST as well as most other instruments.
Note
The
all_world2pix()
uses a vectorized implementation of the method of consecutive approximations (seeNotes
section below) in which it iterates over all input poits regardless until the required accuracy has been reached for all input points. In some cases it may be possible that almost all points have reached the required accuracy but there are only a few of input data points left for which additional iterations may be needed (this depends mostly on the characteristics of the geometric distortions for a given instrument). In this situation it may be advantageous to setadaptive
=True
in which caseall_world2pix()
will continue iterating only over the points that have not yet converged to the required accuracy. However, for the HST’s ACS/WFC detector, which has the strongest distortions of all HST instruments, testing has shown that enabling this option would lead to a about 10-30 % penalty in computational time (depending on specifics of the image, geometric distortions, and number of input points to be converted). Therefore, for HST instruments, it is recommended to setadaptive
=False
. The only danger in getting this setting wrong will be a performance penalty.Note
When
detect_divergence
isTrue
,all_world2pix()
will automatically switch to the adaptive algorithm once divergence has been detected.detect_divergence (bool, optional (Default = True)) –
Specifies whether to perform a more detailed analysis of the convergence to a solution. Normally
all_world2pix()
may not achieve the required accuracy if either thetolerance
ormaxiter
arguments are too low. However, it may happen that for some geometric distortions the conditions of convergence for the the method of consecutive approximations used byall_world2pix()
may not be satisfied, in which case consecutive approximations to the solution will diverge regardless of thetolerance
ormaxiter
settings.When
detect_divergence
isFalse
, these divergent points will be detected as not having achieved the required accuracy (without further details). In addition, ifadaptive
isFalse
then the algorithm will not know that the solution (for specific points) is diverging and will continue iterating and trying to “improve” diverging solutions. This may result in NaN or Inf values in the return results (in addition to a performance penalties). Even whendetect_divergence
isFalse
,all_world2pix()
, at the end of the iterative process, will identify invalid results (NaN or Inf) as “diverging” solutions and will raiseNoConvergence
unless thequiet
parameter is set toTrue
.When
detect_divergence
isTrue
,all_world2pix()
will detect points for which current correction to the coordinates is larger than the correction applied during the previous iteration if the requested accuracy has not yet been achieved. In this case, ifadaptive
isTrue
, these points will be excluded from further iterations and ifadaptive
isFalse
,all_world2pix()
will automatically switch to the adaptive algorithm.Note
When absolute tolerance has been achieved, small increases in current corrections may be possible due to rounding errors (when
adaptive
isFalse
) and such increases will be ignored.Note
Setting
detect_divergence
toTrue
will incurr about 5-10% performance penalty (in our testing on ACS/WFC images). Because the benefits of enabling this feature outweigh the small performance penalty, it is recommended to setdetect_divergence
toTrue
, unless extensive testing of the distortion models for images from specific instruments show a good stability of the numerical method for a wide range of coordinates (even outside the image itself).Note
Indices of the diverging inverse solutions will be reported in the
divergent
attribute of the raisedNoConvergence
object.quiet (bool, optional (Default = False)) – Do not throw
NoConvergence
exceptions when the method does not converge to a solution with the required accuracy within a specified number of maximum iterations set bymaxiter
parameter. Instead, simply return the found solution.
- Raises:
NoConvergence – The method does not converge to a solution with the required accuracy within a specified number of maximum iterations set by the
maxiter
parameter.
Notes
Inputs can either be (RA, Dec, origin) or (RADec, origin) where RA and Dec are 1-D arrays/lists of coordinates and RADec is an array/list of pairs of coordinates.
Using the method of consecutive approximations we iterate starting with the initial approximation, which is computed using the non-distorion-aware
wcs_world2pix()
(or equivalent).The
all_world2pix()
function uses a vectorized implementation of the method of consecutive approximations and therefore it is highly efficient (>30x) when all data points that need to be converted from sky coordinates to image coordinates are passed at once. Therefore, it is advisable, whenever possible, to pass as input a long array of all points that need to be converted toall_world2pix()
instead of callingall_world2pix()
for each data point. Also see the note to theadaptive
parameter.Examples
>>> import stwcs >>> from astropy.io import fits >>> hdulist = fits.open('j94f05bgq_flt.fits') >>> w = stwcs.wcsutil.HSTWCS(hdulist, ext=('sci',1)) >>> hdulist.close()
>>> ra, dec = w.all_pix2world([1,2,3],[1,1,1],1); print(ra); print(dec) [ 5.52645241 5.52649277 5.52653313] [-72.05171776 -72.05171295 -72.05170814] >>> radec = w.all_pix2world([[1,1],[2,1],[3,1]],1); print(radec) [[ 5.52645241 -72.05171776] [ 5.52649277 -72.05171295] [ 5.52653313 -72.05170814]] >>> x, y = w.all_world2pix(ra,dec,1) >>> print(x) [ 1.00000233 2.00000232 3.00000233] >>> print(y) [ 0.99999997 0.99999997 0.99999998] >>> xy = w.all_world2pix(radec,1) >>> print(xy) [[ 1.00000233 0.99999997] [ 2.00000232 0.99999997] [ 3.00000233 0.99999998]] >>> xy = w.all_world2pix(radec,1, maxiter=3, tolerance=1.0e-10, quiet=False) NoConvergence: 'HSTWCS.all_world2pix' failed to converge to requested accuracy after 3 iterations.
>>> Now try to use some diverging data: >>> divradec = w.all_pix2world([[1.0,1.0],[10000.0,50000.0], [3.0,1.0]],1); print(divradec) [[ 5.52645241 -72.05171776] [ 7.15979392 -70.81405561] [ 5.52653313 -72.05170814]]
>>> try: >>> xy = w.all_world2pix(divradec,1, maxiter=20, tolerance=1.0e-4, adaptive=False, detect_divergence=True, quiet=False) >>> except stwcs.wcsutil.hstwcs.NoConvergence as e: >>> print("Indices of diverging points: {}".format(e.divergent)) >>> print("Indices of poorly converging points: {}".format(e.failed2converge)) >>> print("Best solution: {}".format(e.best_solution)) >>> print("Achieved accuracy: {}".format(e.accuracy)) >>> raise e Indices of diverging points: [1] Indices of poorly converging points: None Best solution: [[ 1.00006219e+00 9.99999288e-01] [ -1.99440907e+06 1.44308548e+06] [ 3.00006257e+00 9.99999316e-01]] Achieved accuracy: [[ 5.98554253e-05 6.79918148e-07] [ 8.59514088e+11 6.61703754e+11] [ 6.02334592e-05 6.59713067e-07]] Traceback (innermost last): File "<console>", line 8, in <module> NoConvergence: 'HSTWCS.all_world2pix' failed to converge to the requested accuracy. After 5 iterations, the solution is diverging at least for one input point.
>>> try: >>> xy = w.all_world2pix(divradec,1, maxiter=20, tolerance=1.0e-4, adaptive=False, detect_divergence=False, quiet=False) >>> except stwcs.wcsutil.hstwcs.NoConvergence as e: >>> print("Indices of diverging points: {}".format(e.divergent)) >>> print("Indices of poorly converging points: {}".format(e.failed2converge)) >>> print("Best solution: {}".format(e.best_solution)) >>> print("Achieved accuracy: {}".format(e.accuracy)) >>> raise e Indices of diverging points: [1] Indices of poorly converging points: None Best solution: [[ 1. 1.] [ nan nan] [ 3. 1.]] Achieved accuracy: [[ 0. 0.] [ nan nan] [ 0. 0.]] Traceback (innermost last): File "<console>", line 8, in <module> NoConvergence: 'HSTWCS.all_world2pix' failed to converge to the requested accuracy. After 20 iterations, the solution is diverging at least for one input point.
- property naxis1
- property naxis2
- pc2cd()
- printwcs()
Print the basic WCS keywords.
- readIDCCoeffs(header)
Reads in first order IDCTAB coefficients if present in the header
- readModel(update=False, header=None)
Reads distortion model from IDCTAB.
If IDCTAB is not found (‘N/A’, “”, or not found on disk), then if SIP coefficients and first order IDCTAB coefficients are present in the header, restore the idcmodel from the header. If not - assign None to self.idcmodel.
- Parameters:
header (
astropy.io.fits.Header
) – fits extension headerupdate (bool (False)) – if True - record the following IDCTAB quantities as header keywords: CX10, CX11, CY10, CY11, IDCSCALE, IDCTHETA, IDCXREF, IDCYREF, IDCV2REF, IDCV3REF
- readModelFromIDCTAB(header=None, update=False)
Read distortion model from idc table.
- Parameters:
header (
astropy.io.fits.Header
) – fits extension headerupdate (booln (False)) – if True - save teh following as header keywords: CX10, CX11, CY10, CY11, IDCSCALE, IDCTHETA, IDCXREF, IDCYREF, IDCV2REF, IDCV3REF
- resetLTV()
Reset LTV values for polarizer data
The polarizer field is smaller than the detector field. The distortion coefficients are defined for the entire polarizer field and the LTV values are set as with subarray data. This may also be true for other special filters. This is a special case when the observation is considered a subarray in terms of detector field but a full frame in terms of distortion model. To avoid shifting the distortion coefficients the LTV values are reset to 0.
- setInstrSpecKw(prim_hdr=None, ext_hdr=None)
Populate the instrument specific attributes:
These can be in different headers but each instrument class has knowledge of where to look for them.
- Parameters:
prim_hdr (
astropy.io.fits.Header
) – primary headerext_hdr (
astropy.io.fits.Header
) – extension header
- setOrient()
Computes ORIENTAT from the CD matrix
- setPscale()
Calculates the plate scale from the CD matrix
- updatePscale(scale)
Updates the CD matrix with a new plate scale
- wcs2header(sip2hdr=False, idc2hdr=True, wcskey=None, relax=False)
Create a
astropy.io.fits.Header
object from WCS keywords.If the original header had a CD matrix, return a CD matrix, otherwise return a PC matrix.
- Parameters:
sip2hdr (bool) – If True - include SIP coefficients