macro and flowevent maker to run part of the code in root
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 30 Jan 2009 13:12:48 +0000 (13:12 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 30 Jan 2009 13:12:48 +0000 (13:12 +0000)
PWG2/FLOW/other/code/FlowEventSimpleMaker.cxx [new file with mode: 0644]
PWG2/FLOW/other/code/FlowEventSimpleMaker.h [new file with mode: 0644]
PWG2/FLOW/other/runFlowAnalysis.C [new file with mode: 0644]

diff --git a/PWG2/FLOW/other/code/FlowEventSimpleMaker.cxx b/PWG2/FLOW/other/code/FlowEventSimpleMaker.cxx
new file mode 100644 (file)
index 0000000..b1ba8d9
--- /dev/null
@@ -0,0 +1,150 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+/* $Id */
+
+#include "Riostream.h"
+#include "FlowEventSimpleMaker.h"
+#include "AliFlowEventSimple.h"
+#include "AliFlowTrackSimple.h"
+#include "TTree.h"
+#include "TParticle.h"
+#include "AliFlowTrackSimpleCuts.h"
+
+
+// FlowEventSimpleMaker:
+// Class to fill the AliFlowEventSimple
+// with AliFlowTrackSimple objects
+// ouside the AliRoot Framework
+// Has fill methods for TTree, 
+
+ClassImp(FlowEventSimpleMaker)
+//----------------------------------------------------------------------- 
+FlowEventSimpleMaker::FlowEventSimpleMaker()
+{
+  //constructor
+}
+
+//-----------------------------------------------------------------------   
+FlowEventSimpleMaker::~FlowEventSimpleMaker()
+{
+  //destructor
+}
+
+//-----------------------------------------------------------------------   
+AliFlowEventSimple* FlowEventSimpleMaker::FillTracks(TTree* anInput, AliFlowTrackSimpleCuts* intCuts, AliFlowTrackSimpleCuts* diffCuts)
+{
+  //fills the event from a TTree of kinematic.root files
+  
+  // number of times to use the same particle (trick to introduce nonflow)
+  Int_t iLoops = 1;
+  
+  //flags for particles passing int. and diff. flow cuts
+  Bool_t bPassedIntFlowCuts  = kFALSE;
+  Bool_t bPassedDiffFlowCuts = kFALSE;
+  
+  //track cut values
+  Double_t dPtMaxInt  = intCuts->GetPtMax();
+  Double_t dPtMinInt  = intCuts->GetPtMin();
+  Double_t dEtaMaxInt = intCuts->GetEtaMax();
+  Double_t dEtaMinInt = intCuts->GetEtaMin();
+  Double_t dPhiMaxInt = intCuts->GetPhiMax();
+  Double_t dPhiMinInt = intCuts->GetPhiMin();
+  Int_t iPIDInt       = intCuts->GetPID();
+  
+  Double_t dPtMaxDiff  = diffCuts->GetPtMax();
+  Double_t dPtMinDiff  = diffCuts->GetPtMin();
+  Double_t dEtaMaxDiff = diffCuts->GetEtaMax();
+  Double_t dEtaMinDiff = diffCuts->GetEtaMin();
+  Double_t dPhiMaxDiff = diffCuts->GetPhiMax();
+  Double_t dPhiMinDiff = diffCuts->GetPhiMin();
+  Int_t iPIDDiff       = diffCuts->GetPID();
+  
+  Int_t iNumberOfInputTracks = anInput->GetEntries() ;
+  //cerr<<"iNumberOfInputTracks = "<<iNumberOfInputTracks<<endl;
+  TParticle* pParticle = new TParticle();
+  anInput->SetBranchAddress("Particles",&pParticle);  
+  //  AliFlowEventSimple* pEvent = new AliFlowEventSimple(iNumberOfInputTracks);
+  AliFlowEventSimple* pEvent = new AliFlowEventSimple(10);
+  //cerr<<pEvent<<" pEvent "<<endl;
+  
+  Int_t iN = iNumberOfInputTracks; // additional variable to artificially fix the number of tracks
+  //  Int_t iN = 576; //multiplicity for chi=1.5
+  //  Int_t iN = 256; //multiplicity for chi=1
+  //  Int_t iN = 164; //multiplicity for chi=0.8
+  
+  Int_t iGoodTracks = 0;
+  Int_t itrkN = 0;
+  Int_t iSelParticlesDiff = 0;
+  Int_t iSelParticlesInt = 0;
+  
+  while (itrkN < iNumberOfInputTracks) {
+    anInput->GetEntry(itrkN);   //get input particle
+    //checking the cuts for int. and diff. flow
+    if (pParticle->Pt() > dPtMinInt && pParticle->Pt() < dPtMaxInt &&
+       pParticle->Eta() > dEtaMinInt && pParticle->Eta() < dEtaMaxInt &&
+       pParticle->Phi() > dPhiMinInt && pParticle->Phi() < dPhiMaxInt &&
+       TMath::Abs(pParticle->GetPdgCode()) == iPIDInt) { 
+      bPassedIntFlowCuts = kTRUE; 
+    } 
+    
+    if (pParticle->Pt() > dPtMinDiff && pParticle->Pt() < dPtMaxDiff &&
+       pParticle->Eta() > dEtaMinDiff && pParticle->Eta() < dEtaMaxDiff &&
+       pParticle->Phi() > dPhiMinDiff && pParticle->Phi() < dPhiMaxDiff &&
+       TMath::Abs(pParticle->GetPdgCode()) == iPIDDiff){ 
+      bPassedDiffFlowCuts = kTRUE; 
+    }
+    
+    if (bPassedIntFlowCuts || bPassedDiffFlowCuts) {
+      for(Int_t d=0;d<iLoops;d++) {
+       AliFlowTrackSimple* pTrack = new AliFlowTrackSimple();
+       pTrack->SetPt(pParticle->Pt());
+       pTrack->SetEta(pParticle->Eta());
+       pTrack->SetPhi(pParticle->Phi());
+       
+       //marking the particles used for int. flow:
+       if(bPassedIntFlowCuts && iSelParticlesInt < iN*iLoops) {  
+         pTrack->SetForIntegratedFlow(kTRUE);
+         iSelParticlesInt++;
+       }
+       //marking the particles used for diff. flow:
+       if(bPassedDiffFlowCuts) {
+         pTrack->SetForDifferentialFlow(kTRUE);
+         iSelParticlesDiff++;
+       }
+       //adding a particles which were used either for int. or diff. flow to the list
+       pEvent->TrackCollection()->Add(pTrack);
+       iGoodTracks++;
+      }//end of for(Int_t d=0;d<iLoops;d++)
+    }//end of if(bPassedIntFlowCuts || bPassedDiffFlowCuts) 
+    itrkN++;  
+    bPassedIntFlowCuts  = kFALSE;
+    bPassedDiffFlowCuts = kFALSE;
+  }//end of while (itrkN < iNumberOfInputTracks)
+  
+  pEvent->SetEventNSelTracksIntFlow(iSelParticlesInt);  
+  pEvent->SetNumberOfTracks(iGoodTracks);//tracks used either for int. or for diff. flow
+
+  cout<<" iGoodTracks = "<<iGoodTracks<<endl;
+  cout<<" # of selected tracks for int. flow  = "<<iSelParticlesInt<<endl;
+  cout<<" # of selected tracks for diff. flow = "<<iSelParticlesDiff<<endl;  
+
+  delete pParticle;
+  return pEvent;
+}
+
+
+
+
diff --git a/PWG2/FLOW/other/code/FlowEventSimpleMaker.h b/PWG2/FLOW/other/code/FlowEventSimpleMaker.h
new file mode 100644 (file)
index 0000000..177c088
--- /dev/null
@@ -0,0 +1,36 @@
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+#ifndef FlowEventSimpleMaker_H
+#define FlowEventSimpleMaker_H
+
+#include "AliFlowEventSimple.h"  
+#include "AliFlowTrackSimpleCuts.h"
+
+class TTree;
+
+// FlowEventSimpleMaker:
+// Class to fill the AliFlowEventSimple with AliFlowTrackSimple objects
+// This works also outside the AliRoot Framework
+          
+class FlowEventSimpleMaker {
+
+ public:
+
+  FlowEventSimpleMaker();             //constructor
+  virtual  ~FlowEventSimpleMaker();   //destructor
+  
+  //TTree
+  AliFlowEventSimple* FillTracks(TTree* anInput, AliFlowTrackSimpleCuts* intCuts, AliFlowTrackSimpleCuts* diffCuts);   //use own cut class
+    
+ private:
+  FlowEventSimpleMaker(const FlowEventSimpleMaker& anAnalysis);            //copy constructor
+  FlowEventSimpleMaker& operator=(const FlowEventSimpleMaker& anAnalysis); //assignment operator
+          
+  ClassDef(FlowEventSimpleMaker,0)    // macro for rootcint
+};
+      
+#endif
+
diff --git a/PWG2/FLOW/other/runFlowAnalysis.C b/PWG2/FLOW/other/runFlowAnalysis.C
new file mode 100644 (file)
index 0000000..5fd6848
--- /dev/null
@@ -0,0 +1,471 @@
+#include "TStopwatch.h"
+#include "TObjArray"
+#include "Riostream.h"
+
+//RUN SETTINGS
+//flow analysis method can be: (set to kTRUE or kFALSE)
+Bool_t SP    = kTRUE;
+Bool_t LYZ1  = kTRUE;
+Bool_t LYZ2  = kFALSE;  
+Bool_t LYZEP = kFALSE; 
+Bool_t GFC   = kTRUE;
+Bool_t QC    = kTRUE;
+Bool_t FQD   = kTRUE;
+Bool_t MCEP  = kFALSE; //does not work yet 24/12/08
+
+Int_t offset = 0;
+
+//CUT SETTINGS
+//integrated selection
+Double_t ptMaxInt  = 10.;
+Double_t ptMinInt  = 0.;
+Double_t etaMaxInt = 1.;
+Double_t etaMinInt = -1.;
+Double_t phiMaxInt = 7.5;
+Double_t phiMinInt = 0.;
+Int_t PIDInt       = 211;
+
+//differential selection
+Double_t ptMaxDiff  = 10.;
+Double_t ptMinDiff  = 0.;
+Double_t etaMaxDiff = 1.;
+Double_t etaMinDiff = -1.;
+Double_t phiMaxDiff = 7.5;
+Double_t phiMinDiff = 0.;
+Int_t PIDDiff       = 211;
+
+int runFlowAnalysis(Int_t aRuns = 100, const char* 
+                         //                      dir="/data/alice1/kolk/KineOnly3/")
+                         dir="/Users/snelling/alice_data/KineOnly3/")
+{
+  TStopwatch timer;
+  timer.Start();
+  
+  if (LYZ1 && LYZ2) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1."<<endl; exit(); }
+  
+  if (LYZ2 && LYZEP) {cout<<"WARNING: you cannot run LYZ2 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
+  
+  if (LYZ1 && LYZEP) {cout<<"WARNING: you cannot run LYZ1 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
+
+
+  cout<<endl;
+  cout<<" ---- BEGIN ANALYSIS ---- "<<endl;
+  cout<<endl;
+  
+  //  gSystem->AddIncludePath("-I$ALICE_ROOT/include");
+  gSystem->AddIncludePath("-I$ROOTSYS/include");
+    
+  // constants  
+  gROOT->LoadMacro("code/AliFlowCommonConstants.cxx+");
+  gROOT->LoadMacro("code/AliFlowLYZConstants.cxx+");
+  gROOT->LoadMacro("code/AliFlowCumuConstants.cxx+");
+
+  // flow event
+  gROOT->LoadMacro("code/AliFlowVector.cxx+"); 
+  gROOT->LoadMacro("code/AliFlowTrackSimple.cxx+");    
+  gROOT->LoadMacro("code/AliFlowEventSimple.cxx+");
+
+  // cuts
+  gROOT->LoadMacro("code/AliFlowTrackSimpleCuts.cxx+");    
+
+  // output histosgrams
+  gROOT->LoadMacro("code/AliFlowCommonHist.cxx+");
+  gROOT->LoadMacro("code/AliFlowCommonHistResults.cxx+");
+  gROOT->LoadMacro("code/AliFlowLYZHist1.cxx+");
+  gROOT->LoadMacro("code/AliFlowLYZHist2.cxx+");
+
+  // functions needed for various methods
+  gROOT->LoadMacro("code/AliCumulantsFunctions.cxx+");
+  gROOT->LoadMacro("code/AliQCumulantsFunctions.cxx+");
+  gROOT->LoadMacro("code/AliFittingFunctionsForQDistribution.cxx+");
+  gROOT->LoadMacro("code/AliFlowLYZEventPlane.cxx+");
+
+  // Flow Analysis code for various methods
+  gROOT->LoadMacro("code/AliFlowAnalysisWithMCEventPlane.cxx+"); 
+  gROOT->LoadMacro("code/AliFlowAnalysisWithScalarProduct.cxx+");
+  gROOT->LoadMacro("code/AliFlowAnalysisWithLYZEventPlane.cxx+");
+  gROOT->LoadMacro("code/AliFlowAnalysisWithLeeYangZeros.cxx+");
+  gROOT->LoadMacro("code/AliFlowAnalysisWithCumulants.cxx+");
+  gROOT->LoadMacro("code/AliFlowAnalysisWithQCumulants.cxx+"); 
+  gROOT->LoadMacro("code/AliFittingQDistribution.cxx+");
+
+  // Method to fill the FlowEvent
+  gROOT->LoadMacro("code/FlowEventSimpleMaker.cxx+");   
+
+  //load needed libraries
+  gSystem->Load("libTree.so");
+
+  cout << "finished loading macros!" << endl;  
+  //  gSystem->Load("libANALYSIS.so");
+  //  gSystem->Load("libPWG2flow.so");
+
+  //------------------------------------------------------------------------
+  //cuts
+  AliFlowTrackSimpleCuts* cutsInt = new AliFlowTrackSimpleCuts();
+  cutsInt->SetPtMax(ptMaxInt);
+  cutsInt->SetPtMin(ptMinInt);
+  cutsInt->SetEtaMax(etaMaxInt);
+  cutsInt->SetEtaMin(etaMinInt);
+  cutsInt->SetPhiMax(phiMaxInt);
+  cutsInt->SetPhiMin(phiMinInt);
+  cutsInt->SetPID(PIDInt);
+
+  AliFlowTrackSimpleCuts* cutsDiff = new AliFlowTrackSimpleCuts();
+  cutsDiff->SetPtMax(ptMaxDiff);
+  cutsDiff->SetPtMin(ptMinDiff);
+  cutsDiff->SetEtaMax(etaMaxDiff);
+  cutsDiff->SetEtaMin(etaMinDiff);
+  cutsDiff->SetPhiMax(phiMaxDiff);
+  cutsDiff->SetPhiMin(phiMinDiff);
+  cutsDiff->SetPID(PIDDiff);
+
+  //------------------------------------------------------------------------
+  //flow event
+  FlowEventSimpleMaker* fEventMaker = new FlowEventSimpleMaker(); 
+  
+  AliFlowAnalysisWithQCumulants    *qc    = NULL;
+  AliFlowAnalysisWithCumulants     *gfc   = NULL;
+  AliFittingQDistribution          *fqd   = NULL;
+  AliFlowAnalysisWithLeeYangZeros  *lyz1  = NULL;
+  AliFlowAnalysisWithLeeYangZeros  *lyz2  = NULL;
+  AliFlowAnalysisWithLYZEventPlane *lyzep = NULL;
+  AliFlowAnalysisWithScalarProduct *sp    = NULL;
+  AliFlowAnalysisWithMCEventPlane  *mcep  = NULL;
+   
+  //flow methods:
+  //MCEP = monte carlo event plane
+  AliFlowAnalysisWithMCEventPlane *mcep = new AliFlowAnalysisWithMCEventPlane();
+  if (MCEP)
+    {
+      mcep->Init();
+    }
+
+  //QC = Q-cumulants  
+  if(QC)
+    { 
+      AliFlowAnalysisWithQCumulants* qc = new AliFlowAnalysisWithQCumulants();
+      qc->CreateOutputObjects();
+    }
+  
+  //GFC = Generating Function Cumulants 
+  if(GFC)
+    {
+      AliFlowAnalysisWithCumulants* gfc = new AliFlowAnalysisWithCumulants();
+      gfc->CreateOutputObjects();
+    }
+  
+  //FQD = Fitting q-distribution 
+  AliFittingQDistribution* fqd = new AliFittingQDistribution();
+  if(FQD)
+    {
+      fqd->CreateOutputObjects();
+    }
+  
+  //LYZ1 = Lee-Yang Zeroes first run
+  AliFlowAnalysisWithLeeYangZeros* lyz1 = new AliFlowAnalysisWithLeeYangZeros();
+  if(LYZ1)
+    {
+      lyz1->SetFirstRun(kTRUE);
+      lyz1->SetUseSum(kTRUE);
+      lyz1->Init();
+    }
+
+  //LYZ2 = Lee-Yang Zeroes second run
+  AliFlowAnalysisWithLeeYangZeros* lyz2 = new AliFlowAnalysisWithLeeYangZeros();
+  if(LYZ2)
+    {
+      // read the input file from the first run 
+      TString inputFileNameLYZ2 = "outputLYZ1analysis.root" ;
+      TFile* inputFileLYZ2 = new TFile(inputFileNameLYZ2.Data(),"READ");
+      if(!inputFileLYZ2 || inputFileLYZ2->IsZombie()) { 
+       cerr << " ERROR: NO First Run file... " << endl ; }
+      else { 
+       TList* inputListLYZ2 = (TList*)inputFileLYZ2->Get("cobjLYZ1");  
+       if (!inputListLYZ2) {cout<<"list is NULL pointer!"<<endl;
+       }
+       else {
+         cout<<"LYZ2 input file/list read..."<<endl;
+         lyz2->SetFirstRunList(inputListLYZ2);
+         lyz2->SetFirstRun(kFALSE);
+         lyz2->SetUseSum(kTRUE);
+         lyz2->Init();
+       }
+      }
+    }
+
+ //LYZEP = Lee-Yang Zeroes event plane
+  AliFlowLYZEventPlane* ep = new AliFlowLYZEventPlane() ;
+  AliFlowAnalysisWithLYZEventPlane* lyzep = new AliFlowAnalysisWithLYZEventPlane();
+  if(LYZEP)
+    {
+      // read the input file from the second lyz run 
+      TString inputFileNameLYZEP = "outputLYZ2analysis.root" ;
+      TFile* inputFileLYZEP = new TFile(inputFileNameLYZEP.Data(),"READ");
+      if(!inputFileLYZEP || inputFileLYZEP->IsZombie()) { 
+       cerr << " ERROR: NO Second Run file... " << endl ; }
+      else { 
+       TList* inputListLYZEP = (TList*)inputFileLYZEP->Get("cobjLYZ2");  
+       if (!inputListLYZEP) {cout<<"list is NULL pointer!"<<endl;
+       }
+       else {
+         cout<<"LYZEP input file/list read..."<<endl;
+         ep   ->SetSecondRunList(inputListLYZEP);
+         lyzep->SetSecondRunList(inputListLYZEP);
+         ep   ->Init();
+         lyzep->Init();
+       }
+      }
+    }
+   
+  //SP = Scalar Product 
+  AliFlowAnalysisWithScalarProduct* sp = new AliFlowAnalysisWithScalarProduct();
+  if(SP)
+    {
+      sp->Init();
+    }
+
+
+
+  //------------------------------------------------------------------------
+  
+  //standard code
+  Int_t fCount = 0;
+  TString execDir(gSystem->pwd());
+  TString targetDir(dir);
+  TSystemDirectory* baseDir = new TSystemDirectory(".", dir);  
+  TList* dirList          = baseDir->GetListOfFiles();
+  if (!dirList) {
+    cout << endl << "No input files in: " << targetDir.Data() << endl;
+    break;
+  }
+  Int_t nDirs             = dirList->GetEntries();
+  cout<<endl;
+  cout<<"Int_t nDirs = "<<nDirs<<endl;
+  gSystem->cd(execDir);
+  
+  for(Int_t iDir=0;iDir<nDirs;++iDir)
+    {
+      TSystemFile* presentDir = (TSystemFile*)dirList->At(iDir);
+      if(!presentDir || !presentDir->IsDirectory() || 
+        strcmp(presentDir->GetName(), ".") == 0 || 
+        strcmp(presentDir->GetName(), "..") == 0) 
+       {
+         cout << endl; 
+         cout << "Directory (" << iDir << "): " << presentDir->GetName() << 
+           " - Skipping ... " << endl;
+         continue ;   
+       }
+      
+      if(offset > 0) { --offset ; continue ; }
+      if((aRuns > 0) && (fCount >= aRuns)) { break ; }
+      
+      TString presentDirName(dir); // aDataDir
+      presentDirName += presentDir->GetName();
+      presentDirName += "/";
+      //cerr<<" presentDirName = "<<presentDirName<<endl;
+      
+      TString fileName = presentDirName; 
+      fileName += "galice.root";
+      Long_t *id, *size, *flags, *modtime;
+      if(gSystem->GetPathInfo(fileName.Data(),id,size,flags,modtime)) 
+       { 
+         cout << " File : " << fileName << " does NOT exist ! - Skipping ... " 
+              << endl; 
+         continue; 
+       }
+      cout << endl ; cout << "Directory (" << iDir << "): " << presentDirName 
+          << " ... " << endl;
+      
+      //loop (simulations in the present dir) 
+      TSystemDirectory* evtsDir = new TSystemDirectory(".", 
+                                                      presentDirName.Data());
+      TList* fileList      = evtsDir->GetListOfFiles();
+      Int_t nFiles                 = fileList->GetEntries();
+      cout<<" Int_t nFiles = "<<nFiles<<endl;
+      gSystem->cd(execDir);      
+      for(Int_t iFiles=0; iFiles<nFiles; ++iFiles)
+       {
+         TSystemFile* presentFile = (TSystemFile*) fileList->At(iFiles);
+         TString presentFileName(presentDirName);
+         presentFileName += presentFile->GetName();
+         
+         if(!(presentFileName.Contains("Kinematics") && 
+              presentFileName.Contains("root"))) { continue ; }
+         
+         cout << " found: " << presentFileName.Data() << endl; 
+         
+         TFile* kineFile = new TFile(presentFileName.Data(), "READ"); 
+         // kineFile->ls();
+         Int_t nEvts = kineFile->GetNkeys() ; 
+         cout << "  . found: " << nEvts << " KineTree(s) in " << 
+           presentFileName.Data() << endl;
+         TList* kineEventsList = (TList*)kineFile->GetListOfKeys(); 
+         TTree* kTree;
+         TIter next(kineEventsList); 
+         TKey* key;
+         
+         //loop over the events
+         while( key=(TKey *)next() ) 
+           {
+             TDirectory* tDir = (TDirectory*)key->ReadObj();
+             if(!tDir) break;
+             
+             TString evtDir(tDir->GetName()); 
+             cout << "  . . found: " << tDir->GetName() << endl;
+             
+             kTree = (TTree *)tDir->Get("TreeK");
+             if(!kTree) break;
+             
+             Int_t nPart = kTree->GetEntries();
+             cout << "  . . . kTree " << fCount << " has " << nPart << 
+               " particles " << endl; 
+             
+             
+             
+             //-----------------------------------------------------------
+             //fill and save the flow event
+             
+             AliFlowEventSimple* fEvent = fEventMaker->FillTracks(kTree, cutsInt, cutsDiff);
+             
+             //pass the flow event to flow methods for analysis 
+             //MCEP
+             if (MCEP)
+               {
+                 //mcep->Make(fEvent,fEP);  //fix fEP
+                 cout<<"  --> MCEP analysis..."<<endl;
+               }
+             //QC
+             if(QC)
+               {  
+                 qc->Make(fEvent);
+                 cout<<"  --> QC analysis..."<<endl;
+               }
+             //GFC
+             if(GFC)
+               {
+                 gfc->Make(fEvent);
+                 cout<<"  --> GFC analysis..."<<endl;
+               }
+             //FQD 
+             if(FQD)
+               {
+                 fqd->Make(fEvent);
+                 cout<<"  --> FQD analysis..."<<endl;
+               }
+             //LYZ1
+             if(LYZ1)
+               {
+                 lyz1->Make(fEvent);
+                 cout<<"  --> LYZ1 analysis..."<<endl;
+               }
+             //LYZ2
+             if(LYZ2)
+               {
+                 lyz2->Make(fEvent);
+                 cout<<"  --> LYZ2 analysis..."<<endl;
+               }
+             //LYZEP
+             if(LYZEP)
+               {
+                 lyzep->Make(fEvent,ep);
+                 cout<<"  --> LYZEP analysis..."<<endl;
+               }
+             //SP
+             if(SP)
+               {
+                 sp->Make(fEvent);
+                 cout<<"  --> SP analysis..."<<endl;
+               }
+
+
+
+             //-----------------------------------------------------------
+             
+             fCount++;
+             delete kTree;
+             delete fEvent;
+           }
+         delete kineFile ;
+       }
+      delete evtsDir ;
+    }
+
+  //--------------------------------------------------------------
+  //calculating and storing the final results of flow analysis
+  //MCEP
+  if (MCEP)
+    {
+      mcep->Finish();
+      TString *outputFileNameMCEP = new TString("outputMCEPanalysis.root");
+      //mcep->WriteHistograms(outputFileNameMCEP); //add this method to MCEP classes
+      delete outputFileNameMCEP;
+    }
+  //QC
+  if(QC)
+    {
+      qc->Finish();
+      TString *outputFileNameQC = new TString("outputQCanalysis.root");
+      qc->WriteHistograms(outputFileNameQC);
+      delete outputFileNameQC;
+    }
+  //GFC
+  if(GFC)
+    {
+      gfc->Finish();
+      TString *outputFileNameGFC = new TString("outputGFCanalysis.root");
+      gfc->WriteHistograms(outputFileNameGFC);
+      delete outputFileNameGFC;
+    }
+  //FQD
+  if(FQD)
+    {
+      fqd->Finish();
+      TString *outputFileNameFQD = new TString("outputFQDanalysis.root");
+      fqd->WriteHistograms(outputFileNameFQD);
+      delete outputFileNameFQD;
+    }
+  //LYZ1
+  if(LYZ1)
+    {
+      lyz1->Finish();
+      TString *outputFileNameLYZ1 = new TString("outputLYZ1analysis.root");
+      lyz1->WriteHistograms(outputFileNameLYZ1);
+      delete outputFileNameLYZ1;
+    }
+  //LYZ2
+  if(LYZ2)
+    {
+      lyz2->Finish();
+      TString *outputFileNameLYZ2 = new TString("outputLYZ2analysis.root");
+      lyz2->WriteHistograms(outputFileNameLYZ2);
+      delete outputFileNameLYZ2;
+    }
+  //LYZEP
+  if(LYZEP)
+    {
+      lyzep->Finish();
+      TString *outputFileNameLYZEP = new TString("outputLYZEPanalysis.root");
+      lyzep->WriteHistograms(outputFileNameLYZEP);
+      delete outputFileNameLYZEP;
+    }
+  //SP
+  if(SP)
+    {
+      sp->Finish();
+      TString *outputFileNameSP = new TString("outputSPanalysis.root");
+      sp->WriteHistograms(outputFileNameSP);
+      delete outputFileNameSP;
+    }
+
+
+
+  //--------------------------------------------------------------
+  
+  cout << endl;
+  cout << " Finished ... " << endl;
+  cout << endl;
+  
+  timer.Stop();
+  cout << endl;
+  timer.Print();
+}