// -------------------------------------------------------------------------
// -----                TpcAsciiGiBUUGenerator source file             -----
// -----          Created 03/04/12  by S. Doerheim                     -----
// -------------------------------------------------------------------------
#include "TpcAsciiGiBUUGenerator.h"

#include "FairPrimaryGenerator.h"

#include "TDatabasePDG.h"
#include "TRandom3.h"
#include "TMath.h"
#include <iostream>
#include <string>

using std::cout;
using std::endl;

// -----   Default constructor   ------------------------------------------
TpcAsciiGiBUUGenerator::TpcAsciiGiBUUGenerator()
  :FairGenerator(),
   fInputFile(NULL),
   filter_on_PDG(false),
   fFileName(""),
   pdgToFilter(0),
   evNb(0),
   evNbFile(0),
   targetThickness(5),
   vetoR(2),
   fzoff(0)
{
}
// ------------------------------------------------------------------------


// -----   Standard constructor   -----------------------------------------
TpcAsciiGiBUUGenerator::TpcAsciiGiBUUGenerator(const char* fileName)
  :FairGenerator(),
   fInputFile(0),
   fFileName(fileName),
   filter_on_PDG(false),
   pdgToFilter(0),
   targetThickness(5),
   vetoR(2),
   fzoff(0)
{
  cout << "-I TpcAsciiGiBUUGenerator: Opening input file " << fileName << endl;
  fInputFile = new ifstream(fFileName);
  if ( ! fInputFile->is_open() ) {
    Fatal("TpcAsciiGiBUUGenerator","Cannot open input file.");
  }
  evNb=0;
  evNbFile=0;
  fPDG=TDatabasePDG::Instance();
}
// ------------------------------------------------------------------------


// -----   Destructor   ---------------------------------------------------
TpcAsciiGiBUUGenerator::~TpcAsciiGiBUUGenerator()
{
  CloseInput();
}
// ------------------------------------------------------------------------

// -----   Public method ReadEvent   --------------------------------------
Bool_t TpcAsciiGiBUUGenerator::ReadEvent(FairPrimaryGenerator* primGen)
{
  std::vector<TVector3> mom;
  std::vector<int> pdgs;
  std::vector<double> energies;

  Bool_t foundParticle=false;
  //std::cout<<"Searching for pdgID:"<<pdgToFilter<<"\n";
  while(true){
    if(ReadFileEvent(mom,pdgs,energies,foundParticle)){
      //std::cout<<".";
      if(filter_on_PDG&&!foundParticle){
        evNbFile++;
        continue;
      }else{
        double vx=0,vy=0,vz=0;
        GenerateVertex(vx,vy,vz);
        evNbFile++;
        evNb++;
       // cout << "-I TpcAsciiGiBUUGenerator: Event " << evNb <<", evInFile "<<evNbFile
       //        <<" multiplicity " << mom.size() << "\n";
        for(Int_t itrack=0; itrack<mom.size(); itrack++) {

          primGen->AddTrack(pdgs[itrack], mom[itrack].X(), mom[itrack].Y(),mom[itrack].Z(), vx, vy, vz,-1,true,energies[itrack]);  
        }
        return kTRUE;
      }
    }else{
      return kFALSE;
    }   
  }
  return kFALSE;
}
// ------------------------------------------------------------------------

