// -------------------------------------------------------------------------
// -----                        AnaEvtFilter source file                  -----
// -----          Created 13/04/12  by P. Russotto                    -----
// -------------------------------------------------------------------------

#include "AnaEvtFilter.h"
#include "FairLogger.h"
#include "TClonesArray.h"

#include <iostream>
#include <stdio.h>
#include <sstream>
#include <istream>
#include <ostream>
#include <TString.h>
#include <string>


#include "TMath.h"

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


AnaEvtFilter::AnaEvtFilter()
  : FairTask(),
    //fTLEventContOutName("TLEVENT"),
    fSaveOutputToTree(kFALSE),
    fAsyeosRootGlobalOutName("TROOTGLOBAL")
    //fSaveAsyeosRootGlobalToTree(kFALSE)
{
                   //init calibration
}


// ---------------------------------------------------------------------
AnaEvtFilter::AnaEvtFilter(const char* name, Int_t iVerbose)
  : FairTask(name, iVerbose),
    //fRootCHIEventCopy(0),
    //fRootMBallEventCopy(0),
    //fRootATOFEventCopy(0),
    fRootCHIGlobalCopy(0),  
    fRootLANDEventCopy(0),
    fRootATOFGlobalCopy(0),
    fRootMBallGlobalCopy(0),
    fRootCHITSCopy(0), 
    fRootMBSTSCopy(0), 
    fINFOEVENTCopy(0),
    //fGoodEvent(kFALSE),
    fAsyeosRootGlobal(0),

    fSaveOutputToTree(kFALSE),
    fAsyeosRootGlobalOutName("TROOTGLOBAL")
    //fSaveAsyeosRootGlobalToTree(kFALSE)
{
  fVerbose = 1; 
}
// ---------------------------------------------------------------------
AnaEvtFilter::~AnaEvtFilter() 
{ 
}
// ---------------------------------------------------------------------
InitStatus AnaEvtFilter::Init()
{ 
  
 
  fLogger->Info(MESSAGE_ORIGIN," AnaEvtFilter::Init()------------------Start "); 
  FairRootManager* ioman = FairRootManager::Instance();
  if ( ! ioman ) Fatal("Init", "No FairRootManager");

  fRootCHIGlobalCopy   = (TClonesArray*) ioman->GetObject("CHIGLOBALCLONE");
  fRootLANDEventCopy   = (TClonesArray*) ioman->GetObject("LANDEVENTCLONE");
  fRootATOFGlobalCopy  = (TClonesArray*) ioman->GetObject("ATOFGLOBALCLONE");
  fRootMBallGlobalCopy = (TClonesArray*) ioman->GetObject("MBALLGLOBALCLONE");

  fRootCHITSCopy       = (TClonesArray*) ioman->GetObject("CHITSCLONE");
  fRootMBSTSCopy       = (TClonesArray*) ioman->GetObject("MBSTSCLONE");
  fINFOEVENTCopy       = (TClonesArray*) ioman->GetObject("INFOEVENT");

  CHIGlob      = (TRootCHIGlobal*)       fRootCHIGlobalCopy;     // fCHIGlobal
  CalEvt       = (TRootLANDEvent*)       fRootLANDEventCopy;     // fLANDEvent                           
  ATOFGlob     = (TRootATOFglobalEvent*) fRootATOFGlobalCopy;    // fATOFGlobal
  MBallGlob    = (TRootMBallGlobal*)     fRootMBallGlobalCopy;   // fMBallGlobal

  CHITS        = (TRootTS*)              fRootCHITSCopy;         // fCHITS
  MBSTS        = (TRootTS*)              fRootMBSTSCopy;         // fMBSTS 
  INFOEVENT    = (TRootINFOEvent*)       fINFOEVENTCopy;         // fINFOEVENT

  //fLEvent= new TLEvent();      
  //ioman->Register(fTLEventContOutName, "TL event", fLEvent, fSaveOutputToTree);


  fAsyeosRootGlobal= new TRootGlobal();      
  ioman->Register(fAsyeosRootGlobalOutName, "TRootGlobal event", fAsyeosRootGlobal, fSaveOutputToTree);

  ///SK int imax=(int)(1*h1->GetEntries()/ratio);
  ///SK cout << imax << endl;
                   
  fGoodEvent          = false;
  fBWin0              = false;
  fBWin1              = false;
  fBWin2              = false;
  fBWin3              = false;
  
  fZBound             =-1;
  fZrat_x10           =-1;
  
  fPhirpATOFCHI       = -1000.;
  fDphirpATOFCHI12    = -1000.;
  fPhirpCHI           = -1000.;
  fDphirpCHI12        = -1000.;
  fPhirpMB            = -1000.;
  fDrpATOFCHI_MB      = -1000.;
  fDphiCHIMB          = -1000.;
  fDphiATOFCHI        = -1000.;
  
  fMult               = -10;
  fMultChiAtof        = -10;
  
  /// Counters:
  fEventCounter       = 0; 
  fGoodEventCounter   = 0;
  
  fGoodCHICounter     = 0;
  fGoodATOFCounter    = 0;
  fGoodRPCHIMBCounter = 0;
  fGoodTrigCounter    = 0;
  fGoodMBCounter      = 0;
  
  fBWin0Counter       = 0;
  fBWin1Counter       = 0;
  fBWin2Counter       = 0;
  fBWin3Counter       = 0;
  
}
// ---------------------------------------------------------------------
void AnaEvtFilter::SetParTask()
{
}
// ---------------------------------------------------------------------
void AnaEvtFilter::Finish()
{ 
  if (fVerbose > 1){
    float noEvt = float (fEventCounter);
    float goodEvt = float(fGoodEventCounter);
    cout << "[AnaEvtFilter::Finish:] goodEvt factor =";
    cout << 100.*float(fGoodEventCounter)/noEvt << " %";
    cout << " (" << fGoodEventCounter << " of " << fEventCounter <<")";
    cout << endl;
    
    cout << " - goodCHI:    \t" << 100.* fGoodCHICounter     /noEvt << " [%] (" << fGoodCHICounter     << " of " << fEventCounter << ")" << endl;
    cout << " - goodATOF:   \t" << 100.* fGoodATOFCounter    /noEvt << " [%] (" << fGoodATOFCounter    << " of " << fEventCounter << ")" << endl;
    cout << " - goodRPCHIMB:\t" << 100.* fGoodRPCHIMBCounter /noEvt << " [%] (" << fGoodRPCHIMBCounter << " of " << fEventCounter << ")" << endl;
    cout << " - goodTrig:   \t" << 100.* fGoodTrigCounter    /noEvt << " [%] (" << fGoodTrigCounter    << " of " << fEventCounter << ")" << endl;
    //cout << " - goodMB:     \t" << 100.* fGoodMBCounter      /noEvt << " [%] (" << fGoodMBCounter      << " of " << fEventCounter << ")" << endl;
    
    cout << "[AnaEvtFilter::Finish:] b windows: "<< endl;
    cout << " - bWin=0 :\t" << 100.* fBWin0Counter /noEvt << " [%] (" << fBWin0Counter << " of " << fEventCounter << ")" << endl;
    cout << " - bWin=1 :\t" << 100.* fBWin1Counter /noEvt << " [%] (" << fBWin1Counter << " of " << fEventCounter<< ")" << endl;
    cout << " - bWin=2 :\t" << 100.* fBWin2Counter /noEvt << " [%] (" << fBWin2Counter << " of " << fEventCounter<< ")" << endl;
    cout << " - bWin=3 :\t" << 100.* fBWin3Counter /noEvt << " [%] (" << fBWin3Counter << " of " << fEventCounter << ")" << endl;
    cout << endl;
    
  }
}
// ---------------------------------------------------------------------
void AnaEvtFilter::Reset()
{
}
// ---------------------------------------------------------------------
void AnaEvtFilter::Exec(Option_t* opt)
{
  fEventCounter++;
  fAsyeosRootGlobal->Clean();
  
  if (fVerbose > 5){
    cout << "Before: " << *fAsyeosRootGlobal;
  }
  
  /// Event selection:
  IsGoodEvent();
  
  if ( fGoodEvent )
  {
    fGoodEventCounter++;
    fAsyeosRootGlobal->good_evt      = fGoodEvent;
    fAsyeosRootGlobal->run           = INFOEVENT->run;
    fAsyeosRootGlobal->evt           = INFOEVENT->ievt;
    
    /// [1/2] Impact parameter window and multiplicity:
    GetImpactParamWindow();
    
    fAsyeosRootGlobal->b_win0        = fBWin0;
    fAsyeosRootGlobal->b_win1        = fBWin1;
    fAsyeosRootGlobal->b_win2        = fBWin2;
    fAsyeosRootGlobal->b_win3        = fBWin3;
    fAsyeosRootGlobal->mult          = fMult;        
    fAsyeosRootGlobal->mult_chi_atof = fMultChiAtof;       
    
    /// [2/2] Reaction plae determination
    GetReacionPlane();
    
    fAsyeosRootGlobal->Zbound            = fZBound;
    fAsyeosRootGlobal->Zrat_x10          = fZrat_x10;
    fAsyeosRootGlobal->rp                = fPhirpATOFCHI;
    fAsyeosRootGlobal->drp_subevts       = fDphirpATOFCHI12;
    fAsyeosRootGlobal->rp_chi            = fPhirpCHI;
    fAsyeosRootGlobal->drp_chi_subevts   = fDphirpCHI12;    
    fAsyeosRootGlobal->rp_mball          = fPhirpMB;
    fAsyeosRootGlobal->drp_atofchi_mball = fDrpATOFCHI_MB;  
    fAsyeosRootGlobal->drp_chi_mball     = fDphiCHIMB;      
    fAsyeosRootGlobal->drp_atof_chi      = fDphiATOFCHI;    
    
    if (fVerbose > 5){
      cout << " OK" << endl;
      cout << "After:  ";
      cout << *fAsyeosRootGlobal << endl;
    }
  }else
  {    
    if (fVerbose > 5){
      cout << " BAD" << endl;
    }
  }
}

