def rotate_bvecs(in_bvec, in_matrix):
    """
    Rotates the input bvec file accordingly with a list of matrices.

    .. note:: the input affine matrix transforms points in the destination
      image to their corresponding coordinates in the original image.
      Therefore, this matrix should be inverted first, as we want to know
      the target position of :math:`\\vec{r}`.

    """
    import os
    import numpy as np
    name, fext = os.path.splitext(os.path.basename(in_bvec))
    if fext == '.gz':
        name, _ = os.path.splitext(name)
    out_file = os.path.abspath('%s_rotated.bvec' % name)
    bvecs = np.loadtxt(in_bvec).T
    new_bvecs = []
    if len(bvecs) != len(in_matrix):
        raise RuntimeError(
            'Number of b-vectors (%d) and rotation matrices (%d) should match.'
             % (len(bvecs), len(in_matrix)))
    for bvec, mat in zip(bvecs, in_matrix):
        if np.all(bvec == 0.0):
            new_bvecs.append(bvec)
        else:
            invrot = np.linalg.inv(np.loadtxt(mat))[:3, :3]
            newbvec = invrot.dot(bvec)
            new_bvecs.append(newbvec / np.linalg.norm(newbvec))
    np.savetxt(out_file, np.array(new_bvecs).T, fmt=b'%0.15f')
    return out_file
