def update_refchip_with_shift(chip_wcs, wcslin, fitgeom='rscale',
                              rot=0.0, scale=1.0, xsh=0.0, ysh=0.0,
                              fit=None, xrms=None, yrms=None):
    """ Compute the matrix for the scale and rotation correction

    Parameters
    ----------
    chip_wcs: wcs object
        HST of the input image
    wcslin: wcs object
        Reference WCS from which the offsets/rotations are determined
    fitgeom: str
        NOT USED
    rot : float
        Amount of rotation measured in fit to be applied.
        [Default=0.0]
    scale : float
        Amount of scale change measured in fit to be applied.
        [Default=1.0]
    xsh : float
        Offset in X pixels from defined tangent plane to be applied to image.
        [Default=0.0]
    ysh : float
        Offset in Y pixels from defined tangent plane to be applied to image.
        [Default=0.0]
    fit : arr
        Linear coefficients for fit
        [Default = None]
    xrms : float
        RMS of fit in RA (in decimal degrees) that will be recorded as
        CRDER1 in WCS and header
        [Default = None]
    yrms : float
        RMS of fit in Dec (in decimal degrees) that will be recorded as
        CRDER2 in WCS and header
        [Default = None]
        """
    # compute the matrix for the scale and rotation correction
    if fit is None:
        fit = linearfit.buildFitMatrix(rot, scale)

    shift = np.asarray([xsh, ysh]) - np.dot(wcslin.wcs.crpix, fit) + wcslin.wcs.crpix

    fit = np.linalg.inv(fit).T

    cwcs = chip_wcs.deepcopy()
    cd_eye = np.eye(chip_wcs.wcs.cd.shape[0], dtype=ndfloat128)
    zero_shift = np.zeros(2, dtype=ndfloat128)

    naxis1, naxis2 = chip_wcs.pixel_shape

    # estimate precision necessary for iterative processes:
    maxiter = 100
    crpix2corners = np.dstack([i.flatten() for i in np.meshgrid(
        [1, naxis1], [1, naxis2])])[0] - chip_wcs.wcs.crpix
    maxUerr = 1.0e-5 / np.amax(np.linalg.norm(crpix2corners, axis=1))

    # 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(20.0, (chip_wcs.wcs.crpix[0] - 1.0) / 100.0,
                      (naxis1 - chip_wcs.wcs.crpix[0]) / 100.0))
    hy = max(1.0, min(20.0, (chip_wcs.wcs.crpix[1] - 1.0) / 100.0,
                      (naxis2 - chip_wcs.wcs.crpix[1]) / 100.0))

    # compute new CRVAL for the image WCS:
    crpixinref = wcslin.wcs_world2pix(
        chip_wcs.wcs_pix2world([chip_wcs.wcs.crpix],1),1)
    crpixinref = np.dot(fit, (crpixinref - shift).T).T
    chip_wcs.wcs.crval = wcslin.wcs_pix2world(crpixinref, 1)[0]
    chip_wcs.wcs.set()

    # initial approximation for CD matrix of the image WCS:
    (U, u) = linearize(cwcs, chip_wcs, wcslin, chip_wcs.wcs.crpix,
                       fit, shift, hx=hx, hy=hy)
    err0 = np.amax(np.abs(U-cd_eye)).astype(np.float64)
    chip_wcs.wcs.cd = np.dot(chip_wcs.wcs.cd.astype(ndfloat128), U).astype(np.float64)
    chip_wcs.wcs.set()

    # NOTE: initial solution is the exact mathematical solution (modulo numeric
    # differentiation). However, e.g., due to rounding errors, approximate
    # numerical differentiation, the solution may be improved by performing
    # several iterations. The next step will try to perform
    # fixed-point iterations to "improve" the solution
    # but this is not really required.

    # Perform fixed-point iterations to improve the approximation
    # for CD matrix of the image WCS (actually for the U matrix).
    for i in range(maxiter):
        (U, u) = linearize(chip_wcs, chip_wcs, wcslin, chip_wcs.wcs.crpix,
                           cd_eye, zero_shift, hx=hx, hy=hy)
        err = np.amax(np.abs(U-cd_eye)).astype(np.float64)
        if err > err0:
            break
        chip_wcs.wcs.cd = np.dot(chip_wcs.wcs.cd, U).astype(np.float64)
        chip_wcs.wcs.set()
        if err < maxUerr:
            break
        err0 = err

    if xrms is not None:
        chip_wcs.wcs.crder = np.array([xrms,yrms])