# Path: openmodes/constants.py



import numpy as np
import scipy.special
from scipy.special import factorial
from openmodes.constants import eta_0

    P_lm = np.zeros((len(ct), len(l), m.shape[1]))
    dP_lm = np.zeros_like(P_lm)

    for point_count, ct_n in enumerate(ct):

        # associated Legendre function and its derivative
        P_lmp, dP_lmp = scipy.special.lpmn(max_l, max_l, ct_n)
        P_lmp = P_lmp.T
        dP_lmp = dP_lmp.T

        # Comes from http://dlmf.nist.gov/14.9, eq 14.9.13
        # Calculate negative values of m from positive
        # Does some relation apply to derivatives??
        P_neg = (-1)**m_pos*factorial(l-m_pos)/factorial(l+m_pos)
        P_lmn = P_neg*P_lmp
        dP_lmn = P_neg*dP_lmp

        # combine positive and negative P_lmn
        P_lm[point_count] = np.hstack((P_lmp,  P_lmn[:, :0:-1]))
        dP_lm[point_count] = np.hstack((dP_lmp, dP_lmn[:, :0:-1]))

    # theta derivative of P_lm(cos\theta)
    tau_lm = -st[:, None, None]*dP_lm
    pi_lm = P_lm*m/st[:, None, None]

    return (r, theta, phi, r_hat, theta_hat, phi_hat, P_lm, exp_imp,
            tau_lm, pi_lm)



# Path: openmodes/constants.py



# Path: openmodes/multipole.py
import numpy as np
import scipy.special
from scipy.special import factorial
from openmodes.constants import eta_0


    P_lm = np.zeros((len(ct), len(l), m.shape[1]))
    dP_lm = np.zeros_like(P_lm)

    for point_count, ct_n in enumerate(ct):

        # associated Legendre function and its derivative
        P_lmp, dP_lmp = scipy.special.lpmn(max_l, max_l, ct_n)
        P_lmp = P_lmp.T
        dP_lmp = dP_lmp.T

        # Comes from http://dlmf.nist.gov/14.9, eq 14.9.13
        # Calculate negative values of m from positive
        # Does some relation apply to derivatives??
        P_neg = (-1)**m_pos*factorial(l-m_pos)/factorial(l+m_pos)
        P_lmn = P_neg*P_lmp
        dP_lmn = P_neg*dP_lmp

        # combine positive and negative P_lmn
        P_lm[point_count] = np.hstack((P_lmp,  P_lmn[:, :0:-1]))
        dP_lm[point_count] = np.hstack((dP_lmp, dP_lmn[:, :0:-1]))

    # theta derivative of P_lm(cos\theta)
    tau_lm = -st[:, None, None]*dP_lm
    pi_lm = P_lm*m/st[:, None, None]

    return (r, theta, phi, r_hat, theta_hat, phi_hat, P_lm, exp_imp,
            tau_lm, pi_lm)



def spherical_multipoles(max_l, k, points, current, current_M, eta=eta_0,