Bool_t TpcAsciiGiBUUGenerator::ReadFileEvent(std::vector<TVector3> &mom, 
                     std::vector<int> &pdgs,
                     std::vector<double> &energies,
                     Bool_t &found){

  mom.clear();
  pdgs.clear();
  energies.clear();
  // Check for input file
  if ( ! fInputFile->is_open() ) {
    cout << "-E TpcAsciiGiBUUGenerator: Input file not open!" << endl;
    return kFALSE;
  }

    // If end of input file is reached : close it and abort run
  if ( fInputFile->eof() ) {
    cout << "-I TpcAsciiGiBUUGenerator: End of input file reached " << endl;
    CloseInput();
    return kFALSE;
  }

  // Define event variable to be read from file
  Int_t ntracks = 0, dummyInt = 0;
  Double_t dummyDouble = 0., weight=0.;

  // Define track variables to be read from file
  Int_t pdgID = 0;
  Double_t px = 0., py = 0., pz = 0., E=0., m=0.;
  
   // Read event header line from input file
  *fInputFile >> ntracks >> dummyInt >> weight>> dummyDouble >>  dummyDouble >>  dummyDouble;
  
  // Loop over tracks in the current event
  for (Int_t itrack=0; itrack<ntracks; itrack++) {
    // Read PID and momentum from file
    *fInputFile >> pdgID >> px >> py >> pz >> E >> m;
    if(pdgID==60221){
      continue;
    }
    double pdgMass= fPDG->GetParticle(pdgID)->Mass();
    if(pdgID==pdgToFilter){
      found=kTRUE;
    }
    TVector3 p(px,py,pz);
    E=TMath::Sqrt(m*m+p.Mag2());
    mom.push_back(p);
    pdgs.push_back(pdgID);
    energies.push_back(E);
  }
  return kTRUE;
}

// -----   Private method CloseInput   ------------------------------------
void TpcAsciiGiBUUGenerator::CloseInput()
{
  if ( fInputFile ) {
    if ( fInputFile->is_open() ) {
      cout << "-I TpcAsciiGiBUUGenerator: Closing input file "
           << fFileName << endl;
      fInputFile->close();
    }
    delete fInputFile;
    fInputFile = NULL;
  }
}
// ------------------------------------------------------------------------

void TpcAsciiGiBUUGenerator::GenerateVertex(double& vx,double& vy,double& vz){
  //TODO configurable target
  double v_phi=gRandom->Uniform(0,2*TMath::Pi());
  double v_r=gRandom->Uniform(0,vetoR);
  TVector2 a;
  a.SetMagPhi(v_r,v_phi);
  vx=a.X();  vy=a.Y();
  vz=gRandom->Uniform(-targetThickness/2,targetThickness/2)+fzoff;
}

bool TpcAsciiGiBUUGenerator::SetStartEvent(unsigned int eventNbr){
  std::cout<<"********************************** TpcAsciiGiBUUGenerator::SetStartEvent( "<<eventNbr<<") **********************************************************"<<std::endl;
  if(!filter_on_PDG){
    std::cout<<"no filtering"<<std::endl;
    while(evNbFile<eventNbr){
      if(!fInputFile->eof()){
        // Define event variable to be read from file
        Int_t ntracks = 0, dummyInt = 0;
        Double_t dummyDouble = 0.;

        // Read event header line from input file
        *fInputFile >> ntracks >> dummyInt >> dummyDouble>> dummyDouble >>  dummyDouble >>  dummyDouble;
        std::string buffer="";
        if(ntracks>200){
          std::cerr<<"number of tracks:"<<ntracks<<";"<<evNbFile<<std::endl;
          throw;
        }
        for(int i=0;i<ntracks;++i){
          std::getline(*fInputFile,buffer);
        }
        std::getline(*fInputFile,buffer);
        ++evNb;
        ++evNbFile;
      }else{
        return false;
      }
    }
  }else{
    std::cout<<"filtering for events with pdgID= "<<pdgToFilter<<std::endl;

    while(evNb<eventNbr){
      std::vector<TVector3> mom;
      std::vector<int> pdgs;
      std::vector<double> energies;
      Bool_t foundParticle=false;
      if(ReadFileEvent(mom,pdgs,energies,foundParticle)){
        if(!foundParticle){
          ++evNbFile;
        //  cout<<".";
          continue;
        }else{
          ++evNb;
          ++evNbFile;
          //cout<<";";
        }
      }
    }
  }

  return true;
}

ClassImp(TpcAsciiGiBUUGenerator)