Convert lon / lat (radians) for the spherical triangulation into x,y,z
    on the unit sphere
173362
def lonlat2xyz(lon, lat):
    

    lons = np.array(lon)
    lats = np.array(lat)

    xs = np.cos(lats) * np.cos(lons)
    ys = np.cos(lats) * np.sin(lons)
    zs = np.sin(lats)

    return xs, ys, zs
837177
def _prepare_coords(lons1, lats1, lons2, lats2):
    
    lons1 = numpy.radians(lons1)
    lats1 = numpy.radians(lats1)
    assert lons1.shape == lats1.shape
    lons2 = numpy.radians(lons2)
    lats2 = numpy.radians(lats2)
    assert lons2.shape == lats2.shape
    return lons1, lats1, lons2, lats2
582283
function(lon, lat) {
        var gLatLng;
        if(this.sphericalMercator) {
            var lonlat = this.inverseMercator(lon, lat);
            gLatLng = new GLatLng(lonlat.lat, lonlat.lon);
        } else {
            gLatLng = new GLatLng(lat, lon);
        }
        return gLatLng;
    }
117859
def split_coords_3d(seq):
    
    lons, lats, depths = [], [], []
    for i, el in enumerate(seq):
        if i % 3 == 0:
            lons.append(valid.longitude(el))
        elif i % 3 == 1:
            lats.append(valid.latitude(el))
        elif i % 3 == 2:
            depths.append(valid.depth(el))
    return list(zip(lons, lats, depths))
1194797
def interpolate_lat_lon(points, t, max_dt=1):
    '''
    
    '''
    # find the enclosing points in sorted list
    if (t <= points[0][0]) or (t >= points[-1][0]):
        if t <= points[0][0]:
            dt = abs((points[0][0] - t).total_seconds())
        else:
            dt = (t - points[-1][0]).total_seconds()
        if dt > max_dt:
            raise ValueError(
                "time t not in scope of gpx file by {} seconds".format(dt))
        else:
            print(
                "time t not in scope of gpx file by {} seconds, extrapolating...".format(dt))

        if t < points[0][0]:
            before = points[0]
            after = points[1]
        else:
            before = points[-2]
            after = points[-1]
        bearing = compute_bearing(before[1], before[2], after[1], after[2])

        if t == points[0][0]:
            x = points[0]
            return (x[1], x[2], bearing, x[3])

        if t == points[-1][0]:
            x = points[-1]
            return (x[1], x[2], bearing, x[3])
    else:
        for i, point in enumerate(points):
            if t < point[0]:
                if i > 0:
                    before = points[i - 1]
                else:
                    before = points[i]
                after = points[i]
                break

    # weight based on time
    weight = (t - before[0]).total_seconds() / \
        (after[0] - before[0]).total_seconds()

    # simple linear interpolation in case points are not the same
    if before[1] == after[1]:
        lat = before[1]
    else:
        lat = before[1] - weight * before[1] + weight * after[1]

    if before[2] == after[2]:
        lon = before[2]
    else:
        lon = before[2] - weight * before[2] + weight * after[2]

    # camera angle
    bearing = compute_bearing(before[1], before[2], after[1], after[2])

    # altitude
    if before[3] is not None:
        ele = before[3] - weight * before[3] + weight * after[3]
    else:
        ele = None

    return lat, lon, bearing, ele
505101
def interpolate_cubic(self, lons, lats, data):
        
        return self.interpolate(lons, lats, data, order=3)
959316
def lonlat_to_screen(self, lon, lat):
        
        if type(lon) == list:
            lon = np.array(lon)
        if type(lat) == list:
            lat = np.array(lat)

        lat_rad = np.radians(lat)
        n = 2.0 ** self.zoom
        xtile = (lon + 180.0) / 360.0 * n
        ytile = (1.0 - np.log(np.tan(lat_rad) + (1 / np.cos(lat_rad))) / math.pi) / 2.0 * n
        x = (xtile * TILE_SIZE).astype(int)
        y = (SCREEN_H - ytile * TILE_SIZE).astype(int)
        return x, y
266100
def sphere_average(dset,x,y,z,radius=1):
    ''''''
    return_list = []
    if isinstance(dset,basestring):
        dset = [dset]
    for d in dset:
        return_list += [float(a) for a in subprocess.check_output(['3dmaskave','-q','-dball',str(x),str(y),str(z),str(radius),d],stderr=subprocess.PIPE).split()]
    return return_list
819716
def get_sky_coords(self):
        
        if self._ipix is None:
            theta, phi = hp.pix2ang(
                self._nside, list(range(self._npix)), self._nest)
        else:
            theta, phi = hp.pix2ang(self._nside, self._ipix, self._nest)

        lat = np.degrees((np.pi / 2) - theta)
        lon = np.degrees(phi)
        return np.vstack([lon, lat]).T
179152
def _spherical_to_cartesian(cls, coord, center):
        
        r, theta, phi, r_dot, theta_dot, phi_dot = coord
        x = r * cos(phi) * cos(theta)
        y = r * cos(phi) * sin(theta)
        z = r * sin(phi)

        vx = r_dot * x / r - y * theta_dot - z * phi_dot * cos(theta)
        vy = r_dot * y / r + x * theta_dot - z * phi_dot * sin(theta)
        vz = r_dot * z / r + r * phi_dot * cos(phi)

        return np.array([x, y, z, vx, vy, vz], dtype=float)
