<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def sphere_points_from_angles_and_tilt(angles, tilted_axis):
 """ For a given tilt of the rotational axis `tilted_axis`, compute the points on a unit sphere that correspond to the distribution `angles` along the great circle about this axis. Parameters angles: 1d ndarray The angles that will be distributed on the great circle. tilted_axis: list of length 3 The tilted axis of rotation that determines the great circle. Notes ----- The reference axis is always [0,1,0]. `theta` is the azimuthal angle measured down from the y-axis. `phi` is the polar angle in the x-z plane measured from z towards x. """

    assert len(angles.shape) == 1
    # Normalize tilted axis.
    tilted_axis = norm_vec(tilted_axis)
    [u, v, w] = tilted_axis
    # Initial distribution of points about great circle (x-z).
    newang = np.zeros((angles.shape[0], 3), dtype=float)
    # We subtract angles[0], because in step (a) we want that
    # newang[0]==[0,0,1]. This only works if we actually start
    # at that point.
    newang[:, 0] = np.sin(angles-angles[0])
    newang[:, 2] = np.cos(angles-angles[0])

    # Compute rotational angles w.r.t. [0,1,0].
    # - Draw a unit sphere with the y-axis pointing up and the
    #   z-axis pointing right
    # - The rotation of `tilted_axis` can be described by two
    #   separate rotations. We will use these two angles:
    #   (a) Rotation from y=1 within the y-z plane: theta
    #       This is the rotation that is critical for data
    #       reconstruction. If this angle is zero, then we
    #       have a rotational axis in the imaging plane. If
    #       this angle is PI/2, then our sinogram consists
    #       of a rotating image and 3D reconstruction is
    #       impossible. This angle is counted from the y-axis
    #       onto the x-z plane.
    #   (b) Rotation in the x-z plane: phi
    #       This angle is responsible for matching up the angles
    #       with the correct sinogram images. If this angle is zero,
    #       then the projection of the rotational axis onto the
    #       x-y plane is aligned with the y-axis. If this angle is
    #       PI/2, then the axis and its projection onto the x-y
    #       plane are identical. This angle is counted from the
    #       positive z-axis towards the positive x-axis. By default,
    #       angles[0] is the point that touches the great circle
    #       that lies in the x-z plane. angles[1] is the next point
    #       towards the x-axis if phi==0.

    # (a) This angle is the azimuthal angle theta measured from the
    #     y-axis.
    theta = np.arccos(v)

    # (b) This is the polar angle measured in the x-z plane starting
    #     at the x-axis and measured towards the positive z-axis.
    if np.allclose(u, 0) and np.allclose(w, 0):
        # Avoid flipping the axis of rotation due to numerical
        # errors during its computation.
        phi = 0
    else:
        phi = np.arctan2(u, w)

    # Determine the projection points on the unit sphere.
    # The resulting circle meets the x-z-plane at phi, and
    # is tilted by theta w.r.t. the y-axis.

    # (a) Create a tilted data set. This is achieved in 3 steps.

    # a1) Determine radius of tilted circle and get the centered
    #     circle with a smaller radius.
    rtilt = np.cos(theta)
    newang *= rtilt

    # a2) Rotate this circle about the x-axis by theta
    #     (right-handed/counter-clockwise/basic/elemental rotation)
    Rx = np.array([
        [1,             0,              0],
        [0, np.cos(theta), -np.sin(theta)],
        [0, np.sin(theta),  np.cos(theta)]
    ])
    for ii in range(newang.shape[0]):
        newang[ii] = np.dot(Rx, newang[ii])

    # a3) Shift newang such that newang[0] is located at (0,0,1)
    newang = newang - (newang[0] - np.array([0, 0, 1])).reshape(1, 3)

    # (b) Rotate the entire thing with phi about the y-axis
    #     (right-handed/counter-clockwise/basic/elemental rotation)
    Ry = np.array([
        [+np.cos(phi), 0, np.sin(phi)],
        [0,            1,           0],
        [-np.sin(phi), 0, np.cos(phi)]
    ])

    for jj in range(newang.shape[0]):
        newang[jj] = np.dot(Ry, newang[jj])

    # For visualiztaion:
    # import matplotlib.pylab as plt
    # from mpl_toolkits.mplot3d import Axes3D
    # from matplotlib.patches import FancyArrowPatch
    # from mpl_toolkits.mplot3d import proj3d
    #
    # class Arrow3D(FancyArrowPatch):
    #     def __init__(self, xs, ys, zs, *args, **kwargs):
    #         FancyArrowPatch.__init__(self, (0,0), (0,0), *args, **kwargs)
    #         self._verts3d = xs, ys, zs
    #
    #     def draw(self, renderer):
    #         xs3d, ys3d, zs3d = self._verts3d
    #         xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, renderer.M)
    #         self.set_positions((xs[0],ys[0]),(xs[1],ys[1]))
    #         FancyArrowPatch.draw(self, renderer)
    #
    # fig = plt.figure(figsize=(10,10))
    # ax = fig.add_subplot(111, projection='3d')
    # for vec in newang:
    #     u,v,w = vec
    #     a = Arrow3D([0,u],[0,v],[0,w],
    #                 mutation_scale=20, lw=1, arrowstyle="-|>")
    #     ax.add_artist(a)
    #
    # radius=1
    # ax.set_xlabel('X')
    # ax.set_ylabel('Y')
    # ax.set_zlabel('Z')
    # ax.set_xlim(-radius*1.5, radius*1.5)
    # ax.set_ylim(-radius*1.5, radius*1.5)
    # ax.set_zlim(-radius*1.5, radius*1.5)
    # plt.tight_layout()
    # plt.show()

    return newang