/**
  \file
  Operations involving two vectors

  \author Lukas Nellen
  \version $Id$
  \date 29 Jun 2003
*/

#ifndef _utl_OperationsVV_h_
#define _utl_OperationsVV_h_

#include <utl/Vector.h>
#include <cmath>

static const char CVSId_utl_OperationsVV[] =
  "$Id$";


namespace utl {

  inline
  Vector
  operator+(const Vector& l, const Vector& r)
  {
    r.TransformTo(l.fCoordinateSystem);
    return Vector(l.fVector + r.fVector, l.fCoordinateSystem);
  }


  inline
  Vector
  operator-(const Vector& l, const Vector& r)
  {
    r.TransformTo(l.fCoordinateSystem);
    return Vector(l.fVector - r.fVector, l.fCoordinateSystem);
  }


  inline
  Vector
  operator*(const double d, const Vector& v)
  {
    return Vector(d*v.fVector, v.fCoordinateSystem);
  }


  inline
  Vector
  operator*(const Vector& v, const double d)
  {
    return Vector(d*v.fVector, v.fCoordinateSystem);
  }


  inline
  Vector
  operator/(const Vector& v, const double d)
  {
    return Vector(v.fVector/d, v.fCoordinateSystem);
  }


  inline
  double
  operator*(const Vector& l, const Vector& r)
  {
    r.TransformTo(l.fCoordinateSystem);
    return l.fVector * r.fVector;
  }


  inline
  double
  CosAngle(const Vector& l, const Vector& r)
  {
    const double magnitudeA = l.GetMag();
    const double magnitudeB = r.GetMag();
    return (magnitudeA && magnitudeB) ?
      (l * r) / (magnitudeA * magnitudeB) : 1;
  }


  inline
  double
  Angle(const Vector& left, const Vector& right)
  {
    /*const double cosAngle = CosAngle(left, right);
    return (cosAngle >= 1) ? 0 : acos(cosAngle);*/
    // DV: this is more accurate for small and large angles
    const Vector a = 1/left.GetMag() * left;
    const Vector b = 1/right.GetMag() * right;
    const double d2 = (a - b).GetMag2();
    if (d2 <= 2)
      return 2 * std::asin(0.5 * std::sqrt(d2));
    else
      return 2 * std::acos(0.5 * (a + b).GetMag());
  }

}


#endif

// Configure (x)emacs for this file ...
// Local Variables:
// mode: c++
// compile-command: "make -C .. -k"
// End: