code: def get_PD_direction(X1_prime, X2_prime, X3_prime, PD):
      
      n = len(X1_prime) - 1
      X1 = X1_prime[0] - X1_prime[n]
      X2 = X2_prime[0] - X2_prime[n]
      X3 = X3_prime[0] - X3_prime[n]
      R= numpy.array([X1, X2, X3])
      #print 'R (reference vector for PD direction)', R
      dot = numpy.dot(PD, R) # dot product of reference vector and the principal axis of the V matrix
      #print 'dot (dot of PD and R)', dot
      if dot < -1:
          dot = -1
      elif dot > 1:
          dot = 1
      if numpy.arccos(dot) > old_div(numpy.pi, 2.):
          #print 'numpy.arccos(dot) {} > numpy.pi / 2. {}'.format(numpy.arccos(dot), numpy.pi / 2)
          #print 'correcting PD direction'
          PD = -1. * numpy.array(PD)
      #print 'PD after get PD direction', PD
      return PD
docs: takes arrays of X1_prime, X2_prime, X3_prime, and the PD.  
      checks that the PD vector direction is correct