Precession as a rotation matrix
420234
def precesion(date):  # pragma: no cover
    
    zeta, theta, z = np.deg2rad(_precesion(date))
    return rot3(zeta) @ rot2(-theta) @ rot3(z)
354300
def orbit_y(self, angle):
        '''

        '''
        
        # Subtract pivot point
        self.position -= self.pivot
        
        # Rotate
        rot = rotation_matrix(-angle, self.b)[:3,:3]
        self.position = np.dot(rot, self.position)
        
        # Add again the pivot point
        self.position += self.pivot
        
        self.a = np.dot(rot, self.a)
        self.b = np.dot(rot, self.b)
        self.c = np.dot(rot, self.c)
372182
def set_rotation_matrices(self):
        """"""
        for i in range(len(self.vertices)):
            self.vertices[i].meta['C'] = rotation_matrix(self.vertices[i].meta['axis'][0], 
                                                         self.vertices[i].meta['axis'][1], 
                                                         self.vertices[i].meta['axis'][2], 
                                                         self.vertices[i].meta['axis_order'],
                                                         degrees=True)
            # Todo: invert this by applying angle operations in reverse order
            self.vertices[i].meta['Cinv'] = np.linalg.inv(self.vertices[i].meta['C'])
445508
def rotation_matrix_3D(u, th):
    
    # normalize the axis:
    u = normalize(u)
    # We use the Euler-Rodrigues formula;
    # see https://en.wikipedia.org/wiki/Euler-Rodrigues_formula
    a = math.cos(0.5 * th)
    s = math.sin(0.5 * th)
    (b, c, d) = -s * u
    (a2, b2, c2, d2) = (a*a, b*b, c*c, d*d)
    (bc, ad, ac, ab, bd, cd) = (b*c, a*d, a*c, a*b, b*d, c*d)
    return np.array([[a2 + b2 - c2 - d2, 2*(bc + ad),         2*(bd - ac)],
                     [2*(bc - ad),       a2 + c2 - b2 - d2,   2*(cd + ab)],
                     [2*(bd + ac),       2*(cd - ab),         a2 + d2 - b2 - c2]])
470058
def get_matrix(theta, phi, psi, dx, dy, dz,
                   matrix=np.zeros((4, 4), dtype=DTYPE),
                   angles=np.zeros(3, dtype=DTYPE)):
        
        # NB!: matrix and angles by default are being overwritten on each call
        # thus, only created once at compile time.

        angles[0] = theta
        angles[1] = phi
        angles[2] = psi

        cx, cy, cz = np.cos(angles)
        sx, sy, sz = np.sin(angles)

        rotation = matrix[:3, :3]
        rotation.flat = (cx * cz - sx * cy * sz,
                         cx * sz + sx * cy * cz, sx * sy,
                         -sx * cz - cx * cy * sz,
                         -sx * sz + cx * cy * cz, cx * sy,
                         sy * sz,
                         -sy * cz, cy)

        # Translation component
        matrix[:3, 3] = dx, dy, dz
        matrix[3, 3] = 1.
        return matrix
483074
def rotation_matrix_to_quaternion(rotation_matrix):
    """"""
    invert = (np.linalg.det(rotation_matrix) < 0)
    if invert:
        factor = -1
    else:
        factor = 1
    c2 = 0.25*(factor*np.trace(rotation_matrix) + 1)
    if c2 < 0:
        #print c2
        c2 = 0.0
    c = np.sqrt(c2)
    r2 = 0.5*(1 + factor*np.diagonal(rotation_matrix)) - c2
    #print "check", r2.sum()+c2
    r = np.zeros(3, float)
    for index, r2_comp in enumerate(r2):
        if r2_comp < 0:
            continue
        else:
            row, col = off_diagonals[index]
            if (rotation_matrix[row, col] - rotation_matrix[col, row] < 0):
                r[index] = -np.sqrt(r2_comp)
            else:
                r[index] = +np.sqrt(r2_comp)
    return factor, np.array([c, r[0], r[1], r[2]], float)
1371016
def qrot(vector, quaternion):
    
    t = 2 * np.cross(quaternion[1:], vector)
    v_rot = vector + quaternion[0] * t + np.cross(quaternion[1:], t)
    return v_rot
868346
def _rotate(self, angles, dj_matrix, gm=None, r0=None, omega=None):
        """"""
        if dj_matrix is None:
            dj_matrix = _shtools.djpi2(self.lmax + 1)

        # The coefficients need to be 4pi normalized with csphase = 1
        coeffs = _shtools.SHRotateRealCoef(
            self.to_array(normalization='4pi', csphase=1), angles, dj_matrix)

        # Convert 4pi normalized coefficients to the same normalization
        # as the unrotated coefficients.
        if self.normalization != '4pi' or self.csphase != 1:
            temp = _convert(coeffs, normalization_in='4pi', csphase=1,
                            normalization_out=self.normalization,
                            csphase_out=self.csphase)
            return SHGravCoeffs.from_array(
                temp, normalization=self.normalization,
                csphase=self.csphase, copy=False, gm=gm, r0=r0, omega=omega)
        else:
            return SHGravCoeffs.from_array(coeffs, gm=gm, r0=r0, omega=omega,
                                           copy=False)
