def METHOD_NAME(self, x, y, mu_r, mu_t, parity, phi_G, center_x=0, center_y=0):
    """
    :param x: x-coord (in angles)
    :param y: y-coord (in angles)
    :param mu_r: radial magnification
    :param mu_t: tangential magnification
    :param parity: parity of the side of the macromodel. Either +1 (positive parity) or -1 (negative parity)
    :param phi_G: shear orientation angle (relative to the x-axis)
    :return: hessian matrix (in angles)
    """
    # positive parity case
    if parity == 1:
        gamma = (1.0 / mu_t - 1.0 / mu_r) * 0.5
        kappa = 1 - gamma - 1.0 / mu_r
    # negative parity case
    elif parity == -1:
        gamma = (1.0 / mu_t + 1.0 / mu_r) * 0.5
        kappa = 1 - gamma + 1.0 / mu_r
    else:
        raise ValueError(
            "%f is not a valid value for the parity of the macromodel. Choose either +1 or -1."
            % parity
        )
    # compute the shear along the x and y directions, rotate the vector in the opposite direction than the reference frame (compare with util.rotate)
    gamma1, gamma2 = gamma * np.cos(2 * phi_G), -gamma * np.sin(2 * phi_G)
    f_xx = kappa + gamma1
    f_yy = kappa - gamma1
    f_xy = -gamma2
    return f_xx, f_xy, f_xy, f_yy