Python math.acos() Examples

The following are code examples for showing how to use math.acos(). They are extracted from open source Python projects. You can vote up the examples you like or vote down the exmaples you don't like. You can also save this page to your account.

Example 1
Project: camera_calibration_frontend   Author: groundmelon   File: calibrator.py    (license) View Source Project 7 votes vote down vote up
def _get_skew(corners, board):
    """
    Get skew for given checkerboard detection.
    Scaled to [0,1], which 0 = no skew, 1 = high skew
    Skew is proportional to the divergence of three outside corners from 90 degrees.
    """
    # TODO Using three nearby interior corners might be more robust, outside corners occasionally
    # get mis-detected
    up_left, up_right, down_right, _ = _get_outside_corners(corners, board)

    def angle(a, b, c):
        """
        Return angle between lines ab, bc
        """
        ab = a - b
        cb = c - b
        return math.acos(numpy.dot(ab,cb) / (numpy.linalg.norm(ab) * numpy.linalg.norm(cb)))

    skew = min(1.0, 2. * abs((math.pi / 2.) - angle(up_left, up_right, down_right)))
    return skew 
Example 2
Project: bpy_lambda   Author: bcongdon   File: utils.py    (license) View Source Project 7 votes vote down vote up
def align_bone_z_axis(obj, bone, vec):
    """ Rolls the bone to align its z-axis as closely as possible to
        the given vector.
        Must be in edit mode.
    """
    bone_e = obj.data.edit_bones[bone]

    vec = bone_e.y_axis.cross(vec)
    vec.normalize()

    dot = max(-1.0, min(1.0, bone_e.x_axis.dot(vec)))
    angle = math.acos(dot)

    bone_e.roll += angle

    dot1 = bone_e.x_axis.dot(vec)

    bone_e.roll -= angle * 2

    dot2 = bone_e.x_axis.dot(vec)

    if dot1 > dot2:
        bone_e.roll += angle * 2 
Example 3
Project: pybot   Author: spillai   File: quaternion.py    (license) View Source Project 6 votes vote down vote up
def interpolate(self, other, this_weight):
        q0, q1 = np.roll(self.q, shift=1), np.roll(other.q, shift=1)
        u = 1 - this_weight
        assert(u >= 0 and u <= 1)
        cos_omega = np.dot(q0, q1)

        if cos_omega < 0:
            result = -q0[:]
            cos_omega = -cos_omega
        else:
            result = q0[:]

        cos_omega = min(cos_omega, 1)

        omega = math.acos(cos_omega)
        sin_omega = math.sin(omega)
        a = math.sin((1-u) * omega)/ sin_omega
        b = math.sin(u * omega) / sin_omega

        if abs(sin_omega) < 1e-6:
            # direct linear interpolation for numerically unstable regions
            result = result * this_weight + q1 * u
            result /= math.sqrt(np.dot(result, result))
        else:
            result = result*a + q1*b
        return Quaternion(np.roll(result, shift=-1))

    # To conversions 
Example 4
Project: ssbio   Author: SBRG   File: cpv.py    (MIT License) View Source Project 6 votes vote down vote up
def get_angle_formed_by(p1,p2,p3): # angle formed by three positions in space

    # based on code submitted by Paul Sherwood
    r1 = distance(p1,p2)
    r2 = distance(p2,p3)
    r3 = distance(p1,p3)
    
    small = 1.0e-10
    
    if (r1 + r2 - r3) < small:
        # This seems to happen occasionally for 180 angles 
        theta = math.pi
    else:
        theta = math.acos( (r1*r1 + r2*r2  - r3*r3) / (2.0 * r1*r2) )
    return theta;

#------------------------------------------------------------------------------ 
Example 5
Project: CPNSimulatorGui   Author: chris-kuhr   File: ArcItem.py    (GNU General Public License v3.0) View Source Project 6 votes vote down vote up
def setPolygon(self):
        '''Calculate position and rotation of the arc arrow head.'''
        rotDeg = 0            
        xlength = self.pos1.x() - self.pos2.x()
        ylength = self.pos1.y() - self.pos2.y()
        d = math.sqrt( math.pow( xlength , 2) + math.pow( ylength , 2) )        
        if d > 0:
            beta = math.acos( xlength / d )
            rotDeg = math.degrees( beta ) 
        
        self.arrowPolygonObject.setPolygon( QtGui.QPolygonF( [
                                                 QtCore.QPointF( (self.pos2.x() -10),  (self.pos2.y() +5)), 
                                                 QtCore.QPointF( (self.pos2.x() -10) , (self.pos2.y() -5)), 
                                                 QtCore.QPointF(       self.pos2.x() ,      self.pos2.y())
                                                 ] ) ) 
        
        self.arrowPolygonObject.setBrush( QtGui.QBrush(QtCore.Qt.black) )
        
        """ self.angle()!!!!!!!!!"""
#         self.arcLinePolygon.angle()
#         self.arcLinePolygon.rotate(rotDeg)
#         self.arcLinePolygon.setPos( self.pos2 ) 

    #------------------------------------------------------------------------------------------------ 
Example 6
Project: CPNSimulatorGui   Author: chris-kuhr   File: ArrowGui.py    (GNU General Public License v3.0) View Source Project 6 votes vote down vote up
def setPolygon(self):
        rotDeg = 0            
        xlength = self.pos1.x() - self.pos2.x()
        ylength = self.pos1.y() - self.pos2.y()
        d = math.sqrt( math.pow( xlength , 2) + math.pow( ylength , 2) )        
        if d > 0:
            beta = math.acos( xlength / d )
            rotDeg = math.degrees( beta ) 
        
        self.arcLinePolygon.setPolygon( QtGui.QPolygonF( [
                                                 QtCore.QPointF( (self.pos2.x() -10),  (self.pos2.y() +5)), 
                                                 QtCore.QPointF( (self.pos2.x() -10) , (self.pos2.y() -5)), 
                                                 QtCore.QPointF(       self.pos2.x() ,      self.pos2.y())
                                                 ] ) ) 
        
        self.arcLinePolygon.setBrush( QtGui.QBrush(QtCore.Qt.black) )
        
        """ self.angle()!!!!!!!!!"""
#         self.arcLinePolygon.angle()
#         self.arcLinePolygon.rotate(rotDeg)
#         self.arcLinePolygon.setPos( self.pos2 ) 
    #------------------------------------------------------------------------------------------------ 