972523
def sidereal_time(t):
    """"""

    # Compute the Earth Rotation Angle.  Time argument is UT1.

    theta = earth_rotation_angle(t.ut1)

    # The equinox method.  See Circular 179, Section 2.6.2.
    # Precession-in-RA terms in mean sidereal time taken from third
    # reference, eq. (42), with coefficients in arcseconds.

    t = (t.tdb - T0) / 36525.0
    st =        ( 0.014506 +
        (((( -    0.0000000368   * t
             -    0.000029956  ) * t
             -    0.00000044   ) * t
             +    1.3915817    ) * t
             + 4612.156534     ) * t)

    # Form the Greenwich sidereal time.

    return (st / 54000.0 + theta * 24.0) % 24.0
1239711
def rotate(self, azimuth, axis=None):
        
        target = self._target

        y_axis = self._n_pose[:3, 1].flatten()
        if axis is not None:
            y_axis = axis
        x_rot_mat = transformations.rotation_matrix(azimuth, y_axis, target)
        self._n_pose = x_rot_mat.dot(self._n_pose)

        y_axis = self._pose[:3, 1].flatten()
        if axis is not None:
            y_axis = axis
        x_rot_mat = transformations.rotation_matrix(azimuth, y_axis, target)
        self._pose = x_rot_mat.dot(self._pose)
675228
def quaternion_rotation(quat, vector):
    
    dp = np.dot(quat[1:], vector)
    cos = (2*quat[0]*quat[0] - 1)
    return np.array([
        2 * (quat[0] * (quat[2] * vector[2] - quat[3] * vector[1]) + quat[1] * dp) + cos * vector[0],
        2 * (quat[0] * (quat[3] * vector[0] - quat[1] * vector[2]) + quat[2] * dp) + cos * vector[1],
        2 * (quat[0] * (quat[1] * vector[1] - quat[2] * vector[0]) + quat[3] * dp) + cos * vector[2]
    ], float)
1074996
def get_rotations(self):
        """"""
        if self.centrosymmetric:
            return np.vstack((self.rotations, -self.rotations))
        else:
            return self.rotations
1108538
float[] dtRotate2D(float[] v, float ang) {
        float[] dest = new float[3];
        float c = (float) Math.cos(ang);
        float s = (float) Math.sin(ang);
        dest[0] = v[0] * c - v[2] * s;
        dest[2] = v[0] * s + v[2] * c;
        dest[1] = v[1];
        return dest;
    }
33747
def rotated(self, angle_degrees_ccw):
        """"""
        angle = angle_degrees_ccw / 180.0 * pi
        c, s = cos(angle), sin(angle)
        return self @ PdfMatrix((c, s, -s, c, 0, 0))
911154
def rotation_df(ATT):
    ''''''
    r = Matrix3()
    r.from_euler(radians(ATT.Roll), radians(ATT.Pitch), radians(ATT.Yaw))
    return r
999681
def pt_rotate(pt=(0.0, 0.0), angle=[0.0], center=(0.0, 0.0)):
    '''
    '''
    assert isinstance(pt, tuple)
    l_pt = len(pt)
    assert l_pt > 1
    for i in pt:
        assert isinstance(i, float)
    assert isinstance(angle, list)
    l_angle = len(angle)
    assert l_angle == l_pt-1
    for i in angle:
        assert isinstance(i, float)
        assert abs(i) <= 2*pi
    assert isinstance(center, tuple)
    assert len(center) == l_pt
    for i in center:
        assert isinstance(i, float)

    # Get vector from center to point and use to get relative polar coordinate.
    v_cart = [pt[i] - center[i] for i in range(l_pt)]

    # Length of vector needs to stay constant for new point.
    v_pol_l = [sqrt(v_cart[i]**2 + v_cart[i+1]**2) for i in range(l_angle)]
    v_pol_a = [(atan(v_cart[i+1] / v_cart[i]) if v_cart[i] != 0.0 else pi/2) + pi*int(pt[i] < center[i]) \
               for i in range(l_angle)]

    # Add rotation angle then convert back to cartesian vector.
    n_pol_a = [v_pol_a[i] + angle[i] for i in range(l_angle)]
    n_cart = [v_pol_l[0] * cos(n_pol_a[0])] + [v_pol_l[i] * sin(n_pol_a[i])\
                                               for i in range(l_angle)]

    # Add in the centre offset to get original offset from c.
    r = [n_cart[i] + center[i] for i in range(l_pt)]
    return tuple(r)
609634
def inverse_transform(t):
    
    return [
        quat_rotate(
            tf.transformations.quaternion_inverse(t[1]),
            [-t[0][0], -t[0][1], -t[0][2]]),
        tf.transformations.quaternion_inverse(t[1])
    ]