Produces Ynm(theta,phi)

        theta runs from 0 to pi and has 'nrows' points.

        phi runs from 0 to 2*pi - 2*pi/ncols and has 'ncols' points.
557750
def sph_harmonic_tp(nrows, ncols, n, m):
    
    nuvec = pysphi.ynunm(n, m, n + 1)
    num = nuvec[1:] * (-1) ** m
    num = num[::-1]
    yvec = np.concatenate([num, nuvec])
    
    theta = np.linspace(0.0, np.pi, nrows)
    phi = np.linspace(0.0, 2.0 * np.pi - 2.0 * np.pi / ncols)
    nu = np.array(range(-n, n + 1), dtype=np.complex128)
    
    out = np.zeros([nrows, ncols], dtype=np.complex128)

    for nn in xrange(0, nrows):
        exp_vec = np.exp(1j * nu * theta[nn])
        for mm in xrange(0, ncols):          
            vl = np.sum(yvec * exp_vec)
            out[nn, mm] = 1j ** m * np.exp(1j * m * phi[mm]) * vl

    return out
795457
def rdirichlet(theta, size=1):
    
    gammas = np.vstack([rgamma(theta, 1) for i in xrange(size)])
    if size > 1 and np.size(theta) > 1:
        return (gammas.T / gammas.sum(1))[:-1].T
    elif np.size(theta) > 1:
        return (gammas[0] / gammas[0].sum())[:-1]
    else:
        return 1.
1014793
def rot1(theta):
    
    return np.array([
        [1, 0, 0],
        [0, np.cos(theta), np.sin(theta)],
        [0, -np.sin(theta), np.cos(theta)]
    ])
941878
def arc_distance(theta_1, phi_1,
                       theta_2, phi_2):
    
    temp = np.sin((theta_2-theta_1)/2)**2+np.cos(theta_1)*np.cos(theta_2)*np.sin((phi_2-phi_1)/2)**2
    distance_matrix = 2 * (np.arctan2(np.sqrt(temp),np.sqrt(1-temp)))
    return distance_matrix
731456
def to_cartesian(r, theta, theta_units="radians"):
    
    assert theta_units in ['radians', 'degrees'],\
        "kwarg theta_units must specified in radians or degrees"

    # Convert to radians
    if theta_units == "degrees":
        theta = to_radians(theta)

    theta = to_proper_radians(theta)
    x = r * cos(theta)
    y = r * sin(theta)

    return x, y
1046587
def _theta_E_convert(self, theta_E, w_c, w_t):
        
        if not w_t > w_c:
            w_t, w_c = w_c, w_t
        s_scale_1 = w_c
        s_scale_2 = w_t
        f_x_1, f_y_1 = self.nie.derivatives(1, 0, theta_E=1, e1=0, e2=0, s_scale=s_scale_1)
        f_x_2, f_y_2 = self.nie.derivatives(1, 0, theta_E=1, e1=0, e2=0, s_scale=s_scale_2)
        f_x = f_x_1 - f_x_2
        theta_E_convert = theta_E / f_x
        return theta_E_convert, w_c, w_t
610978
def compute(self, xt1, yt1, xt, yt, theta1t1, theta2t1, theta1, theta2):
    
    dx = xt - xt1
    dy = yt - yt1

    if self.numPoints < self.maxPoints:
      self.dxValues[self.numPoints,0] = dx
      self.dxValues[self.numPoints,1] = dy

      self.thetaValues[self.numPoints,0] = theta1
      self.thetaValues[self.numPoints,1] = theta2

      self.numPoints += 1

      # print >>sys.stderr, "Xt's: ", xt1, yt1, xt, yt, "Delta's: ", dx, dy
      # print >>sys.stderr, "Theta t-1: ", theta1t1, theta2t1, "t:",theta1, theta2

    elif self.numPoints == self.maxPoints:
      print >> sys.stderr,"Max points exceeded, analyzing ",self.maxPoints,"points only"
      self.numPoints += 1
