# Path: attitude/geom/util.py
# def dot(*matrices):
#     return reduce(N.dot, matrices)
# 
# Path: attitude/stereonet.py
# def sph2cart(lat,lon):
#     _ = stereonet_math.sph2cart(lat,lon)
#     #val = N.array(_).flatten()
#     val = N.roll(_,-1)
#     val[:-1] *= -1
#     return val
# 
# Path: attitude/error/axes.py
# def sampling_covariance(fit, **kw):
#     # This is asserted in both Faber and Jolliffe, although the
#     # former expression is ambiguous due to a weirdly-typeset radical
#     ev = fit.eigenvalues
#     cov = 2/(fit.n-1)*ev**2
#     return cov
# 
# def sampling_axes(fit, **kw):
#     """
#     Hyperbolic axis lengths based on sample-size
#     normal statistics.
# 
#     Note: legacy method that simply wraps the
#     `statistical_axes` function.
#     """
#     return statistical_axes(fit, method='sampling', **kw)


import numpy as N
from math import factorial
from scipy.special import gamma, hyp1f1
from scipy.linalg import expm
from itertools import product
from .geom.util import dot
from .stereonet import sph2cart
from .error.axes import sampling_covariance, sampling_axes
from mplstereonet.stereonet_math import _rotate, cart2sph, sph2cart
    Eigenvalues are analogous to |R|/N.
    """
    # Uses eigenvectors of the covariance matrix
    e = fit.hyperbolic_axes #singular_values
    #e = sampling_covariance(fit) # not sure
    e = e[2]**2/e

    kappa = (e-e[2])[:-1]
    kappa /= kappa[-1]
    F = N.sqrt(N.pi)*confluent_hypergeometric_function(*kappa)


    ax = fit.axes

    Z = 1/e
    M = ax
    F = 1/hyp1f1(*1/Z)

    def pdf(coords):
        lon,lat = coords

        I = lat
        D = lon# + N.pi/2
        #D,I = _rotate(N.degrees(D),N.degrees(I),90)

        # Bingham is given in spherical coordinates of inclination
        # and declination in radians

        # From USGS bingham statistics reference


# Path: attitude/geom/util.py
# def dot(*matrices):
#     return reduce(N.dot, matrices)
# 
# Path: attitude/stereonet.py
# def sph2cart(lat,lon):
#     _ = stereonet_math.sph2cart(lat,lon)
#     #val = N.array(_).flatten()
#     val = N.roll(_,-1)
#     val[:-1] *= -1
#     return val
# 
# Path: attitude/error/axes.py
# def sampling_covariance(fit, **kw):
#     # This is asserted in both Faber and Jolliffe, although the
#     # former expression is ambiguous due to a weirdly-typeset radical
#     ev = fit.eigenvalues
#     cov = 2/(fit.n-1)*ev**2
#     return cov
# 
# def sampling_axes(fit, **kw):
#     """
#     Hyperbolic axis lengths based on sample-size
#     normal statistics.
# 
#     Note: legacy method that simply wraps the
#     `statistical_axes` function.
#     """
#     return statistical_axes(fit, method='sampling', **kw)


# Path: attitude/bingham.py
import numpy as N
from math import factorial
from scipy.special import gamma, hyp1f1
from scipy.linalg import expm
from itertools import product
from .geom.util import dot
from .stereonet import sph2cart
from .error.axes import sampling_covariance, sampling_axes
from mplstereonet.stereonet_math import _rotate, cart2sph, sph2cart

    Eigenvalues are analogous to |R|/N.
    """
    # Uses eigenvectors of the covariance matrix
    e = fit.hyperbolic_axes #singular_values
    #e = sampling_covariance(fit) # not sure
    e = e[2]**2/e

    kappa = (e-e[2])[:-1]
    kappa /= kappa[-1]
    F = N.sqrt(N.pi)*confluent_hypergeometric_function(*kappa)


    ax = fit.axes

    Z = 1/e
    M = ax
    F = 1/hyp1f1(*1/Z)

    def pdf(coords):
        lon,lat = coords

        I = lat
        D = lon# + N.pi/2
        #D,I = _rotate(N.degrees(D),N.degrees(I),90)

        # Bingham is given in spherical coordinates of inclination
        # and declination in radians

        # From USGS bingham statistics reference


        xhat = N.array(sph2cart(lon,lat)).T