<commit_before>#!usr/bin/env python
# Filename: math.py
# pylint: disable=C0103
"""
Maths, Geometry, coordinates.
"""

from matplotlib.path import Path
import numpy as np
from scipy import linalg, stats
from sklearn.neighbors import KernelDensity


from .logger import logging

__author__ = "<NAME> and <NAME>"
__copyright__ = "Copyright 2017, Tamas Gal and the KM3NeT collaboration."
__credits__ = []
__license__ = "MIT"
__maintainer__ = "<NAME> and <NAME>"
__email__ = "<EMAIL>"
__status__ = "Development"

log = logging.getLogger(__name__)  # pylint: disable=C0103


def neutrino_to_source_direction(phi, theta, radian=True):
    """Flip the direction.

    Parameters
    ==========
    phi, theta: neutrino direction
    radian: bool [default=True]
        receive + return angles in radian? (if false, use degree)
    """
    phi = np.atleast_1d(phi).copy()
    theta = np.atleast_1d(theta).copy()
    if not radian:
        phi *= np.pi / 180
        theta *= np.pi / 180
    assert np.all(phi <= 2 * np.pi)
    assert np.all(theta <= np.pi)
    azimuth = (phi + np.pi) % (2 * np.pi)
    zenith = np.pi - theta
    if not radian:
        azimuth *= 180 / np.pi
        zenith *= 180 / np.pi
    return azimuth, zenith


def source_to_neutrino_direction(azimuth, zenith, radian=True):
    """Flip the direction.

    Parameters
    ==========
    zenith, azimuth: neutrino origin
    radian: bool [default=True]
        receive + return angles in radian? (if false, use degree)
    """
    azimuth = np.atleast_1d(azimuth).copy()
    zenith = np.atleast_1d(zenith).copy()
    if not radian:
        azimuth *= np.pi / 180
        zenith *= np.pi / 180
    phi = (azimuth - np.pi) % (2 * np.pi)
    theta = np.pi - zenith
    if not radian:
        phi *= 180 / np.pi
        theta *= 180 / np.pi
    return phi, theta


def theta(v):
    """Neutrino direction in polar coordinates.

    Downgoing event: theta = 180deg
    Horizont: 90deg
    Upgoing: theta = 0

    Angles in radians.
    """
    v = np.atleast_2d(v)
    dir_z = v[:, 2]
    return np.arccos(dir_z)


def phi(v):
    """Neutrino direction in polar coordinates.

    ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``.

    Angles in radians.
    """
    v = np.atleast_2d(v)
    dir_x = v[:, 0]
    dir_y = v[:, 1]
    p = np.arctan2(dir_y, dir_x)
    p[p < 0] += 2 * np.pi
    # if len(p) == 1:
    #     return p[0]
    return p


def zenith(v):
    """Return the zenith angle in radians.

    Defined as 'Angle respective to downgoing'.
    Downgoing event: zenith = 0
    Horizont: 90deg
    Upgoing: zenith = 180deg
    """
    return angle_between((0, 0, -1), v)


def azimuth(v):
    """Return the azimuth angle in radians.

    ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``.

    This is the 'normal' azimuth definition -- beware of how you
    define your coordinates. KM3NeT defines azimuth
    differently than e.g. SLALIB, astropy, the AAS.org
    """
    v = np.atleast_2d(v)
    azi = phi(v) - np.pi
    azi[azi < 0] += 2 * np.pi
    if len(azi) == 1:
        return azi[0]
    return azi


def cartesian(phi, theta, radius=1):
    x = radius * np.sin(theta) * np.cos(phi)
    y = radius * np.sin(theta) * np.sin(phi)
    z = radius * np.cos(theta)
    return np.column_stack((x, y, z))


def angle_between(v1, v2):
    """Returns the angle in radians between vectors 'v1' and 'v2'::

    >>> angle_between((1, 0, 0), (0, 1, 0))
    1.5707963267948966
    >>> angle_between((1, 0, 0), (1, 0, 0))
    0.0
    >>> angle_between((1, 0, 0), (-1, 0, 0))
    3.141592653589793

    """
    v1_u = unit_vector(v1)
    v2_u = unit_vector(v2)
    # Don't use `np.dot`, does not work with all shapes
    angle = np.arccos(np.inner(v1_u, v2_u))
    return angle