Example 7
Project: Gaia   Author: splcurran   File: operators.py    (license) View Source Project 6 votes vote down vote up
def eLowDotOperator(stack, z, mode):
	if mode == 1:   # num
		stack.append(utilities.formatNum(math.acos(z)))
	#elif mode == 2:
	elif mode == 3: # str or list
		if len(z) == 0:
			stack.append([])
		else:
			result = ""
			for i in z:
				i = utilities.castToList(i)
				if len(i) >= 2:
					if type(i[1]) == str and type(result) == str:
						result += (i[1] * utilities.castToNumber(i[0]))
					else:
						result = list(result)
						result += [i[1]] * utilities.castToNumber(i[0])
			stack.append(result)
	else:
		monadNotImplemented(mode, '')

# ? 
Example 8
Project: gprime   Author: GenealogyCollective   File: fanchart.py    (license) View Source Project 6 votes vote down vote up
def get_max_width_for_circles(self, rad1, rad2, max_centering_proportion):
        r"""
        (the "r" in the above line is to keep pylint happy)
           __
          /__\ <- compute the line width which is drawable between 2 circles.
         /  _ \   max_centering_proportion : 0, touching the circle1, 1,
        |  |_| |   touching the circle2, 0.5 : middle between the 2 circles
        |      |
         \    /
          \__/
                 basically, max_centering_proportion is
                 max_centering_proportion/nb_lines
        """
        # radius at the center of the 2 circles
        rmid = rad2 - (rad2-rad1)*max_centering_proportion
        return sin(acos(rmid/rad2)) * rad2 * 2 
Example 9
Project: bpy_lambda   Author: bcongdon   File: art2polyarea.py    (license) View Source Project 6 votes vote down vote up
def _Angle(u, v):
    """Return angle between two vectors.

    Args:
      u: (float, float)
      v: (float, float)
    Returns:
      float - angle in radians between u and v, where
        it is +/- depending on sign of ux * vy - uy * vx
    """

    (ux, uy) = u
    (vx, vy) = v
    costheta = (ux * vx + uy * vy) / \
        (math.sqrt(ux ** 2 + uy ** 2) * math.sqrt(vx ** 2 + vy ** 2))
    if costheta > 1.0:
        costheta = 1.0
    if costheta < -1.0:
        costheta = -1.0
    theta = math.acos(costheta)
    if ux * vy - uy * vx < 0.0:
        theta = -theta
    return theta 
Example 10
Project: bpy_lambda   Author: bcongdon   File: triquad.py    (license) View Source Project 6 votes vote down vote up
def Angle(a, b, c, points):
    """Return Angle abc in degrees, in range [0,180),
    where a,b,c are indices into points."""

    u = Sub2(points.pos[c], points.pos[b])
    v = Sub2(points.pos[a], points.pos[b])
    n1 = Length2(u)
    n2 = Length2(v)
    if n1 == 0.0 or n2 == 0.0:
        return 0.0
    else:
        costheta = Dot2(u, v) / (n1 * n2)
        if costheta > 1.0:
            costheta = 1.0
        if costheta < - 1.0:
            costheta = - 1.0
        return math.acos(costheta) * 180.0 / math.pi 
Example 11
Project: bpy_lambda   Author: bcongdon   File: archipack_2d.py    (license) View Source Project 6 votes vote down vote up
def proj_xy(self, t, next=None):
        """
            length of projection of sections at crossing line / circle intersections
            deformation unit vector for profil in xy axis
            so f(x_profile) = position of point in xy plane
        """
        if next is None:
            return self.normal(t).v.normalized(), 1
        v0 = self.normal(1).v.normalized()
        v1 = next.normal(0).v.normalized()
        direction = v0 + v1
        adj = (v0 * self.length) * (v1 * next.length)
        hyp = (self.length * next.length)
        c = min(1, max(-1, adj / hyp))
        size = 1 / cos(0.5 * acos(c))
        return direction.normalized(), min(3, size) 
Example 12
Project: bpy_lambda   Author: bcongdon   File: io_export_md3.py    (license) View Source Project 6 votes vote down vote up
def Encode(self, normal):
        x = normal[0]
        y = normal[1]
        z = normal[2]
        # normalize
        l = math.sqrt((x*x) + (y*y) + (z*z))
        if l == 0:
            return 0
        x = x/l
        y = y/l
        z = z/l

        if (x == 0.0) & (y == 0.0) :
            if z > 0.0:
                return 0
            else:
                return (128 << 8)

        lng = math.acos(z) * 255 / (2 * math.pi)
        lat = math.atan2(y, x) * 255 / (2 * math.pi)
        retval = ((int(lat) & 0xFF) << 8) | (int(lng) & 0xFF)
        return retval 
Example 13
Project: bpy_lambda   Author: bcongdon   File: triquad.py    (license) View Source Project 6 votes vote down vote up
def Angle(a, b, c, points):
    """Return Angle abc in degrees, in range [0,180),
    where a,b,c are indices into points."""

    u = Sub2(points.pos[c], points.pos[b])
    v = Sub2(points.pos[a], points.pos[b])
    n1 = Length2(u)
    n2 = Length2(v)
    if n1 == 0.0 or n2 == 0.0:
        return 0.0
    else:
        costheta = Dot2(u, v) / (n1 * n2)
        if costheta > 1.0:
            costheta = 1.0
        if costheta < - 1.0:
            costheta = - 1.0
        return math.acos(costheta) * 180.0 / math.pi 
Example 14
Project: bpy_lambda   Author: bcongdon   File: delta.py    (license) View Source Project 6 votes vote down vote up
def set_mat(obj, bone_name, matrix):
        """ Sets the bone to have the given transform matrix.
        """
        a = obj.data.edit_bones[bone_name]

        a.head = (0, 0, 0)
        a.tail = (0, 1, 0)

        a.transform(matrix)

        d = acos(a.matrix.to_quaternion().dot(matrix.to_quaternion())) * 2.0

        roll_1 = a.roll + d
        roll_2 = a.roll - d

        a.roll = roll_1
        d1 = a.matrix.to_quaternion().dot(matrix.to_quaternion())
        a.roll = roll_2
        d2 = a.matrix.to_quaternion().dot(matrix.to_quaternion())

        if d1 > d2:
            a.roll = roll_1
        else:
            a.roll = roll_2 