1186889
def real_rotated_Gammas(m,l,phi1,phi2,theta1,theta2,gamma_ml):
    
    
    if m>0:
        ans=(1./sqrt(2))*(rotated_Gamma_ml(m,l,phi1,phi2,theta1,theta2,gamma_ml) + \
                          (-1)**m*rotated_Gamma_ml(-m,l,phi1,phi2,theta1,theta2,gamma_ml))
        return ans.real
    if m==0:
        return rotated_Gamma_ml(0,l,phi1,phi2,theta1,theta2,gamma_ml).real
    if m<0:
        ans=(1./sqrt(2)/complex(0.,1))*(rotated_Gamma_ml(-m,l,phi1,phi2,theta1,theta2,gamma_ml) - \
                                        (-1)**m*rotated_Gamma_ml(m,l,phi1,phi2,theta1,theta2,gamma_ml))
        return ans.real
813993
def r_y(theta, r = Const::R_UNIT)
      r_a = Array.new

      begin
        s, c = Math.sin(theta), Math.cos(theta)
        a00 = c * r[0][0] - s * r[2][0]
        a01 = c * r[0][1] - s * r[2][1]
        a02 = c * r[0][2] - s * r[2][2]
        a20 = s * r[0][0] + c * r[2][0]
        a21 = s * r[0][1] + c * r[2][1]
        a22 = s * r[0][2] + c * r[2][2]
        r_a << [a00, a01, a02]
        r_a << r[1]
        r_a << [a20, a21, a22]
        return r_a
      rescue => e
        raise
      end
    end
1152826
def run(xmin, ymin, xmax, ymax, step, range_, range_x, range_y, t):
    pt = zeros((range_x, range_y, 3))
    ""
    for i in xrange(range_x):
        for j in xrange(range_y):
            pt[i,j,0], pt[i,j,1] = (xmin+step*i)*180/math.pi, (ymin+step*j)*180/math.pi
            for k in xrange(t.shape[0]):
                tmp = 6368.* math.acos( math.cos(xmin+step*i)*math.cos( t[k,0] ) * math.cos((ymin+step*j)-t[k,1])+  math.sin(xmin+step*i)*math.sin(t[k,0]))
                if tmp < range_:
                    pt[i,j,2]+= t[k,2] / (1+tmp)
    return pt
602971
def phi(self):
        
        return np.arctan2(np.sqrt(self.x**2 + self.y**2), self.z)
459067
def _Aterm(N, h, m, k, dW):
    """"""
    sqrt2h = np.sqrt(2.0/h)
    Xk = np.random.normal(0.0, 1.0, (N, m, 1))
    Yk = np.random.normal(0.0, 1.0, (N, m, 1))
    term1 = _dot(Xk, _t(Yk + sqrt2h*dW))
    term2 = _dot(Yk + sqrt2h*dW, _t(Xk))
    return (term1 - term2)/k
623900
def spher2cart(rho, theta, phi):
    """"""
    st = np.sin(theta)
    sp = np.sin(phi)
    ct = np.cos(theta)
    cp = np.cos(phi)
    rhost = rho * st
    x = rhost * cp
    y = rhost * sp
    z = rho * ct
    return np.array([x, y, z])
586597
def source_to_neutrino_direction(azimuth, zenith, radian=True):
    
    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
509415
def nnz(self):
        
        nnz = 0
        for k in range(len(self.snpost)):
            nn = self.snptr[k+1]-self.snptr[k]    
            na = self.relptr[k+1]-self.relptr[k]
            nnz += nn*(nn+1)/2 + nn*na
        return nnz
781682
def rotated_Gamma_ml(m,l,phi1,phi2,theta1,theta2,gamma_ml):
    
    
    rotated_gamma = 0
    
    for ii in range(2*l+1):
        rotated_gamma += Dlmk(l,m,ii-l,phi1,phi2,theta1,theta2).conjugate()*gamma_ml[ii]

    return rotated_gamma
1016985
function(theta) {
            while (theta < 0) {
                theta += mathUtils.TWO_PI;
            }
            return this.sinLUT[((theta * this.rad2deg) + this.quadrant) % this.period];
        }