def innerprod_1d(v1, v2):
    """1d Inner product for vector-of-vectors.

    Example:
    ========
    ```
    v1 = np.array([dir_x, dir_y, dir_z]).T
    v2 = ... # dito
    angle_between_v1_v1 = innerprod_1d(v1, v2)
    ```
    """
    return np.einsum('ij,ij->i', v1, v1)


def unit_vector(vector, **kwargs):
    """Returns the unit vector of the vector."""
    # This also works for a dataframe with columns ['x', 'y', 'z']
    # However, the division operation is picky about the shapes
    # So, remember input vector shape, cast all up to 2d,
    # do the (ugly) conversion, then return unit in same shape as input
    # of course, the numpy-ized version of the input...
    vector = np.array(vector)
    out_shape = vector.shape
    vector = np.atleast_2d(vector)
    unit = vector / linalg.norm(vector, axis=1, **kwargs)[:, None]
    return unit.reshape(out_shape)


def pld3(p1, p2, d2):
    """Calculate the point-line-distance for given point and line."""
    c = np.cross(d2, p2 - p1)
    n1 = scipy.linalg.norm(c)
    n2 = scipy.linalg.norm(d2)
    return n1 / n2


def lpnorm(x, p=2):
    return np.power(np.sum(np.power(x, p)), 1 / p)


def dist(x1, x2):
    np.linalg.norm(np.array(x2) - np.array(x1), axis=1)


def com(points, masses=None):
    """Calculate center of mass for given points.
    If masses is not set, assume equal masses."""
    if masses is None:
        return np.average(points, axis=0)
    else:
        return np.average(points, axis=0, weights=masses)<commit_msg>I changed my mind. sorry for the spam<commit_after>#!usr/bin/env python
# Filename: math.py
# pylint: disable=C0103
"""
Maths, Geometry, coordinates.
"""

from matplotlib.path import Path
import numpy as np
from scipy import linalg, stats
from sklearn.neighbors import KernelDensity


from .logger import logging

__author__ = "<NAME> and <NAME>"
__copyright__ = "Copyright 2017, Tamas Gal and the KM3NeT collaboration."
__credits__ = []
__license__ = "MIT"
__maintainer__ = "<NAME> and <NAME>"
__email__ = "<EMAIL>"
__status__ = "Development"

log = logging.getLogger(__name__)  # pylint: disable=C0103


def neutrino_to_source_direction(phi, theta, radian=True):
    """Flip the direction.

    Parameters
    ==========
    phi, theta: neutrino direction
    radian: bool [default=True]
        receive + return angles in radian? (if false, use degree)
    """
    phi = np.atleast_1d(phi).copy()
    theta = np.atleast_1d(theta).copy()
    if not radian:
        phi *= np.pi / 180
        theta *= np.pi / 180
    assert np.all(phi <= 2 * np.pi)
    assert np.all(theta <= np.pi)
    azimuth = (phi + np.pi) % (2 * np.pi)
    zenith = np.pi - theta
    if not radian:
        azimuth *= 180 / np.pi
        zenith *= 180 / np.pi
    return azimuth, zenith


def source_to_neutrino_direction(azimuth, zenith, radian=True):
    """Flip the direction.

    Parameters
    ==========
    zenith, azimuth: neutrino origin
    radian: bool [default=True]
        receive + return angles in radian? (if false, use degree)
    """
    azimuth = np.atleast_1d(azimuth).copy()
    zenith = np.atleast_1d(zenith).copy()
    if not radian:
        azimuth *= np.pi / 180
        zenith *= np.pi / 180
    phi = (azimuth - np.pi) % (2 * np.pi)
    theta = np.pi - zenith
    if not radian:
        phi *= 180 / np.pi
        theta *= 180 / np.pi
    return phi, theta


