code: def _compute(self, funcTilde, R, z, phi):
          
          Acos, Asin = self._Acos, self._Asin
          N, L, M = Acos.shape    
          r, theta, phi = bovy_coords.cyl_to_spher(R,z,phi)
          
          
     
          PP = lpmn(M-1,L-1,nu.cos(theta))[0].T ##Get the Legendre polynomials
          func_tilde = funcTilde(r, N, L) ## Tilde of the function of interest 
          
          func = nu.zeros((N,L,M), float) ## The function of interest (density or potential)
          
          m = nu.arange(0, M)[nu.newaxis, nu.newaxis, :]
          mcos = nu.cos(m*phi)
          msin = nu.sin(m*phi)
          func = func_tilde[:,:,None]*(Acos[:,:,:]*mcos + Asin[:,:,:]*msin)*PP[None,:,:]
          return func
docs: NAME:
             _compute
          PURPOSE:
             evaluate the NxLxM density or potential
          INPUT:
             funcTidle - must be _rhoTilde or _phiTilde
             R - Cylindrical Galactocentric radius
             z - vertical height
             phi - azimuth
          OUTPUT:
             An NxLxM density or potential at (R,z, phi)
          HISTORY:
             2016-05-18 - Written - Aladdin