Example 15
Project: bpy_lambda   Author: bcongdon   File: utils.py    (license) View Source Project 6 votes vote down vote up
def angle_on_plane(plane, vec1, vec2):
    """ Return the angle between two vectors projected onto a plane.
    """
    plane.normalize()
    vec1 = vec1 - (plane * (vec1.dot(plane)))
    vec2 = vec2 - (plane * (vec2.dot(plane)))
    vec1.normalize()
    vec2.normalize()

    # Determine the angle
    angle = math.acos(max(-1.0, min(1.0, vec1.dot(vec2))))

    if angle < 0.00001:  # close enough to zero that sign doesn't matter
        return angle

    # Determine the sign of the angle
    vec3 = vec2.cross(vec1)
    vec3.normalize()
    sign = vec3.dot(plane)
    if sign >= 0:
        sign = 1
    else:
        sign = -1

    return angle * sign 
Example 16
Project: bpy_lambda   Author: bcongdon   File: utils.py    (license) View Source Project 6 votes vote down vote up
def align_bone_x_axis(obj, bone, vec):
    """ Rolls the bone to align its x-axis as closely as possible to
        the given vector.
        Must be in edit mode.
    """
    bone_e = obj.data.edit_bones[bone]

    vec = vec.cross(bone_e.y_axis)
    vec.normalize()

    dot = max(-1.0, min(1.0, bone_e.z_axis.dot(vec)))
    angle = math.acos(dot)

    bone_e.roll += angle

    dot1 = bone_e.z_axis.dot(vec)

    bone_e.roll -= angle * 2

    dot2 = bone_e.z_axis.dot(vec)

    if dot1 > dot2:
        bone_e.roll += angle * 2 
Example 17
Project: bpy_lambda   Author: bcongdon   File: delta.py    (license) View Source Project 6 votes vote down vote up
def set_mat(obj, bone_name, matrix):
        """ Sets the bone to have the given transform matrix.
        """
        a = obj.data.edit_bones[bone_name]

        a.head = (0, 0, 0)
        a.tail = (0, 1, 0)

        a.transform(matrix)

        d = acos(a.matrix.to_quaternion().dot(matrix.to_quaternion())) * 2.0

        roll_1 = a.roll + d
        roll_2 = a.roll - d

        a.roll = roll_1
        d1 = a.matrix.to_quaternion().dot(matrix.to_quaternion())
        a.roll = roll_2
        d2 = a.matrix.to_quaternion().dot(matrix.to_quaternion())

        if d1 > d2:
            a.roll = roll_1
        else:
            a.roll = roll_2 
Example 18
Project: bpy_lambda   Author: bcongdon   File: utils.py    (license) View Source Project 6 votes vote down vote up
def angle_on_plane(plane, vec1, vec2):
    """ Return the angle between two vectors projected onto a plane.
    """
    plane.normalize()
    vec1 = vec1 - (plane * (vec1.dot(plane)))
    vec2 = vec2 - (plane * (vec2.dot(plane)))
    vec1.normalize()
    vec2.normalize()

    # Determine the angle
    angle = math.acos(max(-1.0, min(1.0, vec1.dot(vec2))))

    if angle < 0.00001:  # close enough to zero that sign doesn't matter
        return angle

    # Determine the sign of the angle
    vec3 = vec2.cross(vec1)
    vec3.normalize()
    sign = vec3.dot(plane)
    if sign >= 0:
        sign = 1
    else:
        sign = -1

    return angle * sign 
Example 19
Project: DHP   Author: YuhangSong   File: move_view_lib.py    (license) View Source Project 6 votes vote down vote up
def get_sph_cor(x,y,z):    #using radian
    lon=0
    if(x==0 and y>0):
        lon=PI/2
    elif(x==0 and y<0):
        lon=3*PI/2
    elif(x==0 and y==0):
        print ("error")
        return
    elif(x>0 and y==0):
        lon=0
    elif(x>0 and y>0):
        lon=math.atan(float(y)/float(x))
    elif(x>0 and y<0):
        lon=2*PI+math.atan(float(y)/float(x))
    elif(x<0 and y==0):
        lon=PI
    elif(x<0 and y>0):
        lon=PI+math.atan(float(y)/float(x))
    elif(x<0 and y<0):
        lon=PI+math.atan(float(y)/float(x))

    lat=PI/2-math.acos(z/1.0)

    return lon,lat 
Example 20
Project: SPIND   Author: LiuLab-CSRC   File: SPIND.py    (license) View Source Project 6 votes vote down vote up
def calc_rotation_matrix(q1, q2, ref_q1, ref_q2):
  ref_nv = np.cross(ref_q1, ref_q2) 
  q_nv = np.cross(q1, q2)
  if min(norm(ref_nv), norm(q_nv)) == 0.:  # avoid 0 degree including angle
    return np.identity(3)
  axis = np.cross(ref_nv, q_nv)
  angle = rad2deg(acos(ref_nv.dot(q_nv) / (norm(ref_nv) * norm(q_nv))))
  R1 = axis_angle_to_rotation_matrix(axis, angle)
  rot_ref_q1, rot_ref_q2 = R1.dot(ref_q1), R1.dot(ref_q2)  # rotate ref_q1,2 plane to q1,2 plane

  cos1 = max(min(q1.dot(rot_ref_q1) / (norm(rot_ref_q1) * norm(q1)), 1.), -1.)  # avoid math domain error
  cos2 = max(min(q2.dot(rot_ref_q2) / (norm(rot_ref_q2) * norm(q2)), 1.), -1.)
  angle1 = rad2deg(acos(cos1))
  angle2 = rad2deg(acos(cos2))
  angle = (angle1 + angle2) / 2.
  axis = np.cross(rot_ref_q1, q1)
  R2 = axis_angle_to_rotation_matrix(axis, angle)

  R = R2.dot(R1)
  return R 
Example 21
Project: cadee   Author: kamerlinlab   File: q_ffld2q.py    (license) View Source Project 6 votes vote down vote up
def calc_angle(ac1, ac2, ac3, fk, th0):
    th0 = math.pi/180.0 * th0   # degrees to radians
    rji1 = ac1[0] - ac2[0]
    rji2 = ac1[1] - ac2[1]
    rji3 = ac1[2] - ac2[2]
    rjk1 = ac3[0] - ac2[0]
    rjk2 = ac3[1] - ac2[1]
    rjk3 = ac3[2] - ac2[2]
    
    # get the angle from the dot product equation ( A*B = |A|*|B|*cos(theta) )
    # where A and B are vectors bewteen atoms (1->2 and 2->3)
    bji2inv = 1./(rji1**2 + rji2**2 + rji3**2)
    bjk2inv = 1./(rjk1**2 + rjk2**2 + rjk3**2)
    bjiinv = math.sqrt(bji2inv)
    bjkinv = math.sqrt(bjk2inv)
    scp = (rji1*rjk1 + rji2*rjk2 + rji3*rjk3)
    scp = scp * bjiinv* bjkinv
    if scp > 1.0:
        scp =  1.0
    elif scp < -1.0:
        scp = -1.0
    theta = math.acos(scp)   # in radians
    dtheta = (theta-th0)  
    en = 0.5*fk*dtheta**2 
    return en,theta*180.0/math.pi    # kcal/mol and degrees 
