//-----------------------------------------------------------------------
// File and Version Information:
//     PacEmcCluster;
//     Class for PacEmcCluster for fast simulation, a stripped down version
//     from Babar's EmcCluster.
//
// Author List:
//     Chih-hsiang Cheng       Caltech   (initial version 2008/06/13)
//
//-----------------------------------------------------------------------

#include "BaBar/BaBar.hh"

//-----------------------
// This Class's Header --
//-----------------------

#include "PacEmc/PacEmcCluster.hh"

//---------------
// C++ Headers --
//---------------

#include <iostream>
#include <iomanip>
#include <math.h>
#include <algorithm>
#include <float.h>

#include <vector>
using std::vector;

//-------------------------------
// Collaborating Class Headers --
//-------------------------------

#include "AbsEnv/AbsEnv.hh"
#include "AbsEnv/TwoCoordIndex.hh"
#include "BaBar/Constants.hh"
#include "BbrGeom/BbrError.hh"
//#include "BbrStdUtils/CollectionUtils.hh"
#include "CLHEP/Vector/LorentzVector.h"
#include "CLHEP/Matrix/Vector.h"
#include "ErrLogger/ErrLog.hh"
#include "G3Data/GTrack.hh"
#include "TrajGeom/TrkLineTraj.hh"

#include "EmcData/EmcClusterLogPos.hh"
#include "EmcData/EmcClusterLiloPos.hh"
#include "EmcData/EmcClusterMoments.hh"
#include "EmcData/EmcClusterEnergySums.hh"
#include "EmcData/EmcClusterDistances.hh"
#include "EmcData/EmcXClMoments.hh"

#include "EmcEnv/EmcEnv.hh"
#include "EmcEnv/EmcStructure.hh"

using std::endl;
using std::ios;
using std::ostream;
using std::setw;

//----------------
// Constructors --
//----------------

PacEmcCluster::PacEmcCluster() :  AbsRecoCalo()
  ,_pacEmcDigis( new std::vector< PacEmcDigi* > )
  ,_nbumps( 0 )
  ,_gtracks( new std::vector< GTrack > )
  ,_energyValid( false )
  ,_energy( 0 )
  ,_whereValid( false )
  ,_where( 0 )
//   ,_theClusEnergySums( 0 )
//   ,_theClusMoments( 0 )
//   ,_theClusXClMoments( 0 )
//   ,_theClusLogPos( 0 )
//   ,_theClusLiloPos( 0 )
//   ,_theClusDistances( 0 )
{
}

//--------------
// Destructor --
//--------------
PacEmcCluster::~PacEmcCluster()
{
  delete _pacEmcDigis;
  delete _gtracks;
  delete _where;

//   if ( _theClusEnergySums != 0) delete _theClusEnergySums;
//   if ( _theClusMoments != 0) delete _theClusMoments;
//   if ( _theClusXClMoments != 0) delete _theClusXClMoments;
//   if ( _theClusLogPos != 0) delete _theClusLogPos;
//   if ( _theClusLiloPos != 0) delete _theClusLiloPos;
//   if ( _theClusDistances != 0) delete _theClusDistances;

}

const std::vector< PacEmcDigi* >*
PacEmcCluster::pacEmcDigis() const {
  return _pacEmcDigis;
}

const std::vector< GTrack >*
PacEmcCluster::gtracks() const {
  return _gtracks;
}

double
PacEmcCluster::energy() const {
  
  if  ( ! _energyValid ) {
    double sum=0;
    
    PacEmcDigiiterator iter= _pacEmcDigis->begin();
    while ( (iter != _pacEmcDigis->end() ) ) {
      sum+= (*iter)->energy();
      iter++;
    }
    
    _energy = sum;
    _energyValid = true;
  }
  
  return _energy;
  
}

double
PacEmcCluster::rawEnergy() const {
  return this->energy();
}


double
PacEmcCluster::theta() const
{
  return where().theta();
}

double
PacEmcCluster::phi() const
{
  return where().phi();
}


double
PacEmcCluster::thetaIndex() const {
  return 0;
}

int
PacEmcCluster::thetaIndexInt() const {
  return 0;
}

