How to use tweakwcs - 10 common examples

To help you get started, we’ve selected a few tweakwcs examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github spacetelescope / tweakwcs / tweakwcs / tpwcs.py View on Github external
**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)
            )
github spacetelescope / tweakwcs / tweakwcs / linearfit.py View on Github external
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)
github spacetelescope / tweakwcs / tweakwcs / tpwcs.py View on Github external
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]
github spacetelescope / tweakwcs / tweakwcs / linearfit.py View on Github external
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
github spacetelescope / tweakwcs / tweakwcs / linearfit.py View on Github external
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
github spacetelescope / tweakwcs / tweakwcs / wcsimage.py View on Github external
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']))
github spacetelescope / tweakwcs / tweakwcs / tpwcs.py View on Github external
# 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'
github spacetelescope / tweakwcs / tweakwcs / tpwcs.py View on Github external
"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()
github spacetelescope / tweakwcs / tweakwcs / wcsimage.py View on Github external
# 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
github spacetelescope / tweakwcs / tweakwcs / tpwcs.py View on Github external
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)