Example 22
Project: zippy   Author: securesystemslab   File: plot_rotation.py    (license) View Source Project 6 votes vote down vote up
def get_spherical_rotatation(p1, p2, width, height, theta_multiplier):
    v1 = get_sphere_mapping(p1[0], p1[1], width, height)
    v2 = get_sphere_mapping(p2[0], p2[1], width, height)

    d = min(max([dot(v1, v2), -1]), 1)

    if abs(d - 1.0) < 0.000001:
        return None

    raxis = norm( cross(v1, v2) )
    rtheta = theta_multiplier * rad2deg * _acos(d)
    #rtheta = 2.0 * rad2deg * _acos(d)

    glPushMatrix()
    glLoadIdentity()
    glRotatef(rtheta, *raxis)
    mat = (c_float*16)()
    glGetFloatv(GL_MODELVIEW_MATRIX, mat)
    glPopMatrix()

    return mat 
Example 23
Project: cvcalib   Author: Algomorph   File: quaternion.py    (license) View Source Project 6 votes vote down vote up
def _get_angle_axis(self):
        lim = 1e-12
        norm = np.linalg.norm(self.q)
        if norm < lim:
            angle = 0
            axis = [0, 0, 0]
        else:
            rnorm = 1.0 / norm
            angle = acos(max(-1, min(1, rnorm * self.q[3])));
            sangle = sin(angle)
            if sangle < lim:
                axis = [0, 0, 0]
            else:
                axis = (rnorm / sangle) * np.array(self.q[0:3])

            angle *= 2

        return (angle, axis) 
Example 24
Project: cvcalib   Author: Algomorph   File: quaternion.py    (license) View Source Project 6 votes vote down vote up
def slerp(qa, qb, t):
        # Calculate angle between them.
        #qa.w * qb.w + qa.x * qb.x + qa.y * qb.y + qa.z * qb.z;
        cosHalfTheta = np.dot(qa.q,qb.q)
        #if qa=qb or qa=-qb then theta = 0 and we can return qa
        if (abs(cosHalfTheta) >= 1.0):
            return Quat(np.copy(qa.q));
        #Calculate temporary values.
        halfTheta = acos(cosHalfTheta);
        sinHalfTheta = sqrt(1.0 - cosHalfTheta*cosHalfTheta);
        
        #if theta = 180 degrees then result is not fully defined
        #we could rotate around any axis normal to qa or qb
        if(abs(sinHalfTheta) < 0.001):
            return Quat(qa.q * 0.5 + qb.q * 0.5);
        
        ratioA = sin((1 - t) * halfTheta) / sinHalfTheta;
        ratioB = sin(t * halfTheta) / sinHalfTheta; 
        #calculate Quaternion for general case.
        return Quat(qa.q * ratioA + qb.q * ratioB); 
Example 25
Project: blender-cnt   Author: bcorso   File: blender_cnt.py    (license) View Source Project 6 votes vote down vote up
def __init__(self, lattice, m, n):
        self.lattice = lattice
        self.m = m
        self.n = n

        # Chiral vector
        self.c = lattice.pos(m,n)
        self.magC = mag(self.c)
        # Translation vector 
        d = gcd(2*n+m,2*m+n)
        self.t = lattice.pos((2*n+m)/d, -(2*m+n)/d);
        self.magT = mag(self.t)
        # Chiral rotation matrix (rotate a1 along x-axis)
        self.theta = acos(norm(self.c)[0]*copysign(1, self.c[1]))
        self.rotM = np.array([
                [cos(self.theta), sin(self.theta)],
                [-sin(self.theta), cos(self.theta)]]).T

        # Calculate atoms and bonds in unit cell
        self._boundsErr = mag(lattice.pos(0,0,0) - lattice.pos(0,0,1))
        self.indices = self._calcIndices(m, n)
        self.atoms = self._calcAtoms(self.indices)
        self.bonds = self._calcBonds(self.indices) 
Example 26
Project: srcsim2017   Author: ZarjRobotics   File: tf_conv.py    (license) View Source Project 6 votes vote down vote up
def compute_distance_to_wall(a, r1, r2):
    """ Given a sighting that determines a distance of r1 to the wall, and then a rotation by angle *a*
        to the left and a sighting of distance r2, returns the angle that we are currently rotated by
        from the perpendicular to the wall and the current distance to the wall as a pair tuple. Angle
        is given in degrees before the rotation from r1 to r2. Rotating to the right by this angle
        should cause us to face the wall directly.
    """
    try:
        if r1 < r2:
            r = r1/r2
            i = 1.0
        elif r1 > r2:
            r = r2/r1
            i = -1.0
        else:
            return 0, r2
        d = radians(a)
        c = cos(d)
        s = sin(d)
        dt = sqrt(1 - c * c * r * r)
        x = c * c * r + s * dt
        return degrees(acos(x)) * i, r2 * x
    except ValueError:
        traceback.print_exc()
        return None, None 
Example 27
Project: FC   Author: JanusWind   File: janus_fc_dat.py    (license) View Source Project 6 votes vote down vote up
def calc_eff_area( self, v ) :

		# Note. #nvn is a vector similar to inflow particle bulk 
		# velocity and ndir is the look direction. 


		# Normalize the particle velocity.

		vn  = calc_arr_norm( v )
		nvn = tuple( [ -c for c in vn ] )

		# Calculate the particle inflow angle (in degrees) relative to
		# the cup normal (i.e., the cup pointing direction).

		psi = acos( calc_arr_dot( self['dir'], nvn ) )*pi/180.
		if ( psi > 90. ) :
			return 0. 
 		
		# Return the effective collecting area corresponding to "psi".

		return interp( psi, self._spec._eff_deg, self._spec._eff_area )

	#-----------------------------------------------------------------------
	# DEFINE THE FUNCTION TO CALCULATE EXPECTED MAXWELLIAN CURRENT.
	#----------------------------------------------------------------------- 