def theta(v):
    """Neutrino direction in polar coordinates.

    Downgoing event: theta = 180deg
    Horizont: 90deg
    Upgoing: theta = 0

    Angles in radians.
    """
    v = np.atleast_2d(v)
    dir_z = v[:, 2]
    return np.arccos(dir_z)


def phi(v):
    """Neutrino direction in polar coordinates.

    ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``.

    Angles in radians.
    """
    v = np.atleast_2d(v)
    dir_x = v[:, 0]
    dir_y = v[:, 1]
    p = np.arctan2(dir_y, dir_x)
    p[p < 0] += 2 * np.pi
    # if len(p) == 1:
    #     return p[0]
    return p


def zenith(v):
    """Return the zenith angle in radians.

    Defined as 'Angle respective to downgoing'.
    Downgoing event: zenith = 0
    Horizont: 90deg
    Upgoing: zenith = 180deg
    """
    return angle_between((0, 0, -1), v)


def azimuth(v):
    """Return the azimuth angle in radians.

    ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``.

    This is the 'normal' azimuth definition -- beware of how you
    define your coordinates. KM3NeT defines azimuth
    differently than e.g. SLALIB, astropy, the AAS.org
    """
    v = np.atleast_2d(v)
    azi = phi(v) - np.pi
    azi[azi < 0] += 2 * np.pi
    if len(azi) == 1:
        return azi[0]
    return azi


def cartesian(phi, theta, radius=1):
    x = radius * np.sin(theta) * np.cos(phi)
    y = radius * np.sin(theta) * np.sin(phi)
    z = radius * np.cos(theta)
    return np.column_stack((x, y, z))


def angle_between(v1, v2):
    """Returns the angle in radians between vectors 'v1' and 'v2'::

    >>> angle_between((1, 0, 0), (0, 1, 0))
    1.5707963267948966
    >>> angle_between((1, 0, 0), (1, 0, 0))
    0.0
    >>> angle_between((1, 0, 0), (-1, 0, 0))
    3.141592653589793

    """
    v1_u = unit_vector(v1)
    v2_u = unit_vector(v2)
    # Don't use `np.dot`, does not work with all shapes
    angle = np.arccos(np.inner(v1_u, v2_u))
    return angle


def innerprod_1d(v1, v2):
    """1d Inner product for vector-of-vectors.

    Example:
    ========
    ```
    v1 = np.array([dir_x, dir_y, dir_z]).T
    v2 = ... # dito
    angle_between_v1_v1 = innerprod_1d(v1, v2)
    ```
    """
    return np.einsum('ij,ij->i', v1, v1)


def unit_vector(vector, **kwargs):
    """Returns the unit vector of the vector."""
    # This also works for a dataframe with columns ['x', 'y', 'z']
    # However, the division operation is picky about the shapes
    # So, remember input vector shape, cast all up to 2d,
    # do the (ugly) conversion, then return unit in same shape as input
    # of course, the numpy-ized version of the input...
    vector = np.array(vector)
    out_shape = vector.shape
    vector = np.atleast_2d(vector)
    unit = vector / linalg.norm(vector, axis=1, **kwargs)[:, None]
    return unit.reshape(out_shape)


def pld3(p1, p2, d2):
    """Calculate the point-line-distance for given point and line."""
    c = np.cross(d2, p2 - p1)
    n1 = scipy.linalg.norm(c)
    n2 = scipy.linalg.norm(d2)
    return n1 / n2


def lpnorm(x, p=2):
    return np.power(np.sum(np.power(x, p)), 1 / p)


def dist(x1, x2):
    np.linalg.norm(x2 - x1, axis=1)


def com(points, masses=None):
    """Calculate center of mass for given points.
    If masses is not set, assume equal masses."""
    if masses is None:
        return np.average(points, axis=0)
    else:
        return np.average(points, axis=0, weights=masses)
{'fraction_non_alphanumeric': 0.14251849042290915, 'fraction_numerical': 0.0333775839180732, 'mean_word_length': 2.6685217391304348, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 33, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2172723', 'n_tokens_mistral': 4209, 'n_tokens_neox': 3836, 'n_words': 1130}