code: def _rotate_to_arbitrary_vector(v,a,inv=False,_dontcutsmall=False):
      
      normv= v/numpy.tile(numpy.sqrt(numpy.sum(v**2.,axis=1)),(3,1)).T
      rotaxis= numpy.cross(normv,a)
      rotaxis/= numpy.tile(numpy.sqrt(numpy.sum(rotaxis**2.,axis=1)),(3,1)).T
      crossmatrix= numpy.empty((len(v),3,3))
      crossmatrix[:,0,:]= numpy.cross(rotaxis,[1,0,0])
      crossmatrix[:,1,:]= numpy.cross(rotaxis,[0,1,0])
      crossmatrix[:,2,:]= numpy.cross(rotaxis,[0,0,1])
      costheta= numpy.dot(normv,a)
      sintheta= numpy.sqrt(1.-costheta**2.)
      if inv: sgn= 1.
      else: sgn= -1.
      out= numpy.tile(costheta,(3,3,1)).T*numpy.tile(numpy.eye(3),(len(v),1,1))\
          +sgn*numpy.tile(sintheta,(3,3,1)).T*crossmatrix\
          +numpy.tile(1.-costheta,(3,3,1)).T\
          *(rotaxis[:,:,numpy.newaxis]*rotaxis[:,numpy.newaxis,:])
      if not _dontcutsmall:
          out[numpy.fabs(costheta-1.) < 10.**-10.]= numpy.eye(3)
          out[numpy.fabs(costheta+1.) < 10.**-10.]= -numpy.eye(3)
      return out
docs: Return a rotation matrix that rotates v to align with unit vector a
          i.e. R . v = |v|\hat{a}