Return the angular distance from the given direction and
    the edge of the projection.
507618
def distance_to_edge(geom, skydir):
    """"""

    # FIXME: We should add a pixel_size property in gammapy.maps
    # FIXME: We should make this into a MapGeom method 
    
    xpix, ypix = skydir.to_pixel(geom.wcs, origin=0)
    deltax = np.array((xpix - geom.center_pix[0]) * geom._cdelt[0],
                      ndmin=1)
    deltay = np.array((ypix - geom.center_pix[1]) * geom._cdelt[1],
                      ndmin=1)

    deltax = np.abs(deltax) - 0.5 * geom.width[0]
    deltay = np.abs(deltay) - 0.5 * geom.width[1]

    m0 = (deltax < 0) & (deltay < 0)
    m1 = (deltax > 0) & (deltay < 0)
    m2 = (deltax < 0) & (deltay > 0)
    m3 = (deltax > 0) & (deltay > 0)
    mx = np.abs(deltax) <= np.abs(deltay)
    my = np.abs(deltay) < np.abs(deltax)

    delta = np.zeros(len(deltax))
    delta[(m0 & mx) | (m3 & my) | m1] = deltax[(m0 & mx) | (m3 & my) | m1]
    delta[(m0 & my) | (m3 & mx) | m2] = deltay[(m0 & my) | (m3 & mx) | m2]
    return delta
692215
def geo_distance(left, right, use_spheroid=None):
    
    op = ops.GeoDistance(left, right, use_spheroid)
    return op.to_expr()
509079
def angle2d(self):
        """"""
        if self.x==0:
            if   self.y<0:  return       math.pi/2.0*3
            elif self.y>0:  return       math.pi/2.0
            else:           return 0
        elif self.y==0:
            if   self.x<0:  return       math.pi
           #elif self.x>0:  return 0
            else:           return 0
        ans = math.atan( self.y / self.x )
        if self.x > 0:
            if self.y>0:    return ans
            else:           return ans + math.pi*2.0
        else:               return ans + math.pi
664572
private static function computeAngleBetween($from, $to) {
        return self::distanceRadians(deg2rad($from['lat']), deg2rad($from['lng']),
                               deg2rad($to['lat']), deg2rad($to['lng']));
    }
87664
def angular_distance(km, lat, lat2=None):
    
    if lat2 is not None:
        # use the largest latitude to compute the angular distance
        lat = max(abs(lat), abs(lat2))
    return km * KM_TO_DEGREES / math.cos(lat * DEGREES_TO_RAD)
803346
def distances(self, points):
        """"""
        return [math.acos(self._pos3d.dot(p.vector)) for p in points]
