code: def random_quaternion(rand=None):
      
      if rand is None:
          rand = np.random.rand(3)
      else:
          assert len(rand) == 3
      r1 = np.sqrt(1.0 - rand[0])
      r2 = np.sqrt(rand[0])
      pi2 = math.pi * 2.0
      t1 = pi2 * rand[1]
      t2 = pi2 * rand[2]
      return np.array([np.cos(t2) * r2, np.sin(t1) * r1,
                       np.cos(t1) * r1, np.sin(t2) * r2])
docs: Return uniform random unit quaternion.
  
      rand: array like or None
          Three independent random variables that are uniformly distributed
          between 0 and 1.
  
      >>> q = random_quaternion()
      >>> np.allclose(1, vector_norm(q))
      True
      >>> q = random_quaternion(np.random.random(3))
      >>> len(q.shape), q.shape[0]==4
      (1, True)