Example 28
Project: kicad_scripts   Author: NilujePerchut   File: td.py    (license) View Source Project 6 votes vote down vote up
def __ComputeCurved(vpercent, w, vec, via, pts, segs):
    """Compute the curves part points"""

    radius = via[1]/2.0

    # Compute the bezier middle points
    req_angle = asin(vpercent/100.0)
    oppside = tan(req_angle)*(radius-(w/sin(req_angle)))
    length = sqrt(radius*radius + oppside*oppside)
    d = req_angle - acos(radius/length)
    vecBC = [vec[0]*cos(d)+vec[1]*sin(d), -vec[0]*sin(d)+vec[1]*cos(d)]
    pointBC = via[0] + wxPoint(int(vecBC[0] * length), int(vecBC[1] * length))
    d = -d
    vecAE = [vec[0]*cos(d)+vec[1]*sin(d), -vec[0]*sin(d)+vec[1]*cos(d)]
    pointAE = via[0] + wxPoint(int(vecAE[0] * length), int(vecAE[1] * length))

    curve1 = __Bezier(pts[1], pointBC, pts[2], n=segs)
    curve2 = __Bezier(pts[4], pointAE, pts[0], n=segs)

    return curve1 + [pts[3]] + curve2 
Example 29
Project: roLabelImg   Author: cgvict   File: canvas.py    (license) View Source Project 6 votes vote down vote up
def getAngle(self, center, p1, p2):
        dx1 = p1.x() - center.x();
        dy1 = p1.y() - center.y();

        dx2 = p2.x() - center.x();
        dy2 = p2.y() - center.y();

        c = math.sqrt(dx1*dx1 + dy1*dy1) * math.sqrt(dx2*dx2 + dy2*dy2)
        if c == 0: return 0
        y = (dx1*dx2+dy1*dy2)/c
        if y>1: return 0
        angle = math.acos(y)

        if (dx1*dy2-dx2*dy1)>0:   
            return angle
        else:
            return -angle 
Example 30
Project: seaglass   Author: seaglass-project   File: misc.py    (license) View Source Project 6 votes vote down vote up
def earth_distance(lat1, lon1, lat2, lon2):
    """ Distance in meters between two points specified in degrees. """
    x1 = calc_rad(lat1) * math.cos(degree_to_radian(lon1)) * math.sin(degree_to_radian(90-lat1))
    x2 = calc_rad(lat2) * math.cos(degree_to_radian(lon2)) * math.sin(degree_to_radian(90-lat2))
    y1 = calc_rad(lat1) * math.sin(degree_to_radian(lon1)) * math.sin(degree_to_radian(90-lat1))
    y2 = calc_rad(lat2) * math.sin(degree_to_radian(lon2)) * math.sin(degree_to_radian(90-lat2))
    z1 = calc_rad(lat1) * math.cos(degree_to_radian(90-lat1))
    z2 = calc_rad(lat2) * math.cos(degree_to_radian(90-lat2))
    a = (x1*x2 + y1*y2 + z1*z2)/pow(calc_rad((lat1+lat2)/2), 2)
    # a should be in [1, -1] but can sometimes fall outside it by
    # a very small amount due to rounding errors in the preceding
    # calculations (this is prone to happen when the argument points
    # are very close together).  Thus we constrain it here.
    if abs(a) > 1:
        a = 1
    elif a < -1:
        a = -1
    return calc_rad((lat1+lat2) / 2) * math.acos(a) 
Example 31
Project: gpvdm   Author: roderickmackenzie   File: spectrum_solargeometry.py    (license) View Source Project 6 votes vote down vote up
def get_zenith(Latitude, Longitude, d, hour, minute, timezone):
    gamma_val = ((2 * math.pi) / 365) * ((d - 1) + (hour - 12) / 24)
    decl_angle = 0.006918 - (0.399912 * math.cos(gamma_val)) + 0.070257 * math.sin(gamma_val) - 0.006758 * math.cos(
        2 * gamma_val) \
                 + 0.000907 * math.sin(2 * gamma_val) - 0.002697 * math.cos(3 * gamma_val) + 0.00148 * math.sin(
        3 * gamma_val)
    eq_time = 229.18 * (
        0.000075 + 0.001868 * math.cos(gamma_val) - 0.032077 * math.sin(gamma_val) - 0.014615 * math.cos(2 * gamma_val)
        - 0.040849 * math.sin(2 * gamma_val))
    time_offset = eq_time - 4 * Longitude + 60*timezone
    true_solar_time = hour * 60 + minute + time_offset
    solar_hour_angle = true_solar_time / 4 - 180

    Z_deg = (180/math.pi)*math.acos((math.sin(Latitude * (math.pi / 180)) * math.sin(decl_angle)) + (
        math.cos(Latitude * (math.pi / 180)) * math.cos(decl_angle) * math.cos(solar_hour_angle * (math.pi / 180))))
    return Z_deg 
Example 32
Project: gpvdm   Author: roderickmackenzie   File: solar_solargeometry.py    (license) View Source Project 6 votes vote down vote up
def get_zenith(Latitude, Longitude, d, hour, minute, timezone):
    gamma_val = ((2 * math.pi) / 365) * ((d - 1) + (hour - 12) / 24)
    decl_angle = 0.006918 - (0.399912 * math.cos(gamma_val)) + 0.070257 * math.sin(gamma_val) - 0.006758 * math.cos(
        2 * gamma_val) \
                 + 0.000907 * math.sin(2 * gamma_val) - 0.002697 * math.cos(3 * gamma_val) + 0.00148 * math.sin(
        3 * gamma_val)
    eq_time = 229.18 * (
        0.000075 + 0.001868 * math.cos(gamma_val) - 0.032077 * math.sin(gamma_val) - 0.014615 * math.cos(2 * gamma_val)
        - 0.040849 * math.sin(2 * gamma_val))
    time_offset = eq_time - 4 * Longitude + 60*timezone
    true_solar_time = hour * 60 + minute + time_offset
    solar_hour_angle = true_solar_time / 4 - 180

    Z_deg = (180/math.pi)*math.acos((math.sin(Latitude * (math.pi / 180)) * math.sin(decl_angle)) + (
        math.cos(Latitude * (math.pi / 180)) * math.cos(decl_angle) * math.cos(solar_hour_angle * (math.pi / 180))))
    return Z_deg 
