// -*- Mode: c++ -*-
#ifndef	__VEC3_H__
#define	__VEC3_H__

// uses old-style iostream include; if using Standard C++ iostream (new style)
// in the same code, change this include to 
// #include <iostream>
// to avoid problems

#include <iostream.h>
#include <assert.h>
#include <math.h>

// Vectors of length 3
// This class, together with overloaded binary finctions,
// allow us to write vector expressions in the form close to 
// mathematical, e.g. 
//   Vec3 v,u,p,r; double s; 
//   r = s*cross(v,u) - dot(v,u)*p.
//
// The following operations are implemented (v,w are vectors, s is a scalar):
// binary:
// v+w; v-w; s*v; v*s; v/s; dot(v,w); cross(v,w); Min(v,w); Max(v,w)
// assignments and self-modifying:
//  v = w; v += w; v-= w; v*= s; v /= s; v.normalize() 
// unary:
//  -v; v.dot(); v.length(); v.dir(); v.l1(); 
// input/output:
// os << v; is >> v;



class Vec3 {
public:
  enum{ X = 0, Y = 1, Z = 2 };

  //default constructor: do nothing
  Vec3( void ) {}
  // copy constructor
  Vec3( const Vec3& w ) { v[X] = w.v[X], v[Y] = w.v[Y],v[Z] = w.v[Z]; }
  // vector from three numbers
  Vec3( double x, double y, double z) { v[X] = x, v[Y] = y; v[Z] = z; }

  // assignment
  Vec3& operator=( const Vec3& w )
    { v[X] = w.v[X]; v[Y] = w.v[Y];v[Z] = w.v[Z]; return *this; }

   // accessors to components
  double& x( void )       { return v[X]; }
  double  x( void ) const { return v[X]; }
  double& y( void )       { return v[Y]; }
  double  y( void ) const { return v[Y]; }
  double& z( void )       { return v[Z]; }
  double  z( void ) const { return v[Z]; }
  double& operator() (const int i)       
    { assert(i >= 0 && i < 3); return v[i]; }
  double  operator() (const int i) const 
    { assert(i >= 0 && i < 3); return v[i]; }
 
 
  // operations combined with assignment 
  Vec3& operator+=( const Vec3& w )
    { v[X] += w.v[X], v[Y] += w.v[Y]; v[Z] += w.v[Z]; return *this; }
  Vec3& operator-=( const Vec3& w )
    { v[X] -= w.v[X], v[Y] -= w.v[Y]; v[Z] -= w.v[Z]; return *this; }
  Vec3& operator*=( double s )
    { v[X] *= s, v[Y] *= s, v[Z] *= s; return *this; }
  Vec3& operator/=( double s )
    { v[X] /= s, v[Y] /= s, v[Z] /= s; return *this; }

  // unary minus 
  Vec3 operator-( void ) const
    { return Vec3( -v[X], -v[Y], -v[Z]); }
   
  // dot product of a vector with itself
  double dot( void ) const
    { return v[X] * v[X] + v[Y] * v[Y] + v[Z] * v[Z]; }

  double length( void ) const { return sqrt( dot() ); }

  //: l1 norm (sum of abs values of the components)
  double l1( void ) const
    { double a = fabs( v[X] ); a += fabs( v[Y] );a += fabs( v[Z] ); return a; }

  // for nonzero vectors, returns the normalized vector
  Vec3 dir( void ) const {
    double l = length();
    if( ( l  == 0) ) return *this;
    else return Vec3(v[X] / l, v[Y]/l, v[Z]/l);
  }
  // make the vector unit length
  // returns 1 if normalization was successful, 0 if the vector was zero vector
  int normalize( void )
    { double l = length(); return (l == 0) ? 0 : ( *this *= 1 / l, 1 ); }


 private:
  double v[3];
};


// binary operations; 

inline Vec3 operator+( const Vec3& v, const Vec3& w ) 
  { return Vec3( v.x() + w.x(), v.y() + w.y(), v.z() + w.z()); }

inline Vec3 operator-( const Vec3& v, const Vec3& w ) 
  { return Vec3( v.x() - w.x(), v.y() - w.y(), v.z() - w.z()); }

inline Vec3 operator*( double s, const Vec3& v ) 
  { return Vec3( s*v.x(), s*v.y(), s*v.z()); }

inline Vec3 operator*(  const Vec3& v, double s ) 
  { return Vec3( s*v.x(), s*v.y(), s*v.z()); }

inline Vec3 operator/(  const Vec3& v, double s ) 
{ assert(s != 0); return Vec3( v.x()/s, v.y()/s, v.z()/s); }

inline double dot( const Vec3& v, const Vec3& w ) { 
  return v.x() * w.x() + v.y() * w.y() + v.z() * w.z();}

inline Vec3 cross( const Vec3& v, const Vec3& w ) { 
  return Vec3(
         v.y() * w.z() - v.z() * w.y(),
         v.z() * w.x() - v.x() * w.z(),
         v.x() * w.y() - v.y() * w.x() ); 
}

  // returns the componentwise minimum; we use Min and Max instead
  // of min and max  to avoid problems with Windows


  // returns the componentwise minimum
inline  Vec3 Min( const Vec3& v, const Vec3& w ) { 
  return Vec3(
    (v.x() < w.x() )? v.x():w.x(),
    (v.y() < w.y() )? v.y():w.y(),
    (v.z() < w.z() )? v.z():w.z());
}

  // returns the componentwise maximum
 
inline  Vec3 Max( const Vec3& v, const Vec3& w ) { 
  return Vec3(
    (v.x() > w.x() )? v.x():w.x(),
    (v.y() > w.y() )? v.y():w.y(),
    (v.z() > w.z() )? v.z():w.z());
}

inline ostream& operator<<(ostream& os, const Vec3& v) 
{ os << v.x() << " " << v.y() << " " << v.z(); return os; }

inline istream& operator>>(istream& is, Vec3& v) { 
  is >> v.x(); 
  is >> v.y();  
  is >> v.z(); 
  return is; 
}
#endif	/* __VEC3_H__ */