392902
def point_to_t(self, point):
        """"""

        def in_range(min, max, val):
            return (min <= val) and (max >= val)

        # Single-precision floats have only 7 significant figures of
        # resolution, so test that we're within 6 sig figs.
        if np.isclose(point, self.start, rtol=0.0, atol=1e-6):
            return 0.0
        elif np.isclose(point, self.end, rtol=0.0, atol=1e-6):
            return 1.0

        if self.rotation != 0.0:
            raise ValueError("Arc.point_to_t() only works on non-rotated Arcs.")

        v = point - self.center
        distance_from_center = sqrt((v.real * v.real) + (v.imag * v.imag))
        min_radius = min(self.radius.real, self.radius.imag)
        max_radius = max(self.radius.real, self.radius.imag)
        if (distance_from_center < min_radius) and not np.isclose(distance_from_center, min_radius):
            return None
        if (distance_from_center > max_radius) and not np.isclose(distance_from_center, max_radius):
            return None

        # x = center_x + radius_x cos(radians(theta + t delta))
        # y = center_y + radius_y sin(radians(theta + t delta))
        #
        # For x:
        # cos(radians(theta + t delta)) = (x - center_x) / radius_x
        # radians(theta + t delta) = acos((x - center_x) / radius_x)
        # theta + t delta = degrees(acos((x - center_x) / radius_x))
        # t_x = (degrees(acos((x - center_x) / radius_x)) - theta) / delta
        #
        # Similarly for y:
        # t_y = (degrees(asin((y - center_y) / radius_y)) - theta) / delta

        x = point.real
        y = point.imag

        #
        # +Y points down!
        #
        # sweep mean clocwise
        # sweep && (delta > 0)
        # !sweep && (delta < 0)
        #
        # -180 <= theta_1 <= 180
        #
        # large_arc && (-360 <= delta <= 360)
        # !large_arc && (-180 < delta < 180)
        #

        end_angle = self.theta + self.delta
        min_angle = min(self.theta, end_angle)
        max_angle = max(self.theta, end_angle)

        acos_arg = (x - self.center.real) / self.radius.real
        if acos_arg > 1.0:
            acos_arg = 1.0
        elif acos_arg < -1.0:
            acos_arg = -1.0

        x_angle_0 = degrees(acos(acos_arg))
        while x_angle_0 < min_angle:
            x_angle_0 += 360.0
        while x_angle_0 > max_angle:
            x_angle_0 -= 360.0

        x_angle_1 = -1.0 * x_angle_0
        while x_angle_1 < min_angle:
            x_angle_1 += 360.0
        while x_angle_1 > max_angle:
            x_angle_1 -= 360.0

        t_x_0 = (x_angle_0 - self.theta) / self.delta
        t_x_1 = (x_angle_1 - self.theta) / self.delta

        asin_arg = (y - self.center.imag) / self.radius.imag
        if asin_arg > 1.0:
            asin_arg = 1.0
        elif asin_arg < -1.0:
            asin_arg = -1.0

        y_angle_0 = degrees(asin(asin_arg))
        while y_angle_0 < min_angle:
            y_angle_0 += 360.0
        while y_angle_0 > max_angle:
            y_angle_0 -= 360.0

        y_angle_1 = 180 - y_angle_0
        while y_angle_1 < min_angle:
            y_angle_1 += 360.0
        while y_angle_1 > max_angle:
            y_angle_1 -= 360.0

        t_y_0 = (y_angle_0 - self.theta) / self.delta
        t_y_1 = (y_angle_1 - self.theta) / self.delta

        t = None
        if np.isclose(t_x_0, t_y_0):
            t = (t_x_0 + t_y_0) / 2.0
        elif np.isclose(t_x_0, t_y_1):
            t= (t_x_0 + t_y_1) / 2.0
        elif np.isclose(t_x_1, t_y_0):
            t = (t_x_1 + t_y_0) / 2.0
        elif np.isclose(t_x_1, t_y_1):
            t = (t_x_1 + t_y_1) / 2.0
        else:
            # Comparing None and float yields a result in python2,
            # but throws TypeError in python3.  This fix (suggested by
            # @CatherineH) explicitly handles and avoids the case where
            # the None-vs-float comparison would have happened below.
            return None

        if (t >= 0.0) and (t <= 1.0):
            return t

        return None
167018
def direction(self):
        
        libfn = utils.get_lib_fn('getDirection%s'%self._libsuffix)
        return libfn(self.pointer)
228397
def _principal_angle(a, B):
  ''''''
  # TODO: check case for t = d-1
  if a.shape[1] == 1:
    return a.T.dot(B)[0,:,0]

  # find normals that maximize distance when projected
  x1 = np.einsum('abc,adc->abd', B, B).dot(a) - a   # b.dot(b.T).dot(a) - a
  x2 = np.einsum('ab,cad->cbd', a.dot(a.T), B) - B  # a.dot(a.T).dot(b) - b
  xx = np.vstack((x1, x2))

  # batch PCA (1st comp. only)
  xx -= xx.mean(axis=1)[:,None]
  c = np.einsum('abc,abd->acd', xx, xx)
  _, vecs = np.linalg.eigh(c)
  fpc = vecs[:,:,-1]
  fpc1 = fpc[:len(x1)]
  fpc2 = fpc[len(x1):]

  # a.dot(fpc1).dot(b.dot(fpc2))
  lhs = a.dot(fpc1.T).T
  rhs = np.einsum('abc,ac->ab', B, fpc2)
  return np.einsum('ij,ij->i', lhs, rhs)