Example 33
Project: 3D-IWGAN   Author: EdwardSmith1884   File: render_model_views.py    (license) View Source Project 6 votes vote down vote up
def camPosToQuaternion(cx, cy, cz):
    camDist = math.sqrt(cx * cx + cy * cy + cz * cz)
    cx = cx / camDist
    cy = cy / camDist
    cz = cz / camDist
    axis = (-cz, 0, cx)
    angle = math.acos(cy)
    a = math.sqrt(2) / 2
    b = math.sqrt(2) / 2
    w1 = axis[0]
    w2 = axis[1]
    w3 = axis[2]
    c = math.cos(angle / 2)
    d = math.sin(angle / 2)
    q1 = a * c - b * d * w1
    q2 = b * c + a * d * w1
    q3 = a * d * w2 + b * d * w3
    q4 = -b * d * w2 + a * d * w3
    return (q1, q2, q3, q4) 
Example 34
Project: pysptools   Author: ctherien   File: dist.py    (license) View Source Project 6 votes vote down vote up
def SAM(s1, s2):
    """
    Computes the spectral angle mapper between two vectors (in radians).

    Parameters:
        s1: `numpy array`
            The first vector.

        s2: `numpy array`
            The second vector.

    Returns: `float`
            The angle between vectors s1 and s2 in radians.
    """
    try:
        s1_norm = math.sqrt(np.dot(s1, s1))
        s2_norm = math.sqrt(np.dot(s2, s2))
        sum_s1_s2 = np.dot(s1, s2)
        angle = math.acos(sum_s1_s2 / (s1_norm * s2_norm))
    except ValueError:
        # python math don't like when acos is called with
        # a value very near to 1
        return 0.0
    return angle 
Example 35
Project: PySCIPOpt   Author: SCIP-Interfaces   File: read_tsplib.py    (license) View Source Project 6 votes vote down vote up
def distGEO(x1,y1,x2,y2):
    print("Implementation is wrong")
    assert False
    PI = 3.141592
    deg = int(x1 + .5)
    min_ = x1 - deg
    lat1 = PI * (deg + 5.*min_/3)/180.
    deg = int(y1 + .5)
    min_ = y1 - deg
    long1 = PI * (deg + 5.*min_/3)/180.
    deg = int(x2 + .5)
    min_ = x2 - deg
    lat2 = PI * (deg + 5.*min_/3)/180.
    deg = int(y2 + .5)
    min_ = y2 - deg
    long2 = PI * (deg + 5.*min_/3)/180.


    RRR = 6378.388
    q1 = math.cos( long1 - long2 );
    q2 = math.cos( lat1 - lat2 );
    q3 = math.cos( lat1 + lat2 );
    return int(RRR * math.acos(.5*((1.+q1)*q2 - (1.-q1)*q3)) + 1.) 
Example 36
Project: PySCIPOpt   Author: SCIP-Interfaces   File: read_tsplib.py    (license) View Source Project 6 votes vote down vote up
def distGEO(x1,y1,x2,y2):
    print("Implementation is wrong")
    assert False
    PI = 3.141592
    deg = int(x1 + .5)
    min_ = x1 - deg
    lat1 = PI * (deg + 5.*min_/3)/180.
    deg = int(y1 + .5)
    min_ = y1 - deg
    long1 = PI * (deg + 5.*min_/3)/180.
    deg = int(x2 + .5)
    min_ = x2 - deg
    lat2 = PI * (deg + 5.*min_/3)/180.
    deg = int(y2 + .5)
    min_ = y2 - deg
    long2 = PI * (deg + 5.*min_/3)/180.


    RRR = 6378.388
    q1 = math.cos( long1 - long2 );
    q2 = math.cos( lat1 - lat2 );
    q3 = math.cos( lat1 + lat2 );
    return int(RRR * math.acos(.5*((1.+q1)*q2 - (1.-q1)*q3)) + 1.) 
Example 37
Project: TauMeta   Author: slspencer   File: patternmath.py    (license) View Source Project 6 votes vote down vote up
def angleOfVector(p1, v, p2):
    """
    Accepts p1,  v,  and p2 of class Point or coordinate pairs
    Returns the angle in radians between the vector v-to-p1 and vector v-to-p2
    """
    #L1=distance(p1, p2)
    #L2=distance(p1, p3)
    #L3=distance(p2, p3)
    #return math.acos((L1**2+L2**2-L3**2)/(2*L1*L2))
    p1 = dPnt(p1)
    p2 = dPnt(p2)
    angle1 = angleOfLine(v, p1)
    angle2 = angleOfLine(v, p2)
    #get the absolute angle
    angle = abs(angle1 - angle2)
    #get the smallest angle of the vector, should not be greater than a straight line
    if angle > pi:
        angle = 2*pi - angle
    return angle 
Example 38
Project: TauMeta   Author: slspencer   File: pmath.py    (license) View Source Project 6 votes vote down vote up
def angleOfVector(p1, v, p2):
    """
    Accepts p1,  v,  and p2 of class Point or coordinate pairs
    Returns the angle in radians between the vector v-to-p1 and vector v-to-p2
    """
    #L1=distance(p1, p2)
    #L2=distance(p1, p3)
    #L3=distance(p2, p3)
    #return math.acos((L1**2+L2**2-L3**2)/(2*L1*L2))
    p1 = dPnt(p1)
    p2 = dPnt(p2)
    angle1 = angleOfLine(v, p1)
    angle2 = angleOfLine(v, p2)
    #get the absolute angle
    angle = abs(angle1 - angle2)
    #get the smallest angle of the vector, should not be greater than a straight line
    if angle > pi:
        angle = 2*pi - angle
    return angle 
Example 39
Project: ted-editor   Author: tarnheld   File: ted-editor.py    (license) View Source Project 6 votes vote down vote up
def dist_sqr_at_l(self, p):
        """return squared distance and arc length at nearest point to p on biarc"""
        p1, t1, cha1 = biarc.point_on_circle(self.p0, self.J, self.t0, p)
        p2, t2, cha2 = biarc.point_on_circle(self.J, self.p1, self.Jt, p)
        c1, k1, a1, l1, c2, k2, a2, l2 = self.circleparameters()
        mind, minl = None, None
        if 0 < t1 < 1:
            mind = la.norm2(p1 - p)
            aa1 = 2 * m.acos(cha1)
            minl = aa1 / abs(k1) if k1 != 0 else t1 * l1
        if 0 < t2 < 1 and (not mind or la.norm2(p2 - p) < mind):
            mind = la.norm2(p2 - p)
            aa2 = 2 * m.acos(cha2)
            minl = l1 + (aa2 / abs(k2) if k2 != 0 else t2 * l2)
        return mind, minl
