Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
**kwargs: optional parameters
Optional parameters for the WCS corrector. `JWSTgWCS` ignores these
arguments (except for storing them in the ``meta`` attribute).
"""
frms = self._wcs.available_frames
if ref_tpwcs is None:
matrix = np.array(matrix, dtype=np.double)
shift = np.array(shift, dtype=np.double)
else:
# compute linear transformation from the tangent plane used for
# alignment to the tangent plane of this wcs:
r, t = _tp2tp(ref_tpwcs, self)
matrix = np.linalg.multi_dot([r, matrix, inv(r)]).astype(np.double)
shift = (np.dot(r, shift) - np.dot(matrix, t) + t).astype(np.double)
# if original WCS did not have tangent-plane corrections, create
# new correction and add it to the WCs pipeline:
if self._tpcorr is None:
self._tpcorr = JWSTgWCS._tpcorr_init(
v2_ref=self._wcsinfo['v2_ref'] / 3600.0,
v3_ref=self._wcsinfo['v3_ref'] / 3600.0,
roll_ref=self._wcsinfo['roll_ref']
)
JWSTgWCS._tpcorr_combine_affines(
self._tpcorr,
matrix,
_ARCSEC2RAD * np.asarray(shift)
)
sxu = np.dot(w, x * u)
syu = np.dot(w, y * u)
sxv = np.dot(w, x * v)
syv = np.dot(w, y * v)
suu = np.dot(w, u * u)
svv = np.dot(w, v * v)
suv = np.dot(w, u * v)
m = np.array([[su, sv, sw], [suu, suv, su], [suv, svv, sv]],
dtype=np.longdouble)
a = np.array([sx, sxu, sxv], dtype=np.longdouble)
b = np.array([sy, syu, syv], dtype=np.longdouble)
try:
inv_m = inv(m)
except np.linalg.LinAlgError:
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
)
p = np.dot(inv_m, a)
q = np.dot(inv_m, b)
if not (np.all(np.isfinite(p)) and np.all(np.isfinite(q))):
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
) # pragma: no cover
# Return the shift, rotation, and scale changes
fit = _build_fit(p, q, 'general')
resids = xy - np.dot(uv, fit['matrix_ld'].T) - fit['shift_ld']
fit['resids'] = resids.astype(np.double)
frms = self._wcs.available_frames
if ref_tpwcs is None:
matrix = np.array(matrix, dtype=np.double)
shift = np.array(shift, dtype=np.double)
else:
# compute linear transformation from the tangent plane used for
# alignment to the tangent plane of this wcs:
r, t = _tp2tp(ref_tpwcs, self)
matrix = np.linalg.multi_dot([r, matrix, inv(r)]).astype(np.double)
shift = (np.dot(r, shift) - np.dot(matrix, t) + t).astype(np.double)
# if original WCS did not have tangent-plane corrections, create
# new correction and add it to the WCs pipeline:
if self._tpcorr is None:
self._tpcorr = JWSTgWCS._tpcorr_init(
v2_ref=self._wcsinfo['v2_ref'] / 3600.0,
v3_ref=self._wcsinfo['v3_ref'] / 3600.0,
roll_ref=self._wcsinfo['roll_ref']
)
JWSTgWCS._tpcorr_combine_affines(
self._tpcorr,
matrix,
_ARCSEC2RAD * np.asarray(shift)
)
self._partial_tpcorr = JWSTgWCS._v2v3_to_tpcorr_from_full(self._tpcorr)
idx_v2v3 = frms.index(self._v23name)
pipeline = deepcopy(self._wcs.pipeline)
pf, pt = pipeline[idx_v2v3]
m = np.array([[su, sv, sw], [suu, suv, su], [suv, svv, sv]],
dtype=np.longdouble)
a = np.array([sx, sxu, sxv], dtype=np.longdouble)
b = np.array([sy, syu, syv], dtype=np.longdouble)
try:
inv_m = inv(m)
except np.linalg.LinAlgError:
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
)
p = np.dot(inv_m, a)
q = np.dot(inv_m, b)
if not (np.all(np.isfinite(p)) and np.all(np.isfinite(q))):
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
) # pragma: no cover
# Return the shift, rotation, and scale changes
fit = _build_fit(p, q, 'general')
resids = xy - np.dot(uv, fit['matrix_ld'].T) - fit['shift_ld']
fit['resids'] = resids.astype(np.double)
_compute_stat(fit, residuals=resids, weights=w)
return fit
syu = np.dot(w, y * u)
sxv = np.dot(w, x * v)
syv = np.dot(w, y * v)
suu = np.dot(w, u * u)
svv = np.dot(w, v * v)
suv = np.dot(w, u * v)
m = np.array([[su, sv, sw], [suu, suv, su], [suv, svv, sv]],
dtype=np.longdouble)
a = np.array([sx, sxu, sxv], dtype=np.longdouble)
b = np.array([sy, syu, syv], dtype=np.longdouble)
try:
inv_m = inv(m)
except np.linalg.LinAlgError:
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
)
p = np.dot(inv_m, a)
q = np.dot(inv_m, b)
if not (np.all(np.isfinite(p)) and np.all(np.isfinite(q))):
raise SingularMatrixError(
"Singular matrix: suspected colinear points."
) # pragma: no cover
# Return the shift, rotation, and scale changes
fit = _build_fit(p, q, 'general')
resids = xy - np.dot(uv, fit['matrix_ld'].T) - fit['shift_ld']
fit['resids'] = resids.astype(np.double)
_compute_stat(fit, residuals=resids, weights=w)
return fit
ref_idx = self._mref_idx
refxy = refxy[ref_idx]
# process weights:
if 'weight' in self._catalog.colnames:
im_weight = np.asarray(self._catalog['weight'])[minput_idx]
else:
im_weight = None
if 'weight' in refcat.catalog.colnames:
ref_weight = np.asarray(refcat.catalog['weight'])[ref_idx]
else:
ref_weight = None
fit = iter_linear_fit(
refxy, im_xyref, ref_weight, im_weight,
fitgeom=fitgeom, nclip=nclip, sigma=sigma, center=None
)
# re-compute shifts for the center at (0, 0):
fit['shift_ld'] += fit['center_ld'] - np.dot(fit['center_ld'], fit['matrix_ld'].T)
fit['shift'] = fit['shift_ld'].astype(np.double)
xy_fit = fit['shift'] + np.dot(im_xyref[fit['fitmask']], fit['matrix'].T)
fit['fit_xy'] = xy_fit
fit['fit_RA'], fit['fit_DEC'] = tanplane_wcs.tanp_to_world(*(xy_fit.T))
log.info("Computed '{:s}' fit for {}:".format(fitgeom, self.name))
if fitgeom == 'shift':
log.info("XSH: {:.6g} YSH: {:.6g}".format(*fit['shift']))
# compute linear transformation from the tangent plane used for
# alignment to the tangent plane of this wcs:
r, t = _tp2tp(ref_tpwcs, self)
matrix = np.linalg.multi_dot([r, matrix, inv(r)]).astype(np.double)
shift = (np.dot(r, shift) - np.dot(matrix, t) + t).astype(np.double)
# if original WCS did not have tangent-plane corrections, create
# new correction and add it to the WCs pipeline:
if self._tpcorr is None:
self._tpcorr = JWSTgWCS._tpcorr_init(
v2_ref=self._wcsinfo['v2_ref'] / 3600.0,
v3_ref=self._wcsinfo['v3_ref'] / 3600.0,
roll_ref=self._wcsinfo['roll_ref']
)
JWSTgWCS._tpcorr_combine_affines(
self._tpcorr,
matrix,
_ARCSEC2RAD * np.asarray(shift)
)
self._partial_tpcorr = JWSTgWCS._v2v3_to_tpcorr_from_full(self._tpcorr)
idx_v2v3 = frms.index(self._v23name)
pipeline = deepcopy(self._wcs.pipeline)
pf, pt = pipeline[idx_v2v3]
pipeline[idx_v2v3] = (pf, deepcopy(self._tpcorr))
frm_v2v3corr = deepcopy(pf)
frm_v2v3corr.name = 'v2v3corr'
pipeline.insert(idx_v2v3 + 1, (frm_v2v3corr, pt))
self._wcs = gwcs.WCS(pipeline, name=self._owcs.name)
self._v23name = 'v2v3corr'
"WCS/TPCorr parameters 'v2ref', 'v3ref', and/or 'roll' "
"differ from the corresponding reference values."
)
self._partial_tpcorr = JWSTgWCS._v2v3_to_tpcorr_from_full(
self._tpcorr
)
else:
self._v23name = 'v2v3'
self._tpcorr = None
self._default_tpcorr = JWSTgWCS._tpcorr_init(
v2_ref=v2_ref / 3600.0,
v3_ref=v3_ref / 3600.0,
roll_ref=roll_ref
)
self._partial_tpcorr = JWSTgWCS._v2v3_to_tpcorr_from_full(
self._default_tpcorr
)
self._update_transformations()
# Find an "optimal" tangent plane to the catalog points based on the
# mean point and then construct a WCS based on the mean point.
# Compute x, y coordinates in this tangent plane based on the
# previously computed WCS and return the set of x, y coordinates and
# "reference WCS".
x, y, z = _S2C(self.catalog['RA'], self.catalog['DEC'])
ra_ref, dec_ref = _C2S(
x.mean(dtype=np.double),
y.mean(dtype=np.double),
z.mean(dtype=np.double)
)
rotm = [planar_rot_3d(np.deg2rad(alpha), 2 - axis)
for axis, alpha in enumerate([ra_ref, dec_ref])]
euler_rot = np.linalg.multi_dot(rotm)
inv_euler_rot = inv(euler_rot)
xr, yr, zr = np.dot(euler_rot, (x, y, z))
x = yr / xr
y = zr / xr
xv, yv = convex_hull(x, y)
if len(xv) == 0:
# no points
raise RuntimeError( # pragma: no cover
"Unexpected error: Contact software developer"
)
elif len(xv) == 1:
# one point. build a small box around it:
x, y = convex_hull(x, y, wcs=None)
tol = 0.5 * self._footprint_tol
meta: dict, None, optional
Dictionary that will be merged to the object's ``meta`` fields.
**kwargs: optional parameters
Optional parameters for the WCS corrector. `FITSWCS` ignores these
arguments (except for storing them in the ``meta`` attribute).
"""
wcs = self._wcs
orig_wcs = wcs.deepcopy()
if ref_tpwcs is None:
ref_tpwcs = FITSWCS(wcs.deepcopy())
naxis1, naxis2 = wcs.pixel_shape
shift = -np.dot(inv(matrix), shift)
# estimate step for numerical differentiation. We need a step
# large enough to avoid rounding errors and small enough to get a
# better precision for numerical differentiation.
# TODO: The logic below should be revised at a later time so that it
# better takes into account the two competing requirements.
hx = max(1.0, min(10, (wcs.wcs.crpix[0] - 1.0) / 100.0,
(naxis1 - wcs.wcs.crpix[0]) / 100.0))
hy = max(1.0, min(10, (wcs.wcs.crpix[1] - 1.0) / 100.0,
(naxis2 - wcs.wcs.crpix[1]) / 100.0))
# compute new CRVAL for the image WCS:
crpix2d = np.atleast_2d(wcs.wcs.crpix)
crval = wcs.wcs_pix2world(crpix2d, 1).ravel()
crpixinref = ref_tpwcs.world_to_tanp(*crval)