#ifndef __JEIGENVALUES2D__
#define __JEIGENVALUES2D__

#include <cmath>
#include <map>
#include <functional>

#include "TMatrixDSymEigen.h"

#include "JGeometry2D/JVector2D.hh"
#include "JGeometry2D/JGeometry2DToolkit.hh"


/**
 * \author mdejong
 */

namespace JGEOMETRY2D {}
namespace JPP { using namespace JGEOMETRY2D; }

namespace JGEOMETRY2D {

  /**
   * Eigen values in 2D.
   */
  struct JEigenValues2D :
    public std::map<double, JVector2D, std::greater<double> >
  {
  public:
    /**
     * Constructor.
     *
     * \param  __begin        begin of data
     * \param  __end          end   of data
     */
    template<class T>
    JEigenValues2D(T __begin,
		   T __end) 
    {
      const JCenter2D center(__begin, __end);
	
      // RMS matrix
	
      TMatrixDSym A(2);

      A = 0.0;
	
      for (T i = __begin; i != __end; ++i) {
	
	const double dx = center.getX() - i->getX();
	const double dy = center.getY() - i->getY();
	  
	A(0,0) += (dx * dx);
	A(0,1) += (dx * dy);
	
	A(1,1) += (dy * dy);
      }
      
      A(1,0) = A(0,1);

      const TMatrixDSymEigen result(A);

      const TVectorD& u = result.GetEigenValues();
      const TMatrixD& V = result.GetEigenVectors();

      for (Int_t i = 0; i != u.GetNoElements(); ++i) {
	(*this)[u[i]] = JVector2D(V(0,i),
				  V(1,i));
      }
    }
  };
}
 
#endif