########################################################################### 
Example 40
Project: archipack   Author: s-leger   File: archipack_2d.py    (license) View Source Project 6 votes vote down vote up
def proj_xy(self, t, next=None):
        """
            length of projection of sections at crossing line / circle intersections
            deformation unit vector for profil in xy axis
            so f(x_profile) = position of point in xy plane
        """
        if next is None:
            return self.normal(t).v.normalized(), 1
        v0 = self.normal(1).v.normalized()
        v1 = next.normal(0).v.normalized()
        direction = v0 + v1
        adj = (v0 * self.length) * (v1 * next.length)
        hyp = (self.length * next.length)
        c = min(1, max(-1, adj / hyp))
        size = 1 / cos(0.5 * acos(c))
        return direction.normalized(), min(3, size) 
Example 41
Project: spockpy   Author: achillesrasquinha   File: gesture.py    (license) View Source Project 6 votes vote down vote up
def _get_defects_count(array, contour, defects, verbose = False):
    ndefects = 0
    
    for i in range(defects.shape[0]):
        s,e,f,_ = defects[i,0]
        beg     = tuple(contour[s][0])
        end     = tuple(contour[e][0])
        far     = tuple(contour[f][0])
        a       = _get_eucledian_distance(beg, end)
        b       = _get_eucledian_distance(beg, far)
        c       = _get_eucledian_distance(end, far)
        angle   = math.acos((b ** 2 + c ** 2 - a ** 2) / (2 * b * c)) * 57
            
        if angle <= 90:
            ndefects = ndefects + 1

            if verbose:
                cv2.circle(array, far, 3, _COLOR_RED, -1)

        if verbose:
            cv2.line(array, beg, end, _COLOR_RED, 1)

    return array, ndefects 
Example 42
Project: Solving-NP-Hard-problems   Author: GauravBh1010tt   File: solver.py    (license) View Source Project 6 votes vote down vote up
def orien(current,initial,others):
    temp=[]
    #print others
    f=0
    for i in others:
        if  i != current:
            #print i,current,initial
            a=length(current,initial)
            b=length(i,initial)
            c=length(current,i)
            d=math.acos((a*a + b*b - c*c) / (2*a*b))
            temp.append((d,f))
        else:
            temp.append((0,f))
        f+=1
    return temp 
Example 43
Project: Solving-NP-Hard-problems   Author: GauravBh1010tt   File: vrp_upload.py    (license) View Source Project 6 votes vote down vote up
def orien(current,initial,others):
    temp=[]
    #print others
    f=0
    for i in others:
        if  i != current:
            #print i,current,initial
            a=length(current,initial)
            b=length(i,initial)
            c=length(current,i)
            d=math.acos((a*a + b*b - c*c) / (2*a*b))
            temp.append((d,f))
        else:
            temp.append((0,f))
        f+=1
    return temp 
Example 44
Project: pybot   Author: spillai   File: quaternion.py    (license) View Source Project 5 votes vote down vote up
def to_angle_axis(self):
        """ Return axis-angle representation """
        q = np.roll(self.q, shift=1)
        halftheta = math.acos(q[0])
        if abs(halftheta) < 1e-12:
            return 0, np.array((0, 0, 1))
        else:
            theta = halftheta * 2
            axis = np.array(q[1:4]) / math.sin(halftheta)
            return theta, axis 
Example 45
Project: pybot   Author: spillai   File: transformations.py    (license) View Source Project 5 votes vote down vote up
def quaternion_slerp(quat0, quat1, fraction, spin=0, shortestpath=True):
    """Return spherical linear interpolation between two quaternions.

    >>> q0 = random_quaternion()
    >>> q1 = random_quaternion()
    >>> q = quaternion_slerp(q0, q1, 0.0)
    >>> numpy.allclose(q, q0)
    True
    >>> q = quaternion_slerp(q0, q1, 1.0, 1)
    >>> numpy.allclose(q, q1)
    True
    >>> q = quaternion_slerp(q0, q1, 0.5)
    >>> angle = math.acos(numpy.dot(q0, q))
    >>> numpy.allclose(2.0, math.acos(numpy.dot(q0, q1)) / angle) or \
        numpy.allclose(2.0, math.acos(-numpy.dot(q0, q1)) / angle)
    True

    """
    q0 = unit_vector(quat0[:4])
    q1 = unit_vector(quat1[:4])
    if fraction == 0.0:
        return q0
    elif fraction == 1.0:
        return q1
    d = numpy.dot(q0, q1)
    if abs(abs(d) - 1.0) < _EPS:
        return q0
    if shortestpath and d < 0.0:
        # invert rotation
        d = -d
        q1 *= -1.0
    angle = math.acos(d) + spin * math.pi
    if abs(angle) < _EPS:
        return q0
    isin = 1.0 / math.sin(angle)
    q0 *= math.sin((1.0 - fraction) * angle) * isin
    q1 *= math.sin(fraction * angle) * isin
    q0 += q1
    return q0 
Example 46
Project: Neural-Networks-for-Inverse-Kinematics   Author: paramrajpura   File: transformations.py    (license) View Source Project 5 votes vote down vote up
def quaternion_slerp(quat0, quat1, fraction, spin=0, shortestpath=True):
    """Return spherical linear interpolation between two quaternions.

    >>> q0 = random_quaternion()
    >>> q1 = random_quaternion()
    >>> q = quaternion_slerp(q0, q1, 0)
    >>> numpy.allclose(q, q0)
    True
    >>> q = quaternion_slerp(q0, q1, 1, 1)
    >>> numpy.allclose(q, q1)
    True
    >>> q = quaternion_slerp(q0, q1, 0.5)
    >>> angle = math.acos(numpy.dot(q0, q))
    >>> numpy.allclose(2, math.acos(numpy.dot(q0, q1)) / angle) or \
        numpy.allclose(2, math.acos(-numpy.dot(q0, q1)) / angle)
    True

    """
    q0 = unit_vector(quat0[:4])
    q1 = unit_vector(quat1[:4])
    if fraction == 0.0:
        return q0
    elif fraction == 1.0:
        return q1
    d = numpy.dot(q0, q1)
    if abs(abs(d) - 1.0) < _EPS:
        return q0
    if shortestpath and d < 0.0:
        # invert rotation
        d = -d
        numpy.negative(q1, q1)
    angle = math.acos(d) + spin * math.pi
    if abs(angle) < _EPS:
        return q0
    isin = 1.0 / math.sin(angle)
    q0 *= math.sin((1.0 - fraction) * angle) * isin
    q1 *= math.sin(fraction * angle) * isin
    q0 += q1
    return q0 