67602
def _expand_coord(self, lat, lon, lmax_calc, degrees):
        """"""
        if self.normalization == '4pi':
            norm = 1
        elif self.normalization == 'schmidt':
            norm = 2
        elif self.normalization == 'unnorm':
            norm = 3
        elif self.normalization == 'ortho':
            norm = 4
        else:
            raise ValueError(
                "Normalization must be '4pi', 'ortho', 'schmidt', or " +
                "'unnorm'. Input value was {:s}"
                .format(repr(self.normalization)))

        if degrees is True:
            latin = lat
            lonin = lon
        else:
            latin = _np.rad2deg(lat)
            lonin = _np.rad2deg(lon)

        if type(lat) is not type(lon):
            raise ValueError('lat and lon must be of the same type. ' +
                             'Input types are {:s} and {:s}'
                             .format(repr(type(lat)), repr(type(lon))))

        if type(lat) is int or type(lat) is float or type(lat) is _np.float_:
            return _shtools.MakeGridPoint(self.coeffs, lat=latin, lon=lonin,
                                          lmax=lmax_calc, norm=norm,
                                          csphase=self.csphase)
        elif type(lat) is _np.ndarray:
            values = _np.empty_like(lat, dtype=float)
            for v, latitude, longitude in _np.nditer([values, latin, lonin],
                                                     op_flags=['readwrite']):
                v[...] = _shtools.MakeGridPoint(self.coeffs, lat=latitude,
                                                lon=longitude,
                                                lmax=lmax_calc, norm=norm,
                                                csphase=self.csphase)
            return values
        elif type(lat) is list:
            values = []
            for latitude, longitude in zip(latin, lonin):
                values.append(
                    _shtools.MakeGridPoint(self.coeffs, lat=latitude,
                                           lon=longitude,
                                           lmax=lmax_calc, norm=norm,
                                           csphase=self.csphase))
            return values
        else:
            raise ValueError('lat and lon must be either an int, float, ' +
                             'ndarray, or list. ' +
                             'Input types are {:s} and {:s}'
                             .format(repr(type(lat)), repr(type(lon))))
457363
function sphericalCoordinates(pos){
    var rad = len(pos);
    var inc = Math.acos(pos.z / rad) / Math.PI;
    var ang = Math.atan2(pos.y, pos.x) / (Math.PI*2);
    return {rad: rad, ang: ang, inc: inc};
  }
903802
def get_pixel_indices(self, lats, lons):
        
        theta = np.radians(90. - lats)
        phi = np.radians(lons)
        return hp.ang2pix(self.nside, theta, phi, self.nest)
1096018
def from_center(self, x=None, y=None, z=None, r=None,
                    theta=None, h=None, reference=None):
        
        coords_to_endpoint = None
        if all([isinstance(i, numbers.Number) for i in (x, y, z)]):
            coords_to_endpoint = self.from_cartesian(x, y, z)

        if all([isinstance(i, numbers.Number) for i in (r, theta, h)]):
            coords_to_endpoint = self.from_polar(r, theta, h)

        coords_to_reference = Vector(0, 0, 0)
        if reference:
            coords_to_reference = self.coordinates(reference)

        return coords_to_reference + coords_to_endpoint
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
593353
def interpolate(self, lat, lon, var):
        

        subset, dims = self.crop(lat, lon, var)

        if np.all([y in dims['lat'] for y in lat]) & \
                np.all([x in dims['lon'] for x in lon]):
                    yn = np.nonzero([y in lat for y in dims['lat']])[0]
                    xn = np.nonzero([x in lon for x in dims['lon']])[0]
                    output = {}
                    for v in subset:
                        # output[v] = subset[v][dn, zn, yn, xn]
                        # Seriously that this is the way to do it?!!??
                        output[v] = subset[v][:, xn][yn]
                    return output

        # The output coordinates shall be created only once.
        points_out = []
        for latn in lat:
            for lonn in lon:
                points_out.append([latn, lonn])
        points_out = np.array(points_out)

        output = {}
        for v in var:
            output[v] = ma.masked_all(
                    (lat.size, lon.size),
                    dtype=subset[v].dtype)

            # The valid data
            idx = np.nonzero(~ma.getmaskarray(subset[v]))

            if idx[0].size > 0:
                points = np.array([
                    dims['lat'][idx[0]], dims['lon'][idx[1]]]).T
                values = subset[v][idx]

                # Interpolate along the dimensions that have more than one
                #   position, otherwise it means that the output is exactly
                #   on that coordinate.
                ind = np.array(
                        [np.unique(points[:, i]).size > 1 for i in
                            range(points.shape[1])])
                assert ind.any()

                values_out = griddata(
                        np.atleast_1d(np.squeeze(points[:, ind])),
                        values,
                        np.atleast_1d(np.squeeze(points_out[:, ind]))
                        )

                # Remap the interpolated value back into a 4D array
                idx = np.isfinite(values_out)
                for [y, x], out in zip(points_out[idx], values_out[idx]):
                    output[v][y==lat, x==lon] = out

        return output
231558
def xyz2lonlat(x, y, z):
    """"""
    lon = xu.rad2deg(xu.arctan2(y, x))
    lat = xu.rad2deg(xu.arctan2(z, xu.sqrt(x**2 + y**2)))
    return lon, lat