# Path: pygeoid/integrals/core.py
# class SphericalKernel(Kernel):
#     """Base class for all spherical kernels.
# 
#     """
# 
#     def _check_spherical_distance(self, spherical_distance):
#         """Check spherical distance.
# 
#         """
#         # if not np.logical_and(spherical_distance >= 0 * u.deg,
#         #                      spherical_distance <= np.pi * u.rad).any():
#         if not Angle(spherical_distance).is_within_bounds(0 * u.deg,
#                                                           180 * u.deg):
#             raise ValueError('spherical_distance must be between 0 and 180 degrees.')
# 
#         return spherical_distance
# 
# Path: pygeoid/integrals/veningmeinesz.py
# class VeningMeineszKernel(SphericalKernel):
#     r"""Vening Meinesz kernel class.
# 
#     """
# 
#     _name = 'Vening Meinesz'
# 
#     @u.quantity_input
#     def kernel(self, spherical_distance: u.deg):
#         r"""Evaluate Vening Meinesz kernel.
# 
#         Parameters
#         ----------
#         spherical_distance : float or array_like of floats
#             Spherical distance, in radians.
#         degrees : bool, optional
#             If True, the spherical distance is given in degrees,
#             otherwise radians.
# 
#         Notes
#         -----
#         The derivative  is the Vening-Meinesz and it  depends
#         on the spherical distance :math:`\psi` by [1]_:
# 
#         .. math ::
#             V\left(\psi\right) =
#             \dfrac{d S\left(\psi\right)}{d\psi} = - \dfrac{\cos{(\psi /
#             2)}}{2\sin^2{(\psi / 2)}} + 8\sin{\psi} -
#             6\cos{(\psi / 2)} - 3\dfrac{1 - \sin{(\psi / 2)}}{\sin{\psi}} +
#             3\sin{\psi}\ln{\left[\sin{(\psi/2)} + \sin^2{(\psi/2)}\right]},
# 
#         which is the derivative of Stokes's kernel with respect to
#         :math:`\psi`.
# 
#         References
#         ----------
#         .. [1] Heiskanen WA, Moritz H (1967) Physical geodesy.
#             Freeman, San Francisco
# 
#         """
# 
#         psi = self._check_spherical_distance(
#             spherical_distance=spherical_distance)
# 
#         cpsi2 = np.cos(psi / 2)
#         spsi2 = np.sin(psi / 2)
#         spsi = np.sin(psi)
# 
#         return -0.5 * cpsi2 / spsi2**2 + 8 * spsi - 6 * cpsi2 -\
#             3 * (1 - spsi2) / spsi + 3 * spsi * np.log(spsi2 + spsi2**2)


