# -*- coding: utf-8 -*-
"""
Created on Tue Oct 14 17:07:23 2014

@author: jussi
"""

def GaV(l,m,r,theta,phi):
    pass

def GbV(l,m,r,theta,phi):
    pass

def YlmV(l,m,theta,phi):
    """ Compute spherical harmonic function Ylm with orders l,m at the point (theta,phi) in spherical coordinates. theta is polar angle from z axis, phi is azimuth.
        """
        
def Plm(l,m,x):
        """ Legendre associated polynomial. Handles negative m. x can be a 1-dim. array.
        """
        # handle negative m via multiplying constant
   if m  < 0:
            C = (-1)**abs(m)*fct(l-abs(m))/fct(l+abs(m))
            m = abs(m)
        else:
            C = 1
        xlen = x.shape[0]
        # init variables
        cx = np.zeros([xlen,l+1])  # store cases from l=0 to l=1
        cx[:,m] = 1
        sqm2 = sqrt(1-x**2)
        # compute l' == m
        D=1  # multiplier
        for j in range(1,m+1):
            cx[:,m] = -cx[:,m]*D*sqm2
            D=D+2
        # if l==m, we are ready
        if m == l:
            return C*cx[:,l]
        # case l' == m+1
        cx[:,m+1] = x*(2*m+1)*cx[:,m]            
        # l'=m+2 to l'=l by recursion
        for j in range(m+2,l+1):
            cx[:,j] = ((2*j-1)*x*cx[:,j-1]+(-j-m+1)*cx[:,j-2])/(j-m)
        return C*cx[:,l]
        
    
def DPlm(self,l,m,x):
        """ Derivative of Legendre associated polynomial Plm. Handles negative m. x can be a 1-dim. array. """
        if m<0:
            C = (-1)**abs(m)*fct(l-abs(m))/fct(l+abs(m))
            m = abs(m)
        else:

            C = 1    
        return C*(m*x*self.Plm(l,m,x)+(l+m)*(l-m+1)*sqrt(1-x**2)*self.Plm(l,m-1,x))/(1-x**2)
        