Example 47
Project: sketch-components   Author: ibhubs   File: geometry.py    (license) View Source Project 5 votes vote down vote up
def __init__(self, arg):
        if isinstance(arg, numbers.Real):
            # We precompute sin and cos for rotations
            self.angle = arg
            self.cos = math.cos(self.angle)
            self.sin = math.sin(self.angle)
        elif isinstance(arg, Point):
            # Point angle is the trigonometric angle of the vector
            # [origin, Point]
            pt = arg
            try:
                self.cos = pt.x / pt.length()
                self.sin = pt.y / pt.length()
            except ZeroDivisionError:
                self.cos = 1
                self.sin = 0

            self.angle = math.acos(self.cos)
            if self.sin < 0:
                self.angle = -self.angle
        else:
            raise TypeError("Angle is defined by a number or a Point") 
Example 48
Project: FlipDotWorker   Author: ArduinoHannover   File: sunrise.py    (GNU General Public License v3.0) View Source Project 5 votes vote down vote up
def __calc(self):
  """
  Perform the actual calculations for sunrise, sunset and
  a number of related quantities.
  
  The results are stored in the instance variables
  sunrise_t, sunset_t and solarnoon_t
  """
  timezone = self.timezone # in hours, east is positive
  longitude= self.long     # in decimal degrees, east is positive
  latitude = self.lat      # in decimal degrees, north is positive

  time  = self.time # percentage past midnight, i.e. noon  is 0.5
  day      = self.day     # daynumber 1=1/1/1900
 
  Jday     =day+2415018.5+time-timezone/24 # Julian day
  Jcent    =(Jday-2451545)/36525    # Julian century

  Manom    = 357.52911+Jcent*(35999.05029-0.0001537*Jcent)
  Mlong    = 280.46646+Jcent*(36000.76983+Jcent*0.0003032)%360
  Eccent   = 0.016708634-Jcent*(0.000042037+0.0001537*Jcent)
  Mobliq   = 23+(26+((21.448-Jcent*(46.815+Jcent*(0.00059-Jcent*0.001813))))/60)/60
  obliq    = Mobliq+0.00256*cos(rad(125.04-1934.136*Jcent))
  vary     = tan(rad(obliq/2))*tan(rad(obliq/2))
  Seqcent  = sin(rad(Manom))*(1.914602-Jcent*(0.004817+0.000014*Jcent))+sin(rad(2*Manom))*(0.019993-0.000101*Jcent)+sin(rad(3*Manom))*0.000289
  Struelong= Mlong+Seqcent
  Sapplong = Struelong-0.00569-0.00478*sin(rad(125.04-1934.136*Jcent))
  declination = deg(asin(sin(rad(obliq))*sin(rad(Sapplong))))
  
  eqtime   = 4*deg(vary*sin(2*rad(Mlong))-2*Eccent*sin(rad(Manom))+4*Eccent*vary*sin(rad(Manom))*cos(2*rad(Mlong))-0.5*vary*vary*sin(4*rad(Mlong))-1.25*Eccent*Eccent*sin(2*rad(Manom)))

  hourangle= deg(acos(cos(rad(90.833))/(cos(rad(latitude))*cos(rad(declination)))-tan(rad(latitude))*tan(rad(declination))))

  self.solarnoon_t=(720-4*longitude-eqtime+timezone*60)/1440
  self.sunrise_t  =self.solarnoon_t-hourangle*4/1440
  self.sunset_t   =self.solarnoon_t+hourangle*4/1440 
Example 49
Project: ssbio   Author: SBRG   File: cpv.py    (MIT License) View Source Project 5 votes vote down vote up
def get_angle(v1,v2): # v1,v2 must be unit vectors
    denom = (math.sqrt(((v1[0]*v1[0]) + (v1[1]*v1[1]) + (v1[2]*v1[2]))) *
                math.sqrt(((v2[0]*v2[0]) + (v2[1]*v2[1]) + (v2[2]*v2[2]))))
    if denom>1e-10:
        result = ( (v1[0]*v2[0]) + (v1[1]*v2[1]) + (v1[2]*v2[2]) ) / denom
    else:
        result = 0.0
    result = math.acos(result)
    return result

#------------------------------------------------------------------------------ 
Example 50
Project: autolab_core   Author: BerkeleyAutomation   File: transformations.py    (Apache License 2.0) View Source Project 5 votes vote down vote up
def quaternion_slerp(quat0, quat1, fraction, spin=0, shortestpath=True):
    """Return spherical linear interpolation between two quaternions.

    >>> q0 = random_quaternion()
    >>> q1 = random_quaternion()
    >>> q = quaternion_slerp(q0, q1, 0.0)
    >>> numpy.allclose(q, q0)
    True
    >>> q = quaternion_slerp(q0, q1, 1.0, 1)
    >>> numpy.allclose(q, q1)
    True
    >>> q = quaternion_slerp(q0, q1, 0.5)
    >>> angle = math.acos(numpy.dot(q0, q))
    >>> numpy.allclose(2.0, math.acos(numpy.dot(q0, q1)) / angle) or \
        numpy.allclose(2.0, math.acos(-numpy.dot(q0, q1)) / angle)
    True

    """
    q0 = unit_vector(quat0[:4])
    q1 = unit_vector(quat1[:4])
    if fraction == 0.0:
        return q0
    elif fraction == 1.0:
        return q1
    d = numpy.dot(q0, q1)
    if abs(abs(d) - 1.0) < _EPS:
        return q0
    if shortestpath and d < 0.0:
        # invert rotation
        d = -d
        q1 *= -1.0
    angle = math.acos(d) + spin * math.pi
    if abs(angle) < _EPS:
        return q0
    isin = 1.0 / math.sin(angle)
    q0 *= math.sin((1.0 - fraction) * angle) * isin
    q1 *= math.sin(fraction * angle) * isin
    q0 += q1
    return q0