import numpy as np
import astropy.units as u
from pygeoid.integrals.core import SphericalKernel
from pygeoid.integrals.veningmeinesz import VeningMeineszKernel
        The derivative of the Stokes's kernel by
        :math:`t = \cos{\psi}`, :math:`-1 \leq t \leq 1` is [1]_:

        .. math ::
            \dfrac{d S\left(t\right)}{dt} = - 8 +
            \dfrac{3\sqrt{2}}{\sqrt{1-t}} +
            \dfrac{1}{\sqrt{2}\left(1-t\right)^{3/2}} +
            \dfrac{3\left(\sqrt{2} - \sqrt{1-t}\right)}
            {\sqrt{2}\left(1-t^2\right)} -
            3\ln{\dfrac{\sqrt{1-t}\left(\sqrt{2} + \sqrt{1-t}\right)}{2}}.

        References
        ----------
        .. [1] Hagiwara Y (1976) A new formula for evaluating the truncation
            error coefficient. Bulletin Géodésique 50:131–135.

        """
        if not np.logical_and(t >= -1, t <= 1).any():
            raise ValueError('t must be between -1 and 1.')

        sq2 = np.sqrt(2)
        sq1t = np.sqrt(1 - t)

        kernel = -8 + 3 * sq2 / sq1t + 1 / (sq2 * sq1t**3)
        kernel += 3 * (sq2 - sq1t) / (sq2 * (1 - t**2))
        kernel -= 3 * np.log(0.5 * (sq1t * (sq2 + sq1t)))

        return kernel



# Path: pygeoid/integrals/core.py
# class SphericalKernel(Kernel):
#     """Base class for all spherical kernels.
# 
#     """
# 
#     def _check_spherical_distance(self, spherical_distance):
#         """Check spherical distance.
# 
#         """
#         # if not np.logical_and(spherical_distance >= 0 * u.deg,
#         #                      spherical_distance <= np.pi * u.rad).any():
#         if not Angle(spherical_distance).is_within_bounds(0 * u.deg,
#                                                           180 * u.deg):
#             raise ValueError('spherical_distance must be between 0 and 180 degrees.')
# 
#         return spherical_distance
# 
# Path: pygeoid/integrals/veningmeinesz.py
# class VeningMeineszKernel(SphericalKernel):
#     r"""Vening Meinesz kernel class.
# 
#     """
# 
#     _name = 'Vening Meinesz'
# 
#     @u.quantity_input
#     def kernel(self, spherical_distance: u.deg):
#         r"""Evaluate Vening Meinesz kernel.
# 
#         Parameters
#         ----------
#         spherical_distance : float or array_like of floats
#             Spherical distance, in radians.
#         degrees : bool, optional
#             If True, the spherical distance is given in degrees,
#             otherwise radians.
# 
#         Notes
#         -----
#         The derivative  is the Vening-Meinesz and it  depends
#         on the spherical distance :math:`\psi` by [1]_:
# 
#         .. math ::
#             V\left(\psi\right) =
#             \dfrac{d S\left(\psi\right)}{d\psi} = - \dfrac{\cos{(\psi /
#             2)}}{2\sin^2{(\psi / 2)}} + 8\sin{\psi} -
#             6\cos{(\psi / 2)} - 3\dfrac{1 - \sin{(\psi / 2)}}{\sin{\psi}} +
#             3\sin{\psi}\ln{\left[\sin{(\psi/2)} + \sin^2{(\psi/2)}\right]},
# 
#         which is the derivative of Stokes's kernel with respect to
#         :math:`\psi`.
# 
#         References
#         ----------
#         .. [1] Heiskanen WA, Moritz H (1967) Physical geodesy.
#             Freeman, San Francisco
# 
#         """
# 
#         psi = self._check_spherical_distance(
#             spherical_distance=spherical_distance)
# 
#         cpsi2 = np.cos(psi / 2)
#         spsi2 = np.sin(psi / 2)
#         spsi = np.sin(psi)
# 
#         return -0.5 * cpsi2 / spsi2**2 + 8 * spsi - 6 * cpsi2 -\
#             3 * (1 - spsi2) / spsi + 3 * spsi * np.log(spsi2 + spsi2**2)


# Path: pygeoid/integrals/stokes.py
import numpy as np
import astropy.units as u
from pygeoid.integrals.core import SphericalKernel
from pygeoid.integrals.veningmeinesz import VeningMeineszKernel

        The derivative of the Stokes's kernel by
        :math:`t = \cos{\psi}`, :math:`-1 \leq t \leq 1` is [1]_:

        .. math ::
            \dfrac{d S\left(t\right)}{dt} = - 8 +
            \dfrac{3\sqrt{2}}{\sqrt{1-t}} +
            \dfrac{1}{\sqrt{2}\left(1-t\right)^{3/2}} +
            \dfrac{3\left(\sqrt{2} - \sqrt{1-t}\right)}
            {\sqrt{2}\left(1-t^2\right)} -
            3\ln{\dfrac{\sqrt{1-t}\left(\sqrt{2} + \sqrt{1-t}\right)}{2}}.

        References
        ----------
        .. [1] Hagiwara Y (1976) A new formula for evaluating the truncation
            error coefficient. Bulletin Géodésique 50:131–135.

        """
        if not np.logical_and(t >= -1, t <= 1).any():
            raise ValueError('t must be between -1 and 1.')

        sq2 = np.sqrt(2)
        sq1t = np.sqrt(1 - t)

        kernel = -8 + 3 * sq2 / sq1t + 1 / (sq2 * sq1t**3)
        kernel += 3 * (sq2 - sq1t) / (sq2 * (1 - t**2))
        kernel -= 3 * np.log(0.5 * (sq1t * (sq2 + sq1t)))

        return kernel



class StokesExtendedKernel(SphericalKernel):