double
PacEmcCluster::phiIndex() const {
  return 0;
}

int
PacEmcCluster::phiIndexInt() const {
  return 0;
}

double
PacEmcCluster::thetaFromIndex() const {
  return 0;
}

double
PacEmcCluster::phiFromIndex() const {
  return 0;
}


double
PacEmcCluster::major_axis() const {
  return 0;
}

const PacEmcDigi*
PacEmcCluster::maxima() const {
  double max=0;
  PacEmcDigi* biggest(0);
  PacEmcDigiiterator iter= _pacEmcDigis->begin();
  while ( (iter!= _pacEmcDigis->end() ) ) {
    if ( max < (*iter)->energy() ) {
      max= (*iter)->energy();
      biggest= (*iter);
    }
  }
  
  return biggest;
}

HepPoint
PacEmcCluster::position( const HepPoint &,
                      const PdtEntry* ) const {
  return this->where();
}

HepPoint
PacEmcCluster::rawPosition( const HepPoint &,
                         const PdtEntry* ) const {
  return this->position();
}

BbrError
PacEmcCluster::errorMatrix( const HepPoint &fromHere,
                         const PdtEntry* particle ) const {
  // Riccardo Faccini 8th April 2000:
  //forward the call toAbsReco in order to share it with BtaDataP (consistency of the micro)

  return errorMatrixOfEmcClusters(energy(),theta(),phi());
}

double
PacEmcCluster::mass() const {
  return this->fourMomentum().mag();
}

HepPoint
PacEmcCluster::where() const {
  if ( !_whereValid )
    {
      delete _where;
      //_where = new HepPoint( algWhere( this ) );
      _where= new HepPoint(gravWhere(this));
      _whereValid = true;
    }
  
  return *_where;
}

double
PacEmcCluster::x() const {
  return where().x();
}

double
PacEmcCluster::y() const {
  return where().y();
}

double
PacEmcCluster::z() const {
  return where().z();
}

//-------------
// Modifiers --
//-------------

void
PacEmcCluster::addDigi(PacEmcDigi* theDigi ) {

  _pacEmcDigis->push_back(theDigi);

}

void
PacEmcCluster::addCluster(PacEmcCluster* cluster) {
  // FIXME
  //_myRefs->push_back((AbsRecoCalo*)cluster);
}

double
PacEmcCluster::distanceToCentre( const HepPoint& aPoint ) const
{
  return  ( where() - aPoint ).mag();
}

double
PacEmcCluster::distanceToCentre( const PacEmcDigi* aDigi ) const
{
  return  ( where() - aDigi->where() ).mag();
}


int
PacEmcCluster::numberOfDigis() const {
  return _pacEmcDigis->size();
}

void
PacEmcCluster::sortDigis() {
  // FIXME
  // Sort digis in cluster with distance from centre of cluster

}

unsigned
PacEmcCluster::nBumps() const {
  return _nbumps;
}

void
PacEmcCluster::setNBumps(unsigned nbumps) {
  _nbumps = nbumps;
}



void
PacEmcCluster::invalidateCache(){
  _energyValid = false;
  _whereValid = false;
  
}


const PacEmcClusterEnergySums&
PacEmcCluster::esums() const
{
  if (_theClusEnergySums == 0)
    _theClusEnergySums = new PacEmcClusterEnergySums( *this );
  return *_theClusEnergySums;
}



//              -----------------------------------------------
//              -- Static Data & Function Member Definitions --
//              -----------------------------------------------

HepPoint
PacEmcCluster::gravWhere( const PacEmcCluster* me ) {
  Hep3Vector aVector(0,0,0);

  PacEmcDigiiterator iter= me->pacEmcDigis()->begin();

  while ( (iter!= me->pacEmcDigis()->end()) ) {
    HepPoint aPoint( (*iter)->where() );
    aVector += Hep3Vector( aPoint.x(), aPoint.y(), aPoint.z() )
      * (*iter)->energy();
  }

  aVector *= 1./me->energy();
  
  return HepPoint( aVector.x(), aVector.y(), aVector.z() );
  
}