// ---------------------------------------------------------------------
Bool_t AnaEvtFilter::IsGoodEvent()
{
  
  /// Flag declaration and initialization:
  fGoodEvent       = false;
  
  bool goodRun      = true;
  bool goodCHI      = false;
  bool goodMB       = false;
  bool goodTrig     = false;
  bool goodRPCHIMB  = false;
  bool goodATOF     = false;

  /// Variables declaration and reading from input data structures:
  int   TPatt         = MBSTS   -> TPatt; 
  
  double phirpCHI     = CHIGlob -> phirp;
  double dphirpCHI12  = CHIGlob -> dphi;
  double Etrans12CHI  = CHIGlob -> Etrans12;
  double ET12CHI      = CHIGlob -> Etrans12;
  double tavecsiCHI   = CHIGlob -> tavecsi;
  double multiCHI     = CHIGlob -> multi;
  
  double multiMB      = MBallGlob->multi;
  double phirpMB      = MBallGlob->phirp;
  
  int    multiATOF    = ATOFGlob->npart;
  double ZMax1ATOF    = ATOFGlob->Zmax1;
  double ATOFZBound   = (double) ATOFGlob->Zbound2;
  double phirpATOF    = (double) ATOFGlob->phirp;
  
  double dphiCHIMB = -1000;
  
  /// Get multiplicity:
  fMult            = multiMB + multiCHI + multiATOF;
  fMultChiAtof     =           multiCHI + multiATOF;
  
  /// FILTERING:
  /// FILTER 1.: GOOD TRIGGER
  if(TPatt==512 || TPatt==1536 ){
     goodTrig = true;
     fGoodTrigCounter++;
   }
  
  // For LAND only TODO: Move that to the LAND Ana
  ///// FILTER 2.: GOOD RUN (TODO: Is this filter needed?)
  //int  curr_run = INFOEVENT->run;
  //if( curr_run >= 1799 && curr_run <= 1808)   goodRun = false;
  //if( curr_run >= 1769 && curr_run <= 1778)   goodRun = false;

  /// FILTER 3.: GOOD CHI
  if( multiCHI>0 && tavecsiCHI>1900 && phirpCHI>-200 && Etrans12CHI>0 && CHIGlob->ZL>-1000 && CHIGlob->ZT>-1000 )    
  {                                                                                                                 
    goodCHI=true;
    fGoodCHICounter++;
  }
  
  /// FILTER 4.: GOOD MB (MicroBall)
  if( multiMB>0  && phirpMB>-200 )
  {
    goodMB = true;
    fGoodMBCounter++;
  }

  /// FILTER 5.: GOOD ATOF
  if( ATOFZBound < 100  && ZMax1ATOF < 50 && ZMax1ATOF > 1 )
  {
    goodATOF = true;
    fGoodATOFCounter++;
  }
  /// FILTER 6.: GOOD dphiCHIMB
  if (goodCHI && goodMB){
    dphiCHIMB      = phirpCHI - phirpMB;
    if( dphiCHIMB > 180 ){
      dphiCHIMB = dphiCHIMB-360;
    }else if ( dphiCHIMB < -180 ){
      dphiCHIMB = dphiCHIMB + 360;
    }
    dphiCHIMB   = TMath::Abs(dphiCHIMB);
  } 
  if(dphiCHIMB>90)
  {
    goodRPCHIMB=true; 
    fGoodRPCHIMBCounter++;
  }
  
  /// Good event final decision:
  if ( fVerbose > 9)
  {
    cout << "[AnaEvtFilter::IsGoodEvent:]";
    cout << "\tgoodTrig="    << goodTrig;
    cout << "\tgoodRPCHIMB=" << goodRPCHIMB;
    cout << "\tgoodCHI="     << goodCHI;
    cout << "\tgoodATOF="    << goodATOF;
    cout << endl;
  }
  if(    goodCHI     && goodATOF
      && goodRPCHIMB && goodTrig
  )
  {
    fGoodEvent = true;
  }
  return fGoodEvent;
}
  
  
// ---------------------------------------------------------------------
Double_t AnaEvtFilter::GetReacionPlane()
{
  
  fPhirpCHI           = -1000.;
  fDphiCHIMB          = -1000.;
  fDphirpCHI12        = -1000.;

  fPhirpATOFCHI       = -1000.;
  fDphirpATOFCHI12    = -1000.;
  fDphiATOFCHI        = -1000.;
  fPhirpMB            = -1000.;
  fDrpATOFCHI_MB      = -1000.;

  /// Variable declaration: 
  double dphirpATOF12 = -1000.; //TODO: Check if this var. needs to be calculated in this function?

  /// Read parameter values from data structures:
  fPhirpCHI         = CHIGlob  -> phirp;
  fDphirpCHI12      = CHIGlob  -> dphi;
  fPhirpMB          = MBallGlob-> phirp;
  
  double QxCHI     = CHIGlob  -> Qx;
  double QyCHI     = CHIGlob  -> Qy;
  double Q1xCHI    = CHIGlob  -> Q1x;
  double Q1yCHI    = CHIGlob  -> Q2y;
  double Q2xCHI    = CHIGlob  -> Q2x;
  double Q2yCHI    = CHIGlob  -> Q2y;
  
  double QxATOF    = ATOFGlob -> Qx;
  double QyATOF    = ATOFGlob -> Qy;
  double Q1xATOF   = ATOFGlob -> Q1x;
  double Q1yATOF   = ATOFGlob -> Q1y;
  double Q2xATOF   = ATOFGlob -> Q2x;
  double Q2yATOF   = ATOFGlob -> Q2y;

  double Qx        =  QxATOF +  QxCHI;
  double Qy        =  QyATOF +  QyCHI;
  double Q1x       = Q1xATOF + Q1xCHI;
  double Q1y       = Q1yATOF + Q1yCHI;
  double Q2x       = Q2xATOF + Q2xCHI;
  double Q2y       = Q2yATOF + Q2yCHI;
  
  double phirpATOF = (double) ATOFGlob->phirp;
  
  /// CALCULATION:
   
  //if (fGoodCHI && fGoodMB){
    fDphiCHIMB      = fPhirpCHI - fPhirpMB;
    if( fDphiCHIMB > 180 )
    {
      fDphiCHIMB = fDphiCHIMB-360;
    }else if ( fDphiCHIMB < -180 ){
      fDphiCHIMB = fDphiCHIMB + 360;
    }
    fDphiCHIMB = TMath::Abs(fDphiCHIMB);
  //}
    
  //if( fGoodCHI && fGoodATOF && phirpATOF > -200)
  if (phirpATOF > -200)
  {
    fPhirpATOFCHI          = atan2(  Qy,  Qx)*180/3.14159;
    double fPhirpATOFCHI1  = atan2( Q1y, Q1x)*180/3.14159;
    double fPhirpATOFCHI2  = atan2( Q2y, Q2x)*180/3.14159;
    
    fDphirpATOFCHI12       = fPhirpATOFCHI1 - fPhirpATOFCHI2;
    
    if( fDphirpATOFCHI12 >  180) fDphirpATOFCHI12 = 360 -fDphirpATOFCHI12;
    if( fDphirpATOFCHI12 < -180) fDphirpATOFCHI12 = 360 +fDphirpATOFCHI12;
    if( fDphirpATOFCHI12 <    0) fDphirpATOFCHI12 =     -fDphirpATOFCHI12;

    double phirpATOF1      = atan2( Q1yATOF, Q1xATOF)*180/3.14159;
    double phirpATOF2      = atan2( Q2yATOF, Q2xATOF)*180/3.14159;
    dphirpATOF12           = phirpATOF1 - phirpATOF2;
    if( dphirpATOF12 >  180)  dphirpATOF12 = 360 -dphirpATOF12;
    if( dphirpATOF12 < -180)  dphirpATOF12 = 360 +dphirpATOF12;
    if( dphirpATOF12 <    0)  dphirpATOF12 =     -dphirpATOF12;

    /// Difference of the RP orientation between CHI and ATOF
    fDphiATOFCHI = fPhirpCHI - phirpATOF;
    if( fDphiATOFCHI > 180 ){
      fDphiATOFCHI = fDphiATOFCHI - 360;
    }else if ( fDphiATOFCHI <-180 ){
      fDphiATOFCHI = fDphiATOFCHI + 360;
    }
    fDphiATOFCHI = TMath::Abs( fDphiATOFCHI );
    
    /// Difference of the RP orientation between ATOF&CHI and MB
    fDrpATOFCHI_MB = fPhirpATOFCHI - fPhirpMB;
    if( fDrpATOFCHI_MB > 180 ){
      fDrpATOFCHI_MB = fDrpATOFCHI_MB - 360;
    }else if ( fDrpATOFCHI_MB <-180 ){
      fDrpATOFCHI_MB = fDrpATOFCHI_MB + 360;
    }
    fDrpATOFCHI_MB = TMath::Abs( fDrpATOFCHI_MB );
  }else{
    /// If the ATOF rp is not determined use rp given by CHIMERA. Code from: FE70dev47lite.C:1234
    fPhirpATOFCHI    = fPhirpCHI;   
    fDphirpATOFCHI12 = fDphirpCHI12;
  }
  
  return fPhirpATOFCHI;
}