1083608
def distance(self, x, y):
        
        if scipy.sparse.issparse(x):
            return numpy.linalg.norm((x-y).toarray().ravel())
        else:
            return numpy.linalg.norm(x-y)
792009
def distance_to(self, *args):
        """"""
        if not len(args) > 0:
            raise ValueError("at least one parameter must be given")

        x = args[0]
        if len(args) > 1:
            unit = args[1]
        else:
            unit = "px"
        u = {"px": (1.,1.), "in": (1.,72.), "cm": (2.54, 72.),
             "mm": (25.4, 72.)}
        f = u[unit][0] / u[unit][1]
        if type(x) is Point:
            return abs(self - x) * f

# from here on, x is a rectangle
# as a safeguard, make a finite copy of it
        r = Rect(x.top_left, x.top_left)
        r = r | x.bottom_right
        if self in r:
            return 0.0
        if self.x > r.x1:
            if self.y >= r.y1:
                return self.distance_to(r.bottom_right, unit)
            elif self.y <= r.y0:
                return self.distance_to(r.top_right, unit)
            else:
                return (self.x - r.x1) * f
        elif r.x0 <= self.x <= r.x1:
            if self.y >= r.y1:
                return (self.y - r.y1) * f
            else:
                return (r.y0 - self.y) * f
        else:
            if self.y >= r.y1:
                return self.distance_to(r.bottom_left, unit)
            elif self.y <= r.y0:
                return self.distance_to(r.top_left, unit)
            else:
                return (r.x0 - self.x) * f
1048504
def distance(self, other):
        """"""
        return math.acos(self._pos3d.dot(other.vector))
36074
function(obj) {
    if (obj.normal) { return obj.isParallelTo(this); }
    var theta = this.direction.angleFrom(obj.direction);
    return (Math.abs(theta) <= Sylvester.precision || Math.abs(theta - Math.PI) <= Sylvester.precision);
  }
1230782
def edge_cluster_angle(edge_dirs, subspaces1, subspaces2):
  ''''''
  QG = edge_dirs / np.linalg.norm(edge_dirs, ord=2, axis=1)[:,None]
  X1 = np.einsum('...ij,...i->...j', subspaces1, QG)
  X2 = np.einsum('...ij,...i->...j', subspaces2, QG)
  # TODO: check the math on this for more cases
  # angles = np.maximum(1-np.sum(X1**2, axis=1), 1-np.sum(X2**2, axis=1))
  C1 = np.linalg.svd(X1[:,:,None], compute_uv=False)
  C2 = np.linalg.svd(X2[:,:,None], compute_uv=False)
  angles = np.maximum(1-C1**2, 1-C2**2)[:,0]
  angles[np.isnan(angles)] = 0.0  # nan when edge length == 0
  return angles
478012
def interp(self):
        """"""
        if self.ndim == 0:
            return 'nearest'
        elif all(interp == self.interp_byaxis[0]
                 for interp in self.interp_byaxis):
            return self.interp_byaxis[0]
        else:
            return self.interp_byaxis
1078447
def projection(self, **kwargs):
        
        viewpoint = kwargs.pop('viewpoint', None)
        if viewpoint is None:
            ndim = self.shape[0]-1
            viewpoint = N.zeros(ndim)
        plane = self.polar_plane(viewpoint)
        return self.slice(plane, **kwargs)
1290350
def get_position(self):
        
        with self._lock:
            return self.x, self.y, self.angle