// ---------------------------------------------------------------------
Int_t AnaEvtFilter::GetImpactParamWindow()
{
  Int_t b_window_bit = 0;  
  
  fBWin0 = false; 
  fBWin1 = false; 
  fBWin2 = false; 
  fBWin3 = false; 
  
  double CHIZL      = CHIGlob->ZL;
  double CHIZT      = CHIGlob->ZT;
  double CHIZBound  = (double) CHIGlob->Zbound;

  double ATOFZL     = (double) ATOFGlob->ZL;
  double ATOFZT     = (double) ATOFGlob->ZT;
  double ATOFZBound = (double) ATOFGlob->Zbound2;
  
  fZBound   = ATOFZBound + CHIZBound;
  double ZT = ATOFZT     + CHIZT;
  double ZL = ATOFZL     + CHIZL;
  
  if( ZT > 0  && ZL != 0) /// ZL should never be equal 0.
  {
    fZrat_x10 = 10.*ZT/ZL; 
  }
  
  //best==0 - better b estimation from root files produced on 26/11/2013 b=3-7 fm  
  if(  16.   <    fZBound  &&   fZBound   < 41.
    &&  0.29 <  fZrat_x10  &&  fZrat_x10  <  0.63)
  {
    fBWin0 = true; 
    fBWin0Counter++; 
    b_window_bit+=1;
    //cout << "B0";
  }

  //best==1 - better b estimation from root files produced on 26/11/2013 b=0-7.5 fm
  if(  0.    <   fZBound   &&   fZBound   < 44.5
    && 0.245 <  fZrat_x10  &&  fZrat_x10  <  2)
  {
    fBWin1 = true;
    fBWin1Counter++;
    b_window_bit+=2;
    //cout << "B1";
  }
  
  //best==2 - better b estimation from root files produced on 26/11/2013 b=0-3 fm
  if(  0   <  fZBound    &&   fZBound   < 18
    && 0.6 <  fZrat_x10  &&  fZrat_x10  <  2)
  {
    fBWin2 = true;
    fBWin2Counter++;
    b_window_bit+=4;
    //cout << "B2";
  }
  
  //best==3 - comp fopi c2 to do
  if( 18.5   < fZBound     &&  fZBound    < 34
    && 0.365 <  fZrat_x10  &&  fZrat_x10  <  0.585)
  { 
    fBWin3 = true;
    fBWin3Counter++;
    b_window_bit+=8;
    //cout << "B3";
  }
  return b_window_bit;
}



ClassImp(AnaEvtFilter)