#include "Riostream.h" //needed as include
#include "TChain.h"
#include "TTree.h"
+#include "TProfile.h"
#include "TFile.h"
#include "TList.h"
#include "AliAnalysisTaskLYZEventPlane.h"
#include "AliFlowEventSimpleMaker.h"
+#include "AliFlowCommonHist.h"
+#include "AliFlowCommonHistResults.h"
#include "AliFlowLYZEventPlane.h"
#include "AliFlowAnalysisWithLYZEventPlane.h"
ClassImp(AliAnalysisTaskLYZEventPlane)
//________________________________________________________________________
-AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane(const char *name) :
+AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane(const char *name, Bool_t on) :
AliAnalysisTask(name, ""),
fESD(NULL),
fAOD(NULL),
fCFManager1(NULL),
fCFManager2(NULL),
fListHistos(NULL),
- fSecondRunFile(NULL)
+ fSecondRunFile(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(on)
{
// Constructor
cout<<"AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane(const char *name)"<<endl;
DefineInput(0, TChain::Class());
DefineInput(1, TList::Class());
// Output slot #0 writes into a TList container
- DefineOutput(0, TList::Class());
+ DefineOutput(0, TList::Class());
+ if(on) {
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, TList::Class()); }
}
//________________________________________________________________________
fCFManager1(NULL),
fCFManager2(NULL),
fListHistos(NULL),
- fSecondRunFile(NULL)
+ fSecondRunFile(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(kFALSE)
{
// Constructor
cout<<"AliAnalysisTaskLYZEventPlane::AliAnalysisTaskLYZEventPlane()"<<endl;
fLyz = new AliFlowAnalysisWithLYZEventPlane() ;
// Get data from input slot
- fSecondRunFile = (TFile*)GetInputData(1);
- cerr<<"fSecondRunFile ("<<fSecondRunFile<<")"<<endl;
- if (fSecondRunFile) cerr<<"fSecondRunFile -> IsOpen() = "<<fSecondRunFile -> IsOpen()<<endl;
-
- fLyzEp -> SetSecondRunFile(fSecondRunFile);
- fLyz -> SetSecondRunFile(fSecondRunFile);
+ TList* pSecondRunList = (TList*)GetInputData(1);
+ if (pSecondRunList) {
+ fLyzEp -> SetSecondRunList(pSecondRunList);
+ fLyz -> SetSecondRunList(pSecondRunList);
+ } else { cout<<"No Second run List!"<<endl; exit(0); }
fLyzEp-> Init();
fLyz-> Init();
}
PostData(0,fListHistos);
+ if (fQA) {
+ PostData(1,fQAInt);
+ PostData(2,fQADiff); }
}
//________________________________________________________________________
void AliAnalysisTaskLYZEventPlane::Terminate(Option_t *)
{
// Called once at the end of the query
- fLyz->Finish();
+ AliFlowAnalysisWithLYZEventPlane* fLyzTerm = new AliFlowAnalysisWithLYZEventPlane() ;
+ fListHistos = (TList*)GetOutputData(0);
+ cout << "histogram list in Terminate" << endl;
+ if (fListHistos) {
+ //Get the common histograms from the output list
+ AliFlowCommonHist *pCommonHist = dynamic_cast<AliFlowCommonHist*>
+ (fListHistos->FindObject("AliFlowCommonHistLYZEP"));
+ AliFlowCommonHistResults *pCommonHistResults = dynamic_cast<AliFlowCommonHistResults*>
+ (fListHistos->FindObject("AliFlowCommonHistResultsLYZEP"));
+
+ TProfile* pHistProR0theta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("First_FlowPro_r0theta_LYZ"));
+
+ TProfile* pHistProFlow = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("FlowPro_VPt_LYZEP"));
+
+ TH1F* pHistQsumforChi = dynamic_cast<TH1F*>
+ (fListHistos->FindObject("Flow_QsumforChi_LYZEP"));
+
+ if (pCommonHist && pCommonHistResults && pHistProR0theta &&
+ pHistProFlow && pHistQsumforChi ) {
+ fLyzTerm->SetCommonHists(pCommonHist);
+ fLyzTerm->SetCommonHistsRes(pCommonHistResults);
+ fLyzTerm->SetFirstr0theta(pHistProR0theta);
+ fLyzTerm->SetHistProFlow(pHistProFlow);
+ fLyzTerm->SetHistQsumforChi(pHistQsumforChi);
+ fLyzTerm->Finish();
+ PostData(0,fListHistos);
+ } else {
+ cout<<"WARNING: Histograms needed to run Finish() are not accessable!"<<endl;
+ }
- TList* fOutHistos = (TList*)GetOutputData(0);
- if (fOutHistos) { fOutHistos->Print(); }
- else { cout<<"hostogram list pointer in Terminate is empty"<<endl; }
+ fListHistos->Print();
+ } else { cout << "histogram list pointer is empty" << endl;}
- //delete fLyz;
- //delete fLyzEp;
- //delete fEventMaker;
+ cout<<".....finished"<<endl;
}
class AliAnalysisTaskLYZEventPlane : public AliAnalysisTask {
public:
AliAnalysisTaskLYZEventPlane();
- AliAnalysisTaskLYZEventPlane(const char *name);
+ AliAnalysisTaskLYZEventPlane(const char *name, Bool_t QAon = kFALSE);
virtual ~AliAnalysisTaskLYZEventPlane();
virtual void ConnectInputData(Option_t *);
AliCFManager* GetCFManager1() const {return this->fCFManager1; }
void SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; }
AliCFManager* GetCFManager2() const {return this->fCFManager2; }
+ void SetQAList1(TList* list) {this->fQAInt = list; }
+ TList* GetQAList1() {return this->fQAInt; }
+ void SetQAList2(TList* list) {this->fQADiff = list; }
+ TList* GetQAList2() {return this->fQADiff; }
+ void SetQAOn(Bool_t kt) {this->fQA = kt; }
+ Bool_t GetQAOn() {return this->fQA; }
private:
AliAnalysisTaskLYZEventPlane(const AliAnalysisTaskLYZEventPlane& aAnalysis);
AliAnalysisTaskLYZEventPlane& operator=(const AliAnalysisTaskLYZEventPlane& aAnalysis);
- AliESDEvent *fESD; //ESD object
- AliAODEvent *fAOD; //AOD object
+ AliESDEvent* fESD; //ESD object
+ AliAODEvent* fAOD; //AOD object
TString fAnalysisType; //string to set the kind of input for the analysis: ESD, AOD or MC
AliFlowLYZEventPlane* fLyzEp; //LYZ EP object
AliFlowAnalysisWithLYZEventPlane* fLyz; //LYZ EP analysis object
AliFlowEventSimpleMaker* fEventMaker; //FlowEventSimple maker object
AliCFManager* fCFManager1; //Correction framework manager
AliCFManager* fCFManager2; //Correction framework manager
- TList* fListHistos; //collection of output hists
- TFile* fSecondRunFile; //output from the second LYZ loop
+ TList* fListHistos; //collection of output hists
+ TFile* fSecondRunFile; //output from the second LYZ loop
+ TList* fQAInt; // QA histogram list
+ TList* fQADiff; // QA histogram list
+
+ Bool_t fQA; // flag to set the filling of the QA hostograms
ClassDef(AliAnalysisTaskLYZEventPlane, 1); // example of analysis
};
#include "TTree.h"
#include "TFile.h"
#include "TList.h"
+#include "TProfile.h"
class AliAnalysisTask;
#include "AliAnalysisTaskLeeYangZeros.h"
#include "AliFlowEventSimpleMaker.h"
+#include "AliFlowCommonHist.h"
+#include "AliFlowCommonHistResults.h"
+#include "AliFlowLYZHist1.h"
+#include "AliFlowLYZHist2.h"
#include "AliFlowAnalysisWithLeeYangZeros.h"
// AliAnalysisTaskLeeYangZeros:
ClassImp(AliAnalysisTaskLeeYangZeros)
//________________________________________________________________________
-AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros(const char *name, Bool_t firstrun) :
+AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros(const char *name, Bool_t firstrun, Bool_t on) :
AliAnalysisTask(name, ""),
fESD(0),
fAOD(0),
fEventMaker(0),
fFirstRunFile(0),
fListHistos(NULL),
- fFirstRunLYZ(firstrun), //set boolean for firstrun to initial value
- fUseSumLYZ(kTRUE) //set boolean for use sum to initial value
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fFirstRunLYZ(firstrun), //set boolean for firstrun to initial value
+ fUseSumLYZ(kTRUE), //set boolean for use sum to initial value
+ fQA(on)
{
// Constructor
cout<<"AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros(const char *name)"<<endl;
if (!firstrun) DefineInput(1, TList::Class()); //for second loop
// Output slot #0 writes into a TList container
DefineOutput(0, TList::Class());
-}
+ if(on) {
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, TList::Class()); }
+}
+
-/*
//________________________________________________________________________
AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros() :
fESD(0),
fEventMaker(0),
fFirstRunFile(0),
fListHistos(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
fFirstRunLYZ(kTRUE), //set boolean for firstrun to initial value
- fUseSumLYZ(kTRUE) //set boolean for use sum to initial value
+ fUseSumLYZ(kTRUE), //set boolean for use sum to initial value
+ fQA(kFALSE)
{
// Constructor
- cout<<"AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros(const char *name)"<<endl;
+ cout<<"AliAnalysisTaskLeeYangZeros::AliAnalysisTaskLeeYangZeros()"<<endl;
}
-*/
-
//________________________________________________________________________
AliAnalysisTaskLeeYangZeros::~AliAnalysisTaskLeeYangZeros()
{
}
-
-
//________________________________________________________________________
void AliAnalysisTaskLeeYangZeros::ConnectInputData(Option_t *)
{
// Get data from input slot 1
if (GetNinputs() == 2) { //if there are two input slots
- fFirstRunFile = (TFile*)GetInputData(1);
- cerr<<"fFirstRunFile ("<<fFirstRunFile<<")"<<endl;
- if (fFirstRunFile) { cerr<<"fFirstRunFile -> IsOpen() = "<<fFirstRunFile -> IsOpen()<<endl;}
- else { cerr<<"fFirstRunFile has a NULL pointer!!"<<endl; exit(0);}
- fLyz -> SetFirstRunFile(fFirstRunFile);
+ TList* pFirstRunList = (TList*)GetInputData(1);
+ if (pFirstRunList) {
+ fLyz -> SetFirstRunList(pFirstRunList);
+ } else { cout<<"No first run List!"<<endl; exit(0); }
}
- fLyz-> Init();
+ fLyz -> Init();
if (fLyz->GetHistList()) {
- fLyz->GetHistList()->Print();
- fListHistos = fLyz->GetHistList();
- fListHistos->Print();
+ fListHistos = fLyz->GetHistList();
+ fListHistos->Print();
}
else {Printf("ERROR: Could not retrieve histogram list"); }
delete fEvent;
}
- //PostData(0,fListHistos); //here for CAF
-
+ PostData(0,fListHistos); //here for CAF
+ if (fQA) {
+ PostData(1,fQAInt);
+ PostData(2,fQADiff); }
}
//________________________________________________________________________
void AliAnalysisTaskLeeYangZeros::Terminate(Option_t *)
{
// Called once at the end of the query
- if (GetNinputs() == 2) {
- cerr<<"fFirstRunFile -> IsOpen() = "<<fFirstRunFile -> IsOpen()<<endl;
- }
-
- fLyz->Finish(); //remove for CAF
- PostData(0,fListHistos); //remove for CAF
-
- //print histogram list:
- TList* fOutListHistos = (TList*)GetOutputData(0);
+ AliFlowAnalysisWithLeeYangZeros* fLyzTerm = new AliFlowAnalysisWithLeeYangZeros() ;
+ fLyzTerm -> SetFirstRun(GetFirstRunLYZ()); //set first run true or false
+ fLyzTerm -> SetUseSum(GetUseSumLYZ()); //set use sum true or false
+
+ fListHistos = (TList*)GetOutputData(0);
cout << "histogram list in Terminate" << endl;
- if (fOutListHistos) {
- //fOutListHistos->Print(); //gives error for secondrun??
+ if (fListHistos) {
+ //Get the common histograms from the output list
+ AliFlowCommonHist *pCommonHist = dynamic_cast<AliFlowCommonHist*>
+ (fListHistos->FindObject("AliFlowCommonHistLYZ"));
+
+ AliFlowCommonHistResults *pCommonHistResults = dynamic_cast<AliFlowCommonHistResults*>
+ (fListHistos->FindObject("AliFlowCommonHistResultsLYZ"));
+
+ TProfile* pHistProR0theta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("First_FlowPro_r0theta_LYZ"));
+
+ TH1F* pHistQsumforChi = dynamic_cast<TH1F*>
+ (fListHistos->FindObject("Flow_QsumforChi_LYZEP"));
+
+ //define histograms for first and second run
+ TProfile* pHistProVtheta = 0;
+ TProfile* pHistProReDenom = 0;
+ TProfile* pHistProImDenom = 0;
+ TProfile* pHistProReDtheta = 0;
+ TProfile* pHistProImDtheta = 0;
+ TProfile* pHistProVeta = 0;
+ TProfile* pHistProVPt = 0;
+
+ AliFlowLYZHist1 *pLYZHist1[5] = {0}; //array of pointers to AliFlowLYZHist1
+ AliFlowLYZHist2 *pLYZHist2[5] = {0}; //array of pointers to AliFlowLYZHist2
+
+ if (GetFirstRunLYZ()) { //for firstrun
+ //Get the histograms from the output list
+ for(Int_t theta = 0;theta<5;theta++){
+ TString name = "AliFlowLYZHist1_";
+ name += theta;
+ pLYZHist1[theta] = dynamic_cast<AliFlowLYZHist1*>
+ (fListHistos->FindObject(name));
+ }
+ pHistProVtheta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("First_FlowPro_Vtheta_LYZ"));
+
+ //Set the histogram pointers and call Finish()
+ if (pCommonHist && pCommonHistResults && pLYZHist1 &&
+ pHistProVtheta && pHistProR0theta && pHistQsumforChi ) {
+ fLyzTerm->SetCommonHists(pCommonHist);
+ fLyzTerm->SetCommonHistsRes(pCommonHistResults);
+ fLyzTerm->SetHist1(pLYZHist1);
+ fLyzTerm->SetHistProVtheta(pHistProVtheta);
+ fLyzTerm->SetHistProR0theta(pHistProR0theta);
+ fLyzTerm->SetHistQsumforChi(pHistQsumforChi);
+ fLyzTerm->Finish();
+ PostData(0,fListHistos);
+ } else {
+ cout<<"WARNING: Histograms needed to run Finish() firstrun are not accessable!"<<endl;
+ }
+ } else { //for second run
+ //Get the histograms from the output list
+ for(Int_t theta = 0;theta<5;theta++){
+ TString name = "AliFlowLYZHist2_";
+ name += theta;
+ pLYZHist2[theta] = dynamic_cast<AliFlowLYZHist2*>
+ (fListHistos->FindObject(name));
+ }
+
+ pHistProReDenom = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_ReDenom_LYZ"));
+ pHistProImDenom = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_ImDenom_LYZ"));
+
+ pHistProReDtheta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_ReDtheta_LYZ"));
+ pHistProImDtheta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_ImDtheta_LYZ"));
+
+ pHistProVeta = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_Veta_LYZ"));
+ pHistProVPt = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("Second_FlowPro_VPt_LYZ"));
+
+ //Set the histogram pointers and call Finish()
+ if (pCommonHist && pCommonHistResults && pLYZHist2 && pHistProR0theta &&
+ pHistProReDenom && pHistProImDenom && pHistProVeta && pHistProVPt) {
+ fLyzTerm->SetCommonHists(pCommonHist);
+ fLyzTerm->SetCommonHistsRes(pCommonHistResults);
+ fLyzTerm->SetHist2(pLYZHist2);
+ fLyzTerm->SetHistProR0theta(pHistProR0theta);
+ fLyzTerm->SetHistProReDenom(pHistProReDenom);
+ fLyzTerm->SetHistProImDenom(pHistProImDenom);
+ fLyzTerm->SetHistProReDtheta(pHistProReDtheta);
+ fLyzTerm->SetHistProImDtheta(pHistProImDtheta);
+ fLyzTerm->SetHistProVeta(pHistProVeta);
+ fLyzTerm->SetHistProVPt(pHistProVPt);
+ fLyzTerm->SetHistQsumforChi(pHistQsumforChi);
+ fLyzTerm->Finish();
+ PostData(0,fListHistos);
+ } else {
+ cout<<"WARNING: Histograms needed to run Finish() secondrun are not accessable!"<<endl;
+ }
+ }
+
+ fListHistos->Print();
}
- else {
- cout << "histgram list pointer is empty" << endl;}
-
- //delete fLyz;
- //delete fEventMaker;
+ else { cout << "histogram list pointer is empty" << endl;}
cout<<".....finished"<<endl;
}
// analysis task for
// Lee Yang Zeroes method
// Author:
-// Naomi van der Kolk (kolk@nikhef.nl)
+// Naomi van der Kolk (kolk@nikhef.nl)
class AliESDEvent;
class AliAODEvent;
class AliAnalysisTaskLeeYangZeros : public AliAnalysisTask {
public:
- //AliAnalysisTaskLeeYangZeros();
- AliAnalysisTaskLeeYangZeros(const char *name = "AliAnalysisTaskLeeYangZeros", Bool_t firstrun = kTRUE);
+ AliAnalysisTaskLeeYangZeros();
+ AliAnalysisTaskLeeYangZeros(const char *name, Bool_t firstrun, Bool_t QAon);
virtual ~AliAnalysisTaskLeeYangZeros();
virtual void ConnectInputData(Option_t *);
AliCFManager* GetCFManager1() {return this->fCFManager1; }
void SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; }
AliCFManager* GetCFManager2() {return this->fCFManager2; }
-
+ void SetQAList1(TList* list) {this->fQAInt = list; }
+ TList* GetQAList1() {return this->fQAInt; }
+ void SetQAList2(TList* list) {this->fQADiff = list; }
+ TList* GetQAList2() {return this->fQADiff; }
+ void SetQAOn(Bool_t kt) {this->fQA = kt; }
+ Bool_t GetQAOn() {return this->fQA; }
private:
AliAnalysisTaskLeeYangZeros(const AliAnalysisTaskLeeYangZeros& aAnalysis);
AliAnalysisTaskLeeYangZeros& operator=(const AliAnalysisTaskLeeYangZeros& aAnalysis);
- AliESDEvent* fESD; //ESD object
- AliAODEvent* fAOD; //AOD object
- TString fAnalysisType; //string to select which kind of input to analyse: ESD, AOD or MC
+ AliESDEvent* fESD; // ESD object
+ AliAODEvent* fAOD; // AOD object
+ TString fAnalysisType; // string to select which kind of input to analyse: ESD, AOD or MC
AliCFManager* fCFManager1; // correction framework manager
AliCFManager* fCFManager2; // correction framework manager
- AliFlowAnalysisWithLeeYangZeros* fLyz; //LYZ analysis object
- AliFlowEventSimpleMaker* fEventMaker; //FlowEventSimple maker object
+ AliFlowAnalysisWithLeeYangZeros* fLyz; // LYZ analysis object
+ AliFlowEventSimpleMaker* fEventMaker; // FlowEventSimple maker object
TFile* fFirstRunFile; // file from the first loop over events
TList* fListHistos; // collection of output
+ TList* fQAInt; // QA histogram list
+ TList* fQADiff; // QA histogram list
//flags
- Bool_t fFirstRunLYZ ; //! flag for lyz analysis
- Bool_t fUseSumLYZ ; //! flag for lyz analysis
+ Bool_t fFirstRunLYZ ; // flag for lyz analysis
+ Bool_t fUseSumLYZ ; // flag for lyz analysis
+ Bool_t fQA; // flag to set the filling of the QA hostograms
- ClassDef(AliAnalysisTaskLeeYangZeros, 1); // example of analysis
+ ClassDef(AliAnalysisTaskLeeYangZeros, 1); //AliAnalysisTaskLeeYangZeros class object
};
#endif
#include "Riostream.h" //needed as include
#include "TChain.h"
+#include "TProfile.h"
#include "TTree.h"
#include "TFile.h" //needed as include
#include "TList.h"
#include "AliAnalysisTaskMCEventPlane.h"
#include "AliFlowEventSimpleMaker.h"
+#include "AliFlowCommonHist.h"
+#include "AliFlowCommonHistResults.h"
#include "AliFlowAnalysisWithMCEventPlane.h"
// AliAnalysisTaskMCEventPlane:
ClassImp(AliAnalysisTaskMCEventPlane)
//________________________________________________________________________
-AliAnalysisTaskMCEventPlane::AliAnalysisTaskMCEventPlane(const char *name) :
+AliAnalysisTaskMCEventPlane::AliAnalysisTaskMCEventPlane(const char *name, Bool_t on) :
AliAnalysisTask(name, ""),
fESD(0),
fAOD(0),
fCFManager2(NULL),
fMc(0),
fEventMaker(0),
- fListHistos(NULL)
+ fListHistos(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(on)
{
// Constructor
cout<<"AliAnalysisTaskMCEventPlane::AliAnalysisTaskMCEventPlane(const char *name)"<<endl;
// Input slot #0 works with a TChain
DefineInput(0, TChain::Class());
// Output slot #0 writes into a TList container
- DefineOutput(0, TList::Class());
+ DefineOutput(0, TList::Class());
+ if(on) {
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, TList::Class()); }
}
//________________________________________________________________________
fCFManager2(NULL),
fMc(0),
fEventMaker(0),
- fListHistos(NULL)
+ fListHistos(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(kFALSE)
{
// Constructor
cout<<"AliAnalysisTaskMCEventPlane::AliAnalysisTaskMCEventPlane()"<<endl;
// Disable all branches and enable only the needed ones
if (fAnalysisType == "MC") {
- cout<<"!!!!!reading MC kinematics only"<<endl;
// we want to process only MC
tree->SetBranchStatus("*", kFALSE);
}
else if (fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" ) {
- cout<<"!!!!!reading the ESD only"<<endl;
tree->SetBranchStatus("*", kFALSE);
tree->SetBranchStatus("Tracks.*", kTRUE);
}
else if (fAnalysisType == "AOD") {
- cout<<"!!!!!reading the AOD only"<<endl;
AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
if (!aodH) {
fMc-> Init();
if (fMc->GetHistList()) {
- fMc->GetHistList()->Print();
+ //fMc->GetHistList()->Print();
fListHistos = fMc->GetHistList();
fListHistos->Print();
}
}
Double_t fRP = hdh->ReactionPlaneAngle();
- cout<<"The reactionPlane is "<<hdh->ReactionPlaneAngle()<<endl;
+ //cout<<"The reactionPlane is "<<hdh->ReactionPlaneAngle()<<endl;
if (fAnalysisType == "MC") {
// analysis
delete fEvent;
}
- //PostData(0,fListHistos); //here for CAF
+ PostData(0,fListHistos); //here for CAF
+ if (fQA) {
+ PostData(1,fQAInt);
+ PostData(2,fQADiff); }
}
+
//________________________________________________________________________
void AliAnalysisTaskMCEventPlane::Terminate(Option_t *)
{
// Called once at the end of the query
- fMc->Finish();
- PostData(0,fListHistos);
-
- //delete fMc;
- //delete fEventMaker;
+ AliFlowAnalysisWithMCEventPlane* fMcTerm = new AliFlowAnalysisWithMCEventPlane() ;
+
+ //Get output data
+ fListHistos = (TList*)GetOutputData(0);
+ cout << "histgram list in Terminate" << endl;
+ if (fListHistos) {
+ //Get the common histograms from the output list
+ AliFlowCommonHist *pCommonHists = dynamic_cast<AliFlowCommonHist*>
+ (fListHistos->FindObject("AliFlowCommonHistMCEP"));
+ AliFlowCommonHistResults *pCommonHistResults =
+ dynamic_cast<AliFlowCommonHistResults*>
+ (fListHistos->FindObject("AliFlowCommonHistResultsMCEP"));
+
+ TProfile *pHistProFlow = dynamic_cast<TProfile*>
+ (fListHistos->FindObject("FlowPro_VPt_MCEP"));
+
+ if (pCommonHists && pCommonHistResults && pHistProFlow) {
+ fMcTerm->SetCommonHists(pCommonHists);
+ fMcTerm->SetCommonHistsRes(pCommonHistResults);
+ fMcTerm->SetHistProFlow(pHistProFlow);
+ fMcTerm->Finish();
+ PostData(0,fListHistos);
+ } else {
+ cout<<"WARNING: Histograms needed to run Finish() are not accessable!"<<endl; }
+
+ fListHistos->Print();
+ } else { cout << "histogram list pointer is empty" << endl;}
+
+ cout<<"...finished."<<endl;
}
public:
AliAnalysisTaskMCEventPlane();
- AliAnalysisTaskMCEventPlane(const char *name);
+ AliAnalysisTaskMCEventPlane(const char *name,Bool_t QAon = kFALSE);
virtual ~AliAnalysisTaskMCEventPlane();
virtual void ConnectInputData(Option_t *);
virtual void Exec(Option_t *option);
virtual void Terminate(Option_t *);
- void SetAnalysisType(TString type) { this->fAnalysisType = type; }
- TString GetAnalysisType() const { return this->fAnalysisType; }
+ void SetAnalysisType(TString type) { this->fAnalysisType = type; }
+ TString GetAnalysisType() const { return this->fAnalysisType; }
void SetCFManager1(AliCFManager* cfmgr) {this->fCFManager1 = cfmgr; }
AliCFManager* GetCFManager1() {return this->fCFManager1; }
void SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; }
AliCFManager* GetCFManager2() {return this->fCFManager2; }
+ void SetQAList1(TList* list) {this->fQAInt = list; }
+ TList* GetQAList1() {return this->fQAInt; }
+ void SetQAList2(TList* list) {this->fQADiff = list; }
+ TList* GetQAList2() {return this->fQADiff; }
+ void SetQAOn(Bool_t kt) {this->fQA = kt; }
+ Bool_t GetQAOn() {return this->fQA; }
private:
AliAnalysisTaskMCEventPlane(const AliAnalysisTaskMCEventPlane& aAnalysis);
AliAnalysisTaskMCEventPlane& operator=(const AliAnalysisTaskMCEventPlane& aAnalysis);
- AliESDEvent *fESD; // ESD object
- AliAODEvent *fAOD; // AOD object
- TString fAnalysisType; // can be MC, ESD or AOD
- AliCFManager* fCFManager1; // correction framework manager
- AliCFManager* fCFManager2; // correction framework manager
- AliFlowAnalysisWithMCEventPlane* fMc; // MC EP analysis object
- AliFlowEventSimpleMaker* fEventMaker; // FlowEventSimple maker object
+ AliESDEvent* fESD; // ESD object
+ AliAODEvent* fAOD; // AOD object
+ TString fAnalysisType; // can be MC, ESD or AOD
+ AliCFManager* fCFManager1; // correction framework manager
+ AliCFManager* fCFManager2; // correction framework manager
+ AliFlowAnalysisWithMCEventPlane* fMc; // MC EP analysis object
+ AliFlowEventSimpleMaker* fEventMaker; // FlowEventSimple maker object
- TList* fListHistos; // collection of output
+ TList* fListHistos; // collection of output
+ TList* fQAInt; // QA histogram list
+ TList* fQADiff; // QA histogram list
- ClassDef(AliAnalysisTaskMCEventPlane, 1); // example of analysis
+ Bool_t fQA; // flag to set the filling of the QA hostograms
+
+
+ ClassDef(AliAnalysisTaskMCEventPlane, 1); // AliAnalysisTaskMCEventPlane class object
};
#endif
#include "Riostream.h" //needed as include
#include "TChain.h"
#include "TTree.h"
-//#include "TFile.h" //needed as include
+#include "TFile.h" //needed as include
#include "TList.h"
ClassImp(AliAnalysisTaskScalarProduct)
//________________________________________________________________________
-AliAnalysisTaskScalarProduct::AliAnalysisTaskScalarProduct(const char *name) :
+AliAnalysisTaskScalarProduct::AliAnalysisTaskScalarProduct(const char *name, Bool_t on) :
AliAnalysisTask(name, ""),
fESD(NULL),
fAOD(NULL),
fAnalysisType("ESD"),
fCFManager1(NULL),
fCFManager2(NULL),
- fListHistos(NULL)
+ fListHistos(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(on)
{
// Constructor
cout<<"AliAnalysisTaskScalarProduct::AliAnalysisTaskScalarProduct(const char *name)"<<endl;
DefineInput(0, TChain::Class());
// Output slot #0 writes into a TList container
DefineOutput(0, TList::Class());
-
+ if(on) {
+ DefineOutput(1, TList::Class());
+ DefineOutput(2, TList::Class()); }
}
//________________________________________________________________________
fAnalysisType("ESD"),
fCFManager1(NULL),
fCFManager2(NULL),
- fListHistos(NULL)
+ fListHistos(NULL),
+ fQAInt(NULL),
+ fQADiff(NULL),
+ fQA(kFALSE)
{
// Constructor
cout<<"AliAnalysisTaskScalarProduct::AliAnalysisTaskScalarProduct()"<<endl;
} else {
// Disable all branches and enable only the needed ones
if (fAnalysisType == "MC") {
- cout<<"!!!!!reading MC kinematics only"<<endl;
// we want to process only MC
tree->SetBranchStatus("*", kFALSE);
}
}
else if (fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1") {
- cout<<"!!!!!reading the ESD only"<<endl;
tree->SetBranchStatus("*", kFALSE);
tree->SetBranchStatus("Tracks.*", kTRUE);
fESD = esdH->GetEvent();
}
else if (fAnalysisType == "AOD") {
- cout<<"!!!!!reading the AOD only"<<endl;
AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
if (!aodH) {
if (fSP->GetHistList()) {
- fSP->GetHistList()->Print();
+ //fSP->GetHistList()->Print();
fListHistos = fSP->GetHistList();
fListHistos->Print();
}
//fListHistos->Print();
PostData(0,fListHistos);
+ if (fQA) {
+ PostData(1,fQAInt);
+ PostData(2,fQADiff); }
}
//________________________________________________________________________
void AliAnalysisTaskScalarProduct::Terminate(Option_t *)
{
- // Called once at the end of the query
+ // Called once at the end of the query -- do not call in case of CAF
// fSP->Finish();
// PostData(0,fListHistos);
fListHistos = (TList*)GetOutputData(0);
cout << "histgram list in Terminate" << endl;
- if (fListHistos)
- {
- fListHistos->Print();
-// AliFlowCommonHist *SDHistClass = dynamic_cast<AliFlowCommonHist*>
-// (fListHistos->FindObject("SP"));
-// if (SDHistClass)
-// {
-// // open a file and write the CLASS to the file
-// TFile *f = new TFile("mytest.root","recreate");
-// SDHistClass->Write();
-// SDHistClass->Print();
-// f->Close();
-// }
-// else
-// {
-// cout << "SD class pointer is NULL" << endl;
-// }
- }
- else
- {
- cout << "histgram list pointer is empty" << endl;
- }
-// delete fSP;
-// delete fEventMaker;
+ if (fListHistos) {
+ fListHistos->Print();
+ }
+ else { cout << "histgram list pointer is empty" << endl; }
+
}
class AliAnalysisTaskScalarProduct : public AliAnalysisTask {
public:
AliAnalysisTaskScalarProduct();
- AliAnalysisTaskScalarProduct(const char *name);
+ AliAnalysisTaskScalarProduct(const char *name,Bool_t QAon = kFALSE);
virtual ~AliAnalysisTaskScalarProduct();
virtual void ConnectInputData(Option_t *);
AliCFManager* GetCFManager1() {return this->fCFManager1; }
void SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; }
AliCFManager* GetCFManager2() {return this->fCFManager2; }
+ void SetQAList1(TList* list) {this->fQAInt = list; }
+ TList* GetQAList1() {return this->fQAInt; }
+ void SetQAList2(TList* list) {this->fQADiff = list; }
+ TList* GetQAList2() {return this->fQADiff; }
+ void SetQAOn(Bool_t kt) {this->fQA = kt; }
+ Bool_t GetQAOn() {return this->fQA; }
private:
AliAnalysisTaskScalarProduct(const AliAnalysisTaskScalarProduct& aAnalysisTask);
AliAnalysisTaskScalarProduct& operator=(const AliAnalysisTaskScalarProduct& aAnalysisTask);
- AliESDEvent *fESD; // ESD object
- AliAODEvent *fAOD; // AOD object
+ AliESDEvent* fESD; // ESD object
+ AliAODEvent* fAOD; // AOD object
AliFlowAnalysisWithScalarProduct* fSP; // analysis object
AliFlowEventSimpleMaker* fEventMaker; // FlowEventSimple maker object
- TString fAnalysisType; // can be MC, ESD or AOD
+ TString fAnalysisType; // can be MC, ESD or AOD
AliCFManager* fCFManager1; // correction framework manager
AliCFManager* fCFManager2; // correction framework manager
- TList *fListHistos; // collection of output
+ TList* fListHistos; // collection of output
+ TList* fQAInt; // QA histogram list
+ TList* fQADiff; // QA histogram list
+
+ Bool_t fQA; // flag to set the filling of the QA hostograms
ClassDef(AliAnalysisTaskScalarProduct, 1); // example of analysis
};
#include "TCanvas.h" //needed as include
#include "TLegend.h" //needed as include
#include "TProfile.h" //needed as include
+#include "TVector2.h"
class TH1F;
//-----------------------------------------------------------------------
AliFlowAnalysisWithLYZEventPlane::AliFlowAnalysisWithLYZEventPlane():
- fSecondRunFile(0),
- fSecondRunFileName(0),
- fHistList(NULL),
- fSecondReDtheta(NULL),
- fSecondImDtheta(NULL),
- fFirstr0theta(NULL),
- fSecondVPt(NULL),
- fHistProFlow(NULL),
- fHistProFlow2(NULL),
- fHistProWr(NULL),
- fHistProWrCorr(NULL),
- fHistFlow(NULL),
- fHistDeltaPhi(NULL),
- fHistDeltaPhi2(NULL),
- fHistDeltaPhihere(NULL),
- fHistPhiEP(NULL),
- fHistPhiEPhere(NULL),
- fHistPhiLYZ(NULL),
- fHistPhiLYZ2(NULL),
- fCommonHists(NULL),
- fCommonHistsRes(NULL),
- fEventNumber(0),
- fQsum(NULL),
- fQ2sum(0)
+ fHistList(NULL),
+ fSecondRunList(NULL),
+ fSecondReDtheta(NULL),
+ fSecondImDtheta(NULL),
+ fFirstr0theta(NULL),
+ fHistProFlow(NULL),
+ fHistProFlow2(NULL),
+ fHistProWr(NULL),
+ fHistProWrCorr(NULL),
+ fHistQsumforChi(NULL),
+ fHistDeltaPhi(NULL),
+ fHistDeltaPhi2(NULL),
+ fHistDeltaPhihere(NULL),
+ fHistPhiEP(NULL),
+ fHistPhiEPhere(NULL),
+ fHistPhiLYZ(NULL),
+ fHistPhiLYZ2(NULL),
+ fCommonHists(NULL),
+ fCommonHistsRes(NULL),
+ fEventNumber(0),
+ fQsum(NULL),
+ fQ2sum(0)
{
// Constructor.
fQsum = new TVector2(); // flow vector sum
fHistList = new TList();
+ fSecondRunList = new TList();
}
//destructor
delete fQsum;
delete fHistList;
+ delete fSecondRunList;
}
cout<<"---Analysis with Lee Yang Zeros Event Plane Method---"<<endl;
//input histograms
- if (fSecondRunFile->IsZombie()){ //check if file exists
- cout << "Error opening file, no input file from LYZ analysis" << endl;
- exit(-1);
- } else if (fSecondRunFile->IsOpen()){
- cout<<"----secondRunFile is open----"<<endl;
- //get List
- TList* list = (TList*)fSecondRunFile->Get("cobj1");
- fSecondVPt = (TProfile*)list->FindObject("Flow_Differential_Pt_LYZ"); //to compare to
- fSecondReDtheta = (TProfile*)list->FindObject("Second_FlowPro_ReDtheta_LYZ");
+ if (fSecondRunList) {
+ fSecondReDtheta = (TProfile*)fSecondRunList->FindObject("Second_FlowPro_ReDtheta_LYZ");
fHistList->Add(fSecondReDtheta);
- fSecondImDtheta = (TProfile*)list->FindObject("Second_FlowPro_ImDtheta_LYZ");
+
+ fSecondImDtheta = (TProfile*)fSecondRunList->FindObject("Second_FlowPro_ImDtheta_LYZ");
fHistList->Add(fSecondImDtheta);
- fFirstr0theta = (TProfile*)list->FindObject("First_FlowPro_r0theta_LYZ");
+
+ fFirstr0theta = (TProfile*)fSecondRunList->FindObject("First_FlowPro_r0theta_LYZ");
fHistList->Add(fFirstr0theta);
- }
+ //warnings
+ if (!fSecondReDtheta) {cout<<"fSecondReDtheta is NULL!"<<endl; }
+ if (!fSecondImDtheta) {cout<<"fSecondImDtheta is NULL!"<<endl; }
+ if (!fFirstr0theta) {cout<<"fFirstr0theta is NULL!"<<endl; }
- fCommonHists = new AliFlowCommonHist("LYZEP");
- fHistList->Add(fCommonHists->GetHistMultOrig());
- fHistList->Add(fCommonHists->GetHistMultInt());
- fHistList->Add(fCommonHists->GetHistMultDiff());
- fHistList->Add(fCommonHists->GetHistPtInt());
- fHistList->Add(fCommonHists->GetHistPtDiff());
- fHistList->Add(fCommonHists->GetHistPhiInt());
- fHistList->Add(fCommonHists->GetHistPhiDiff());
- fHistList->Add(fCommonHists->GetHistEtaInt());
- fHistList->Add(fCommonHists->GetHistEtaDiff());
- fHistList->Add(fCommonHists->GetHistProMeanPtperBin());
- fHistList->Add(fCommonHists->GetHistQ());
-
- fCommonHistsRes = new AliFlowCommonHistResults("LYZEP");
- fHistList->Add(fCommonHistsRes->GetHistDiffFlow());
- fHistList->Add(fCommonHistsRes->GetHistChi());
- fHistList->Add(fCommonHistsRes->GetHistIntFlow());
+ }
+ fCommonHists = new AliFlowCommonHist("AliFlowCommonHistLYZEP");
+ fHistList->Add(fCommonHists);
+ fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResultsLYZEP");
+ fHistList->Add(fCommonHistsRes);
+
Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
Double_t dPtMin = AliFlowCommonConstants::GetPtMin();
Double_t dPtMax = AliFlowCommonConstants::GetPtMax();
fHistProWr->SetYTitle("Wr");
fHistList->Add(fHistProWr);
+ fHistQsumforChi = new TH1F("Flow_QsumforChi_LYZEP","Flow_QsumforChi_LYZEP",3,-1.,2.);
+ fHistQsumforChi->SetXTitle("Qsum.X , Qsum.Y, Q2sum");
+ fHistQsumforChi->SetYTitle("value");
+ fHistList->Add(fHistQsumforChi);
+
fHistDeltaPhi = new TH1F("Flow_DeltaPhi_LYZEP","Flow_DeltaPhi_LYZEP",100,0.,3.14);
fHistDeltaPhi->SetXTitle("DeltaPhi");
fHistDeltaPhi->SetYTitle("Counts");
//get the Q vector from the FlowEvent
AliFlowVector vQ = anEvent->GetQ();
+ if (vQ.X()== 0. && vQ.Y()== 0. ) { cout<<"Q vector is NULL!"<<endl; }
//Weight with the multiplicity
Double_t dQX = 0.;
Double_t dQY = 0.;
dQY = vQ.Y()/vQ.GetMult();
} else {cerr<<"vQ.GetMult() is zero!"<<endl; }
vQ.Set(dQX,dQY);
+ //cout<<"vQ("<<dQX<<","<<dQY<<")"<<endl;
//for chi calculation:
*fQsum += vQ;
+ fHistQsumforChi->SetBinContent(1,fQsum->X());
+ fHistQsumforChi->SetBinContent(2,fQsum->Y());
fQ2sum += vQ.Mod2();
- cout<<"fQ2sum = "<<fQ2sum<<endl;
+ fHistQsumforChi->SetBinContent(3,fQ2sum);
+ //cout<<"fQ2sum = "<<fQ2sum<<endl;
//call AliFlowLYZEventPlane::CalculateRPandW() here!
aLYZEP->CalculateRPandW(vQ);
//plot histograms etc.
cout<<"AliFlowAnalysisWithLYZEventPlane::Terminate()"<<endl;
- //constands:
+
+ //constants:
Double_t dJ01 = 2.405;
Int_t iNtheta = AliFlowLYZConstants::kTheta;
Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
-
+ //set the event number
+ SetEventNumber((int)fCommonHists->GetHistMultOrig()->GetEntries());
+ //cout<<"number of events processed is "<<fEventNumber<<endl;
+
+ //set the sum of Q vectors
+ fQsum->Set(fHistQsumforChi->GetBinContent(1),fHistQsumforChi->GetBinContent(2));
+ SetQ2sum(fHistQsumforChi->GetBinContent(3));
+
//calculate dV the mean of dVtheta
Double_t dVtheta = 0;
Double_t dV = 0;
} //loop over b
-
- cout<<"Making some plots to check the results:"<<endl<<endl;
-
- TCanvas *canvas = new TCanvas("canvas","compare v2 vs pt",800,800);
- canvas->cd();
- fSecondVPt->SetLineColor(3);
- fSecondVPt->SetLineWidth(2);
- fSecondVPt->Draw();
- fHistFlow = fCommonHistsRes->GetHistDiffFlow();
- fHistFlow->Draw("SAME");
- // draw the legend
- TLegend *legend2 = new TLegend(0.6,0.65,0.88,0.85);
- legend2->SetTextFont(72);
- legend2->SetTextSize(0.04);
- legend2->AddEntry(fSecondVPt,"stand. LYZ","lpe");
- legend2->AddEntry(fHistFlow,"new LYZ with calculated errors","lpe");
- legend2->Draw();
-
- TCanvas *canvas5 = new TCanvas("canvas5","phi and Delta Phi",800,600);
- canvas5->Divide(3,1);
- canvas5->cd(1);
- fHistDeltaPhi->Draw();
- canvas5->cd(2);
- fHistPhiLYZ->Draw();
- canvas5->cd(3);
- fHistPhiEP->Draw();
-
cout<<".....finished"<<endl;
}
virtual void Make(AliFlowEventSimple* fEvent, AliFlowLYZEventPlane* fLYZEP);
virtual void Finish();
- // input files
- void SetSecondRunFileName(TString name) { this->fSecondRunFileName = name ; } // Sets input file name
- TString GetSecondRunFileName() const { return this->fSecondRunFileName ; } // Gets output file name
- void SetSecondRunFile(TFile* file) { this->fSecondRunFile = file ; } // Sets first run file
+ void SetEventNumber(Int_t n) { this->fEventNumber = n; }
+ Int_t GetEventNumber() const { return this->fEventNumber; }
+ void SetQ2sum(Double_t d) { this->fQ2sum = d; }
+ Double_t GetQ2sum() { return this->fQ2sum; }
//output
- TList* GetHistList() const {return this->fHistList; }
+ TList* GetHistList() const {return this->fHistList; }
+ AliFlowCommonHist* GetCommonHists() const { return this->fCommonHists; }
+ void SetCommonHists(AliFlowCommonHist* aCommonHist)
+ { this->fCommonHists = aCommonHist; }
+ AliFlowCommonHistResults* GetCommonHistsRes() const
+ { return this->fCommonHistsRes; }
+ void SetCommonHistsRes(AliFlowCommonHistResults* aCommonHistResult)
+ { this->fCommonHistsRes = aCommonHistResult; }
+
+ // !!!!! make getters and setters for all histograms
+ TProfile* GetSecondReDtheta() {return this->fSecondReDtheta; }
+ void SetSecondReDtheta(TProfile* aSecondReDtheta)
+ {this->fSecondReDtheta = aSecondReDtheta; }
+ TProfile* GetSecondImDtheta() {return this->fSecondImDtheta; }
+ void SetSecondImDtheta(TProfile* aSecondImDtheta)
+ {this->fSecondImDtheta = aSecondImDtheta; }
+ TProfile* GetFirstr0theta() {return this->fFirstr0theta; }
+ void SetFirstr0theta(TProfile* aFirstr0theta)
+ {this->fFirstr0theta = aFirstr0theta; }
+ TProfile* GetHistProFlow() {return this->fHistProFlow;}
+ void SetHistProFlow(TProfile* aHistProFlow)
+ {this->fHistProFlow =aHistProFlow; }
+ TProfile* GetHistProFlow2() {return this->fHistProFlow2;}
+ void SetHistProFlow2(TProfile* aHistProFlow2)
+ {this->fHistProFlow2 = aHistProFlow2; }
+ TProfile* GetHistProWr() {return this->fHistProWr; }
+ void SetHistProWr(TProfile* aHistProWr)
+ {this->fHistProWr = aHistProWr; }
+ TProfile* GetHistProWrCorr() {return this->fHistProWrCorr; }
+ void SetHistProWrCorr(TProfile* aHistProWrCorr)
+ {this->fHistProWrCorr = aHistProWrCorr; }
+ TH1F* GetHistQsumforChi() {return this->fHistQsumforChi; }
+ void SetHistQsumforChi(TH1F* aHistQsumforChi)
+ {this->fHistQsumforChi = aHistQsumforChi; }
+ TH1F* GetHistDeltaPhi() {return this->fHistDeltaPhi; }
+ void SetHistDeltaPhi(TH1F* aHistDeltaPhi)
+ {this->fHistDeltaPhi = aHistDeltaPhi; }
+ TH1F* GetHistDeltaPhi2() {return this->fHistDeltaPhi2; }
+ void SetHistDeltaPhi2(TH1F* aHistDeltaPhi2)
+ {this->fHistDeltaPhi2 = aHistDeltaPhi2; }
+ TH1F* GetHistDeltaPhihere() {return this->fHistDeltaPhihere; }
+ void SetHistDeltaPhihere(TH1F* aHistDeltaPhihere)
+ {this->fHistDeltaPhihere = aHistDeltaPhihere; }
+ TH1F* GetHistPhiEP() {return this->fHistPhiEP; }
+ void SetHistPhiEP(TH1F* aHistPhiEP)
+ {this->fHistPhiEP = aHistPhiEP; }
+ TH1F* GetHistPhiEPhere() {return this->fHistPhiEPhere; }
+ void SetHistPhiEPhere(TH1F* aHistPhiEPhere)
+ {this->fHistPhiEPhere = aHistPhiEPhere; }
+ TH1F* GetHistPhiLYZ() {return this->fHistPhiLYZ; }
+ void SetHistPhiLYZ(TH1F* aHistPhiLYZ)
+ {this->fHistPhiLYZ = aHistPhiLYZ; }
+ TH1F* GetHistPhiLYZ2() {return this->fHistPhiLYZ2;}
+ void SetHistPhiLYZ2(TH1F* aHistPhiLYZ2)
+ {this->fHistPhiLYZ2 = aHistPhiLYZ2; }
+ //input
+ void SetSecondRunList(TList* list) { this->fSecondRunList = list; }
+ TList* GetSecondRunList() { return this->fSecondRunList; }
private:
AliFlowAnalysisWithLYZEventPlane(const AliFlowAnalysisWithLYZEventPlane& aAnalysis); // copy constructor
AliFlowAnalysisWithLYZEventPlane& operator=(const AliFlowAnalysisWithLYZEventPlane& aAnalysis); // assignment operator
- TFile* fSecondRunFile ; // pointer to file from second run
- TString fSecondRunFileName;
-
//histograms
TList* fHistList; //list ro hold all histograms
+ TList* fSecondRunList; //list from Second LYZ run output
//input
TProfile* fSecondReDtheta; // input profile
TProfile* fSecondImDtheta; // input profile
TProfile* fFirstr0theta; // input profile
- TProfile* fSecondVPt; // input profile
//output
TProfile* fHistProFlow; //
TProfile* fHistProFlow2; //
TProfile* fHistProWr; //
TProfile* fHistProWrCorr; //
- TH1D* fHistFlow; //
+ TH1F* fHistQsumforChi; //
TH1F* fHistDeltaPhi; //
TH1F* fHistDeltaPhi2; //
TH1F* fHistDeltaPhihere; //
Double_t fQ2sum; // flow vector sum squared
- ClassDef(AliFlowAnalysisWithLYZEventPlane, 0); // lyz analysis
+ ClassDef(AliFlowAnalysisWithLYZEventPlane, 1); // lyz analysis
};
#endif
fDoubleLoop(kFALSE),
fDebug(kFALSE),
fHistList(NULL),
- firstRunFileName(0),
- firstRunFile(NULL),
+ fFirstRunList(NULL),
fHistProVtheta(NULL),
fHistProVeta(NULL),
fHistProVPt(NULL),
fHistProImDenom(NULL),
fHistProReDtheta(NULL),
fHistProImDtheta(NULL),
+ fHistQsumforChi(NULL),
fCommonHists(NULL),
fCommonHistsRes(NULL)
if (fDebug) cout<<"****AliFlowAnalysisWithLeeYangZeros::AliFlowAnalysisWithLeeYangZeros default constructor****"<<endl;
fHistList = new TList();
+ fFirstRunList = new TList();
for(Int_t i = 0;i<5;i++)
{
if (fDebug) cout<<"****~AliFlowAnalysisWithLeeYangZeros****"<<endl;
delete fQsum;
delete fHistList;
+ delete fFirstRunList;
+
}
//-----------------------------------------------------------------------
Double_t dPtMax = AliFlowCommonConstants::GetPtMax();
Double_t dEtaMin = AliFlowCommonConstants::GetEtaMin();
Double_t dEtaMax = AliFlowCommonConstants::GetEtaMax();
+
+ //for control histograms
+ fCommonHists = new AliFlowCommonHist("AliFlowCommonHistLYZ");
+ fHistList->Add(fCommonHists);
+ fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResultsLYZ");
+ fHistList->Add(fCommonHistsRes);
- //for control histograms
- fCommonHists = new AliFlowCommonHist("LYZ");
- //fHistList->Add(fCommonHists->GetHistList());
- fHistList->Add(fCommonHists->GetHistMultOrig());
- fHistList->Add(fCommonHists->GetHistMultInt());
- fHistList->Add(fCommonHists->GetHistMultDiff());
- fHistList->Add(fCommonHists->GetHistPtInt());
- fHistList->Add(fCommonHists->GetHistPtDiff());
- fHistList->Add(fCommonHists->GetHistPhiInt());
- fHistList->Add(fCommonHists->GetHistPhiDiff());
- fHistList->Add(fCommonHists->GetHistEtaInt());
- fHistList->Add(fCommonHists->GetHistEtaDiff());
- fHistList->Add(fCommonHists->GetHistProMeanPtperBin());
- fHistList->Add(fCommonHists->GetHistQ());
- fCommonHistsRes = new AliFlowCommonHistResults("LYZ");
- //fHistList->Add(fCommonHistsRes->GetHistList());
- fHistList->Add(fCommonHistsRes->GetHistDiffFlow());
- fHistList->Add(fCommonHistsRes->GetHistChi());
- fHistList->Add(fCommonHistsRes->GetHistIntFlow());
+ fHistQsumforChi = new TH1F("Flow_QsumforChi_LYZEP","Flow_QsumforChi_LYZEP",3,-1.,2.);
+ fHistQsumforChi->SetXTitle("Qsum.X , Qsum.Y, Q2sum");
+ fHistQsumforChi->SetYTitle("value");
+ fHistList->Add(fHistQsumforChi);
//for first loop over events
if (fFirstRun){
//class AliFlowLYZHist1 defines the histograms: fHistProGtheta, fHistProReGtheta, fHistProImGtheta, fHistProR0theta
for (Int_t theta=0;theta<iNtheta;theta++) {
- fHist1[theta]=new AliFlowLYZHist1(theta);
- //fHistList->Add(fHist1[theta]->GetHistList() );
- fHistList->Add(fHist1[theta]->GetHistGtheta() );
- fHistList->Add(fHist1[theta]->GetHistProReGtheta() );
- fHistList->Add(fHist1[theta]->GetHistProImGtheta() );
+ TString name = "AliFlowLYZHist1_";
+ name += theta;
+ fHist1[theta]=new AliFlowLYZHist1(theta, name);
+ fHistList->Add(fHist1[theta]);
}
-
+
}
//for second loop over events
else {
//class AliFlowLYZHist2 defines the histograms:
for (Int_t theta=0;theta<iNtheta;theta++) {
- fHist2[theta]=new AliFlowLYZHist2(theta);
- //fHistList->Add(fHist2[theta]->GetHistList() );
- fHistList->Add(fHist2[theta]->GetHistProReNumer() );
- fHistList->Add(fHist2[theta]->GetHistProImNumer() );
- fHistList->Add(fHist2[theta]->GetHistProReNumerPt() );
- fHistList->Add(fHist2[theta]->GetHistProImNumerPt() );
- //fHistList->Add(fHist2[theta]->GetHistProReNumer2D() ); //gives error in compilation
- //fHistList->Add(fHist2[theta]->GetHistProImNumer2D() ); //gives error in compilation
+ TString name = "AliFlowLYZHist2_";
+ name += theta;
+ fHist2[theta]=new AliFlowLYZHist2(theta,name);
+ fHistList->Add(fHist2[theta]);
}
-
- //read hists from first run file
- //firstRunFile = new TFile("fof_flowLYZAnal_firstrun.root","READ"); //default is read
- if (firstRunFile->IsZombie()){ //check if file exists
- cout << "Error opening file, run first with fFirstrun = kTRUE" << endl;
- exit(-1);
- } else if (firstRunFile->IsOpen()){
- cout<<"----firstRunFile is open----"<<endl<<endl;
- TList* list = (TList*)firstRunFile->Get("cobj1");
- if (!list) {cout<<"list is NULL pointer!"<<endl;}
- fHistProR0theta = (TProfile*)list->FindObject("First_FlowPro_r0theta_LYZ");
+
+ //read histogram fHistProR0theta from the first run list
+ if (fFirstRunList) {
+ fHistProR0theta = (TProfile*)fFirstRunList->FindObject("First_FlowPro_r0theta_LYZ");
if (!fHistProR0theta) {cout<<"fHistProR0theta has a NULL pointer!"<<endl;}
fHistList->Add(fHistProR0theta);
- }
+ } else { cout<<"list is NULL pointer!"<<endl; }
+
}
//define variables for both runs
Double_t dJ01 = 2.405;
Int_t iNtheta = AliFlowLYZConstants::kTheta;
-
+ //set the event number
+ SetEventNumber((int)fCommonHists->GetHistMultOrig()->GetEntries());
+ //cout<<"number of events processed is "<<fEventNumber<<endl;
+
+ //set the sum of Q vectors
+ fQsum->Set(fHistQsumforChi->GetBinContent(1),fHistQsumforChi->GetBinContent(2));
+ SetQ2sum(fHistQsumforChi->GetBinContent(3));
+
if (fFirstRun){
Double_t dR0 = 0;
Double_t dVtheta = 0;
// recalculate statistical errors on integrated flow
//combining 5 theta angles to 1 relative error BP eq. 89
Double_t dRelErr2comb = 0.;
- Int_t iEvts = fEventNumber;
+ Int_t iEvts = fEventNumber;
for (Int_t theta=0;theta<iNtheta;theta++){
Double_t dTheta = ((double)theta/iNtheta)*TMath::Pi();
Double_t dApluscomb = TMath::Exp((dJ01*dJ01)/(2*dChi*dChi)*
}
dRelErr2comb /= iNtheta;
Double_t dRelErrcomb = TMath::Sqrt(dRelErr2comb);
- cout<<"dRelErrcomb = "<<dRelErrcomb<<endl;
//copy content of profile into TH1D and add error
Double_t dv2pro = dV; //in the case that fv is equal to fV
Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
Double_t dEta, dPt, dReRatio, dVeta, dVPt;
-
-
+
Double_t dR0 = 0.;
Double_t dVtheta = 0.;
Double_t dV = 0.;
//} //loop over theta
cDenom(dReDenom,dImDenom);
-
-
+
//for new method and use by others (only with the sum generating function):
if (fUseSum) {
dR0 = fHistProR0theta->GetBinContent(theta+1);
//fill TH1D
fCommonHistsRes->FillDifferentialFlow(b, dv2pro, dErrdifcomb);
} //loop over bins b
-
-
- //close the first run file
- //firstRunFile->Close(); //gives error when writing to TList
-
-
+
} //secondrun
cout<<"----LYZ analysis finished....----"<<endl<<endl;
dQY = vQ.Y()/vQ.GetMult();
}
vQ.Set(dQX,dQY);
+
//for chi calculation:
*fQsum += vQ;
+ fHistQsumforChi->SetBinContent(1,fQsum->X());
+ fHistQsumforChi->SetBinContent(2,fQsum->Y());
fQ2sum += vQ.Mod2();
+ fHistQsumforChi->SetBinContent(3,fQ2sum);
//cerr<<"fQ2sum = "<<fQ2sum<<endl;
for (Int_t theta=0;theta<iNtheta;theta++)
dQX = vQ.X()/vQ.GetMult();
dQY = vQ.Y()/vQ.GetMult();
}
- vQ.Set(dQX,dQY);
+ vQ.Set(dQX,dQY);
+
//for chi calculation:
*fQsum += vQ;
+ fHistQsumforChi->SetBinContent(1,fQsum->X());
+ fHistQsumforChi->SetBinContent(2,fQsum->Y());
fQ2sum += vQ.Mod2();
+ fHistQsumforChi->SetBinContent(3,fQ2sum);
for (Int_t theta=0;theta<iNtheta;theta++)
{
Bool_t Init(); //defines variables and histograms
Bool_t Make(AliFlowEventSimple* anEvent); //calculates variables and fills histograms
Bool_t Finish(); //saves histograms
-
+
Double_t GetQtheta(AliFlowVector aQ, Double_t aTheta);
- void SetFirstRun(Bool_t kt) { this->fFirstRun = kt ; }
- Bool_t GetFirstRun() const { return this->fFirstRun ; }
- void SetUseSum(Bool_t kt) { this->fUseSum = kt ; }
- Bool_t GetUseSum() const { return this->fUseSum ; }
- void SetDoubleLoop(Bool_t kt) { this->fDoubleLoop = kt ; }
- Bool_t GetDoubleLoop() const { return this->fDoubleLoop ; }
- void SetDebug(Bool_t kt) { this->fDebug = kt ; }
- Bool_t GetDebug() const { return this->fDebug ; }
-
+ void SetFirstRun(Bool_t kt) { this->fFirstRun = kt ; }
+ Bool_t GetFirstRun() const { return this->fFirstRun ; }
+ void SetUseSum(Bool_t kt) { this->fUseSum = kt ; }
+ Bool_t GetUseSum() const { return this->fUseSum ; }
+ void SetDoubleLoop(Bool_t kt) { this->fDoubleLoop = kt ; }
+ Bool_t GetDoubleLoop() const { return this->fDoubleLoop ; }
+ void SetDebug(Bool_t kt) { this->fDebug = kt ; }
+ Bool_t GetDebug() const { return this->fDebug ; }
+
+ void SetEventNumber(Int_t n) { this->fEventNumber = n; }
+ Int_t GetEventNumber() const { return this->fEventNumber; }
+ void SetQ2sum(Double_t d) { this->fQ2sum = d; }
+ Double_t GetQ2sum() { return this->fQ2sum; }
// Output
- TList* GetHistList() const { return this->fHistList ; } // Gets output histogram list
-
- // input for second run
- void SetFirstRunFileName(TString name) { this->firstRunFileName = name ; } // Sets input file name
- TString GetFirstRunFileName() const { return this->firstRunFileName ; } // Gets output file name
- void SetFirstRunFile(TFile* file) { this->firstRunFile = file ; } // Sets first run file
+ TList* GetHistList() const { return this->fHistList ; }
+ AliFlowCommonHist* GetCommonHists() const { return this->fCommonHists; }
+ void SetCommonHists(AliFlowCommonHist* aCommonHist)
+ { this->fCommonHists = aCommonHist; }
+ AliFlowCommonHistResults* GetCommonHistsRes() const
+ { return this->fCommonHistsRes; }
+ void SetCommonHistsRes(AliFlowCommonHistResults* aCommonHistResult)
+ { this->fCommonHistsRes = aCommonHistResult; }
+ //AliFlowLYZHist1* GetHist1() const {return this->fHist1; }
+ void SetHist1(AliFlowLYZHist1* aLYZHist1[])
+ {for (Int_t i=0;i<5;i++) {this->fHist1[i] = aLYZHist1[i];} }
+ //AliFlowLYZHist2* GetHist2() const {return this->fHist2; }
+ void SetHist2(AliFlowLYZHist2* aLYZHist2[])
+ {for (Int_t i=0;i<5;i++) {this->fHist2[i] = aLYZHist2[i];} }
+
+ TProfile* GetHistProVtheta() const {return this->fHistProVtheta; }
+ void SetHistProVtheta(TProfile* aHistProVtheta)
+ { this->fHistProVtheta = aHistProVtheta; }
+ TProfile* GetHistProVeta() const {return this->fHistProVeta; }
+ void SetHistProVeta(TProfile* aHistProVeta)
+ {this->fHistProVeta = aHistProVeta; }
+ TProfile* GetHistProVPt() const {return this->fHistProVPt;}
+ void SetHistProVPt(TProfile* aHistProVPt)
+ {this->fHistProVPt = aHistProVPt; }
+ TProfile* GetHistProR0theta() const {return this->fHistProR0theta; }
+ void SetHistProR0theta(TProfile* aHistProR0theta)
+ {this->fHistProR0theta = aHistProR0theta; }
+ TProfile* GetHistProReDenom() const {return this->fHistProReDenom; }
+ void SetHistProReDenom(TProfile* aHistProReDenom)
+ {this->fHistProReDenom = aHistProReDenom; }
+ TProfile* GetHistProImDenom() const {return this->fHistProImDenom; }
+ void SetHistProImDenom(TProfile* aHistProImDenom)
+ {this->fHistProImDenom = aHistProImDenom; }
+ TProfile* GetHistProReDtheta() const {return this->fHistProReDtheta; }
+ void SetHistProReDtheta(TProfile* aHistProReDtheta)
+ {this->fHistProReDtheta = aHistProReDtheta; }
+ TProfile* GetHistProImDtheta() const {return this->fHistProImDtheta; }
+ void SetHistProImDtheta(TProfile* aHistProImDtheta)
+ {this->fHistProImDtheta = aHistProImDtheta; }
+ TH1F* GetHistQsumforChi() {return this->fHistQsumforChi; }
+ void SetHistQsumforChi(TH1F* aHistQsumforChi)
+ {this->fHistQsumforChi = aHistQsumforChi; }
+
+ void SetFirstRunList(TList* list) { this->fFirstRunList = list; }
+ TList* GetFirstRunList() { return this->fFirstRunList; }
private:
- AliFlowAnalysisWithLeeYangZeros(const AliFlowAnalysisWithLeeYangZeros& aAnalysis); // copy constructor
+ AliFlowAnalysisWithLeeYangZeros(const AliFlowAnalysisWithLeeYangZeros& aAnalysis); // copy constructor
AliFlowAnalysisWithLeeYangZeros& operator=(const AliFlowAnalysisWithLeeYangZeros& aAnalysis); //assignment operator
Bool_t MakeControlHistograms(AliFlowEventSimple* anEvent);
TComplex GetGrtheta(AliFlowEventSimple* anEvent, Double_t aR, Double_t aTheta);
TComplex GetDiffFlow(AliFlowEventSimple* anEvent, Double_t aR, Int_t theta);
- Double_t GetR0(TH1D* fHistGtheta);
-
-
+ Double_t GetR0(TH1D* fHistGtheta);
#ifndef __CINT__
Bool_t fDoubleLoop ; // flag for studying non flow effects
Bool_t fDebug ; // flag for lyz analysis: more print statements
- TList* fHistList; //list to hold all output histograms
- TString firstRunFileName; //
- TFile* firstRunFile; //
-
-
- TProfile* fHistProVtheta; //
- TProfile* fHistProVeta; //
- TProfile* fHistProVPt; //
- TProfile* fHistProR0theta; //
- TProfile* fHistProReDenom; //
- TProfile* fHistProImDenom; //
- TProfile* fHistProReDtheta; //
- TProfile* fHistProImDtheta; //
-
-
+ TList* fHistList; //list to hold all output histograms
+ TList* fFirstRunList; //list from first run output
+
+ TProfile* fHistProVtheta;
+ TProfile* fHistProVeta;
+ TProfile* fHistProVPt;
+ TProfile* fHistProR0theta;
+ TProfile* fHistProReDenom;
+ TProfile* fHistProImDenom;
+ TProfile* fHistProReDtheta;
+ TProfile* fHistProImDtheta;
+ TH1F* fHistQsumforChi; //
+
+
//class AliFlowLYZHist1 defines the histograms: fHistProGtheta, fHistProReGtheta, fHistProImGtheta
AliFlowLYZHist1* fHist1[5]; //
AliFlowCommonHist* fCommonHists; //
AliFlowCommonHistResults* fCommonHistsRes; //
- ClassDef(AliFlowAnalysisWithLeeYangZeros,0) // macro for rootcint
+ ClassDef(AliFlowAnalysisWithLeeYangZeros,1) // macro for rootcint
};
Double_t dPtMin = AliFlowCommonConstants::GetPtMin();
Double_t dPtMax = AliFlowCommonConstants::GetPtMax();
- fCommonHists = new AliFlowCommonHist("MC");
- fHistList->Add(fCommonHists->GetHistMultOrig());
- fHistList->Add(fCommonHists->GetHistMultInt());
- fHistList->Add(fCommonHists->GetHistMultDiff());
- fHistList->Add(fCommonHists->GetHistPtInt());
- fHistList->Add(fCommonHists->GetHistPtDiff());
- fHistList->Add(fCommonHists->GetHistPhiInt());
- fHistList->Add(fCommonHists->GetHistPhiDiff());
- fHistList->Add(fCommonHists->GetHistEtaInt());
- fHistList->Add(fCommonHists->GetHistEtaDiff());
- fHistList->Add(fCommonHists->GetHistProMeanPtperBin());
- fHistList->Add(fCommonHists->GetHistQ());
- fCommonHistsRes = new AliFlowCommonHistResults("MC");
- fHistList->Add(fCommonHistsRes->GetHistDiffFlow());
- fHistList->Add(fCommonHistsRes->GetHistChi());
- fHistList->Add(fCommonHistsRes->GetHistIntFlow());
-
- fHistProFlow = new TProfile("FlowPro_VPt_MC","FlowPro_VPt_MC",iNbinsPt,dPtMin,dPtMax);
+ fCommonHists = new AliFlowCommonHist("AliFlowCommonHistMCEP");
+ fHistList->Add(fCommonHists);
+ fCommonHistsRes = new AliFlowCommonHistResults("AliFlowCommonHistResultsMCEP");
+ fHistList->Add(fCommonHistsRes);
+
+ fHistProFlow = new TProfile("FlowPro_VPt_MCEP","FlowPro_VPt_MCEP",iNbinsPt,dPtMin,dPtMax);
fHistProFlow->SetXTitle("Pt");
fHistProFlow->SetYTitle("v2 (%)");
fHistList->Add(fHistProFlow);
- fHistRP = new TH1F("Flow_RP_MC","Flow_RP_MC",100,0.,3.14);
+ fHistRP = new TH1F("Flow_RP_MCEP","Flow_RP_MCEP",100,0.,3.14);
fHistRP->SetXTitle("Reaction Plane Angle");
fHistRP->SetYTitle("Counts");
fHistList->Add(fHistRP);
*fQsum += vQ;
//cout<<"fQsum.Mod() = "<<fQsum.Mod()<<endl;
fQ2sum += vQ.Mod2();
- cout<<"fQ2sum = "<<fQ2sum<<endl;
+ //cout<<"fQ2sum = "<<fQ2sum<<endl;
fHistRP->Fill(aRP);
if (fDebug) cout<<"AliFlowAnalysisWithMCEventPlane::Terminate()"<<endl;
Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
-
+
+
TH1F* fHistPtDiff = fCommonHists->GetHistPtDiff();
Double_t dV = 0.;
Double_t dErrV = 0.;
Double_t dSum = 0.;
for(Int_t b=0;b<iNbinsPt;b++){
Double_t dv2pro = 0.;
- Double_t dErrdifcomb = 0.; //in case error from profile is correct
+ Double_t dErrdifcomb = 0.;
if(fHistProFlow) {
dv2pro = fHistProFlow->GetBinContent(b);
- dErrdifcomb = fHistProFlow->GetBinError(b);
+ dErrdifcomb = fHistProFlow->GetBinError(b); //in case error from profile is correct
//fill TH1D
fCommonHistsRes->FillDifferentialFlow(b, dv2pro, dErrdifcomb);
if (fHistPtDiff){
//error on integrated flow
dErrV += dYield*dYield*(dErrdifcomb/100)*(dErrdifcomb/100);
}
- }
+ } else { cout<<"fHistProFlow is NULL"<<endl; }
}
- dV /= dSum; //because pt distribution should be normalised
- dErrV /= dSum*dSum;
- dErrV = TMath::Sqrt(dErrV);
+ if (dSum != 0. ) {
+ dV /= dSum; //because pt distribution should be normalised
+ dErrV /= dSum*dSum;
+ dErrV = TMath::Sqrt(dErrV); }
cout<<"dV is "<<dV<<" +- "<<dErrV<<endl;
fCommonHistsRes->FillIntegratedFlow(dV,dErrV);
cout<<".....finished"<<endl;
- }
+}
+
AliFlowAnalysisWithMCEventPlane(); //default constructor
virtual ~AliFlowAnalysisWithMCEventPlane(); //destructor
- void Init(); //defines variables and histograms
- void Make(AliFlowEventSimple* anEvent, Double_t aRP); //calculates variables and fills histograms
- void Finish(); //saves histograms
-
- void SetDebug(Bool_t kt) { this->fDebug = kt ; }
- Bool_t GetDebug() const { return this->fDebug ; }
+ void Init(); //defines variables and histograms
+ void Make(AliFlowEventSimple* anEvent, Double_t aRP); //calculates variables and fills histograms
+ void Finish(); //saves histograms
+
+ void SetDebug(Bool_t kt) { this->fDebug = kt ; }
+ Bool_t GetDebug() const { return this->fDebug ; }
+ void SetEventNumber(Int_t n) { this->fEventNumber = n; }
+ Int_t GetEventNumber() const { return this->fEventNumber; }
// Output
- TList* GetHistList() const { return this->fHistList ; }
-
+ TList* GetHistList() const { return this->fHistList ; }
+ AliFlowCommonHist* GetCommonHists() const { return this->fCommonHists; }
+ void SetCommonHists(AliFlowCommonHist* aCommonHist)
+ { this->fCommonHists = aCommonHist; }
+ AliFlowCommonHistResults* GetCommonHistsRes() const
+ { return this->fCommonHistsRes; }
+ void SetCommonHistsRes(AliFlowCommonHistResults* aCommonHistResult)
+ { this->fCommonHistsRes = aCommonHistResult; }
+ //histograms
+ TProfile* GetHistProFlow() {return this->fHistProFlow; }
+ void SetHistProFlow(TProfile* aHistProFlow)
+ {this->fHistProFlow = aHistProFlow; }
+ TH1F* GetHistRP() {return this->fHistRP; }
+ void SetHistRP(TH1F* aHistRP) {this->fHistRP = aHistRP; }
+
private:
AliFlowAnalysisWithMCEventPlane(const AliFlowAnalysisWithMCEventPlane& aAnalysis); //copy constructor
TProfile* fHistProFlow; //
TH1F* fHistRP; //
- ClassDef(AliFlowAnalysisWithMCEventPlane,0) // macro for rootcint
+ ClassDef(AliFlowAnalysisWithMCEventPlane,1) // macro for rootcint
};
fHistList->Add(fHistProUQ);
fCommonHists = new AliFlowCommonHist("AliFlowCommonHistSP");
- // fHistList->Add(fCommonHists->GetHistList());
-
fHistList->Add(fCommonHists);
- // commented for test writing full object
- // fHistList->Add(fCommonHists->GetHistMultOrig());
- // fHistList->Add(fCommonHists->GetHistMultInt());
- // fHistList->Add(fCommonHists->GetHistMultDiff());
- // fHistList->Add(fCommonHists->GetHistPtInt());
- // fHistList->Add(fCommonHists->GetHistPtDiff());
- // fHistList->Add(fCommonHists->GetHistPhiInt());
- // fHistList->Add(fCommonHists->GetHistPhiDiff());
- // fHistList->Add(fCommonHists->GetHistEtaInt());
- // fHistList->Add(fCommonHists->GetHistEtaDiff());
- // fHistList->Add(fCommonHists->GetHistProMeanPtperBin());
- // fHistList->Add(fCommonHists->GetHistQ());
- // end test
-
+
//fCommonHistsRes = new AliFlowCommonHistResults("SP");
fEventNumber = 0; //set number of events to zero
AliFlowVector vQ = anEvent->GetQ();
//weight by the Multiplicity
- Double_t dQX = vQ.X()/vQ.GetMult();
- Double_t dQY = vQ.Y()/vQ.GetMult();
+ Double_t dQX = 0.;
+ Double_t dQY = 0.;
+ if (iNumberOfTracks!=0.) {
+ dQX = vQ.X()/vQ.GetMult();
+ dQY = vQ.Y()/vQ.GetMult();
+ }
vQ.Set(dQX,dQY);
fHistQ->Fill(vQ.Mod());
#include "TH1D.h" //needed as include
#include "TMath.h" //needed as include
#include "TList.h"
+#include "TBrowser.h"
class TH1F;
class AliFlowVector;
//-----------------------------------------------------------------------
AliFlowCommonHistResults::AliFlowCommonHistResults():
- TObject(),
+ TNamed(),
fHistIntFlow(0),
fHistDiffFlow(0),
fHistChi(0),
//-----------------------------------------------------------------------
- AliFlowCommonHistResults::AliFlowCommonHistResults(TString input):
- TObject(),
+ AliFlowCommonHistResults::AliFlowCommonHistResults(const char *anInput,const char *title):
+ TNamed(anInput,title),
fHistIntFlow(0),
fHistDiffFlow(0),
fHistChi(0),
//integrated flow
name = "Flow_Integrated_";
- name +=input;
+ name += anInput;
fHistIntFlow = new TH1D(name.Data(), name.Data(),1,0.5,1.5);
fHistIntFlow ->SetXTitle("");
fHistIntFlow ->SetYTitle("Integrated Flow value (%)");
//differential flow
name = "Flow_Differential_Pt_";
- name +=input;
+ name += anInput;
fHistDiffFlow = new TH1D(name.Data(), name.Data(),iNbinsPt,dPtMin,dPtMax);
fHistDiffFlow ->SetXTitle("Pt");
fHistDiffFlow ->SetYTitle("v (%)");
//Chi (needed for rebinning later on)
name = "Flow_Chi_";
- name +=input;
+ name += anInput;
fHistChi = new TH1D(name.Data(), name.Data(),1,0.5,1.5);
fHistChi ->SetXTitle("");
fHistChi ->SetYTitle("Chi");
}
+//-----------------------------------------------------------------------
+void AliFlowCommonHistResults::Print(Option_t *option) const
+{
+ // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
+ // ===============================================
+ // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
+ printf( "Class.Print Name = %s, Histogram list:\n",GetName());
+
+ if (fHistList) {
+ fHistList->Print(option);
+ }
+ else
+ {
+ printf( "Empty histogram list \n");
+ }
+}
+
+//-----------------------------------------------------------------------
+ void AliFlowCommonHistResults::Browse(TBrowser *b)
+{
+
+ if (!b) return;
+ if (fHistList) b->Add(fHistList,"AliFlowCommonHistResultsList");
+}
+
+
+
class TH1D;
class TCollection;
class TList;
+class TBrowser;
// AliFlowCommonHistResults:
// Class to organize the common histograms for Flow Analysis
// authors: N. van der Kolk (kolk@nikhef.nl) and A. Bilandzic (anteb@nikhef.nl)
-class AliFlowCommonHistResults : public TObject {
+class AliFlowCommonHistResults : public TNamed {
public:
AliFlowCommonHistResults(); //default constructor
- AliFlowCommonHistResults(TString input); //constructor
+ AliFlowCommonHistResults(const char *name,const char *title = "AliFlowCommonHist"); //constructor
virtual ~AliFlowCommonHistResults(); //destructor
+ Bool_t IsFolder() const {return kTRUE;};
+ void Browse(TBrowser *b);
+
//make fill methods here
Bool_t FillIntegratedFlow(Double_t aV, Double_t anError); //fill fHistIntFlow
Bool_t FillDifferentialFlow(Int_t aBin, Double_t av, Double_t anError); //fill fHistDiffFlow
TList* GetHistList() {return fHistList;} ;
virtual Double_t Merge(TCollection *aList); //merge function
+ void Print(Option_t* option = "") const; //method to print stats
private:
TH1D* fHistChi; //resolution
TList* fHistList; //list to hold all histograms
- ClassDef(AliFlowCommonHistResults,0) // macro for rootcint
+ ClassDef(AliFlowCommonHistResults,1) // macro for rootcint
};
#endif
#include "TProfile.h"
#include "TFile.h"
#include "TComplex.h"
+#include "TList.h"
#include "AliFlowVector.h"
#include "AliFlowLYZConstants.h"
//-----------------------------------------------------------------------
AliFlowLYZEventPlane::AliFlowLYZEventPlane():
- fSecondRunFile(0),
- fSecondRunFileName("noname.ESD"),
+ fSecondRunList(0),
fWR(0),
fPsi(0),
fSecondReDtheta(0),
AliFlowLYZEventPlane::~AliFlowLYZEventPlane()
{
//destructor
-
+ delete fSecondRunList;
}
cout<<"---Lee Yang Zeros Event Plane Method---"<<endl;
//input histograms
- if (fSecondRunFile->IsZombie()){ //check if file exists
- cout << "Error opening file, run first regular LYZ second run" << endl;
- exit(-1);
- } else if (fSecondRunFile->IsOpen()){
- cout<<"----secondRunFile is open----"<<endl;
- //get List
- TList* list = (TList*)fSecondRunFile->Get("cobj1");
- fSecondReDtheta = ( TProfile*)list->FindObject("Second_FlowPro_ReDtheta_LYZ");
- fSecondImDtheta = ( TProfile*)list->FindObject("Second_FlowPro_ImDtheta_LYZ");
- fFirstr0theta = ( TProfile*)list->FindObject("First_FlowPro_r0theta_LYZ");
+ if (fSecondRunList) {
+ fSecondReDtheta = (TProfile*)fSecondRunList->FindObject("Second_FlowPro_ReDtheta_LYZ");
+ fSecondImDtheta = (TProfile*)fSecondRunList->FindObject("Second_FlowPro_ImDtheta_LYZ");
+ fFirstr0theta = (TProfile*)fSecondRunList->FindObject("First_FlowPro_r0theta_LYZ");
+
+ //warnings
+ if (!fSecondReDtheta) {cout<<"fSecondReDtheta is NULL!"<<endl; }
+ if (!fSecondImDtheta) {cout<<"fSecondImDtheta is NULL!"<<endl; }
+ if (!fFirstr0theta) {cout<<"fFirstr0theta is NULL!"<<endl; }
}
+
}
//-----------------------------------------------------------------------
Double_t dSinTerm = 0;
TComplex cDtheta;
TComplex cRatio;
-
- for (Int_t theta=0;theta<iNtheta;theta++) {
- Double_t dTheta = ((float)theta/iNtheta)*TMath::Pi()/2;
- //Calculate Qtheta
- Double_t dQtheta = aQ.X()*cos(2*dTheta)+aQ.Y()*sin(2*dTheta); //get Qtheta from Q vector
+
+ if (aQ.Mod()==0.) { cout<<"Q vector is NULL"<<endl; }
+ else {
+ for (Int_t theta=0;theta<iNtheta;theta++) {
+ Double_t dTheta = ((float)theta/iNtheta)*TMath::Pi()/2;
+ //Calculate Qtheta
+ Double_t dQtheta = aQ.X()*cos(2*dTheta)+aQ.Y()*sin(2*dTheta); //get Qtheta from Q vector
- //get R0
- Double_t dR0 = fFirstr0theta->GetBinContent(theta+1);
- //cerr<<"dR0 = "<<dR0<<endl;
-
- //get Dtheta
- Double_t dReDtheta = fSecondReDtheta->GetBinContent(theta+1);
- //cerr<<"dReDtheta = "<<dReDtheta<<endl;
- Double_t dImDtheta = fSecondImDtheta->GetBinContent(theta+1);
- //cerr<<"dImDtheta = "<<dImDtheta<<endl;
- cDtheta(dReDtheta,dImDtheta);
+ //get R0
+ Double_t dR0 = fFirstr0theta->GetBinContent(theta+1);
+ //cerr<<"dR0 = "<<dR0<<endl;
+
+ //get Dtheta
+ Double_t dReDtheta = fSecondReDtheta->GetBinContent(theta+1);
+ //cerr<<"dReDtheta = "<<dReDtheta<<endl;
+ Double_t dImDtheta = fSecondImDtheta->GetBinContent(theta+1);
+ //cerr<<"dImDtheta = "<<dImDtheta<<endl;
+ cDtheta(dReDtheta,dImDtheta);
- TComplex cExpo(0.,dR0*dQtheta); //Complex number: 0 +(i r0 Qtheta)
- if (cDtheta.Rho()!=0.) { cRatio =(TComplex::Exp(cExpo))/cDtheta; } //(e^(i r0 Qtheta))/Dtheta
- else { cRatio(0.,0.); }
+ TComplex cExpo(0.,dR0*dQtheta); //Complex number: 0 +(i r0 Qtheta)
+ if (cDtheta.Rho()!=0.) { cRatio =(TComplex::Exp(cExpo))/cDtheta; } //(e^(i r0 Qtheta))/Dtheta
+ else { cRatio(0.,0.); }
- //sum over theta
- dCosTerm += cRatio.Re() * TMath::Cos(2*dTheta); //Re{(e^(i r0 Qtheta))/Dtheta } cos(2 theta)
- dSinTerm += cRatio.Re() * TMath::Sin(2*dTheta); //Re{(e^(i r0 Qtheta))/Dtheta } sin(2 theta)
+ //sum over theta
+ dCosTerm += cRatio.Re() * TMath::Cos(2*dTheta); //Re{(e^(i r0 Qtheta))/Dtheta } cos(2 theta)
+ dSinTerm += cRatio.Re() * TMath::Sin(2*dTheta); //Re{(e^(i r0 Qtheta))/Dtheta } sin(2 theta)
- }//loop over theta
+ }//loop over theta
- //average over theta
- dCosTerm /= iNtheta;
- dSinTerm /= iNtheta;
- //cerr<<"cosTerm & sinTerm are: "<<dCosTerm<<" & "<<dSinTerm<<endl;
+ //average over theta
+ dCosTerm /= iNtheta;
+ dSinTerm /= iNtheta;
+ //cerr<<"cosTerm & sinTerm are: "<<dCosTerm<<" & "<<dSinTerm<<endl;
- //calculate fWR
- fWR = TMath::Sqrt(dCosTerm*dCosTerm + dSinTerm*dSinTerm);
+ //calculate fWR
+ fWR = TMath::Sqrt(dCosTerm*dCosTerm + dSinTerm*dSinTerm);
- //calculate fRP
- fPsi = 0.5*TMath::ATan2(dSinTerm,dCosTerm); //takes care of the signs correctly!
- if (fPsi < 0.) { fPsi += TMath::Pi(); } //to shift distribution from (-pi/2 to pi/2) to (0 to pi)
-
+ //calculate fRP
+ fPsi = 0.5*TMath::ATan2(dSinTerm,dCosTerm); //takes care of the signs correctly!
+ if (fPsi < 0.) { fPsi += TMath::Pi(); } //to shift distribution from (-pi/2 to pi/2) to (0 to pi)
+ }
+
}
class AliFlowEventSimple;
class TProfile;
class TFile;
+class TList;
class AliFlowLYZEventPlane {
public:
Double_t GetWR() const {return this->fWR; }
Double_t GetPsi() const {return this->fPsi; }
- // input file
- void SetSecondRunFileName(TString name)
- { this->fSecondRunFileName = name ; } // Sets input file name
- TString GetSecondRunFileName() const
- { return this->fSecondRunFileName ; } // Gets output file name
- void SetSecondRunFile(TFile* file)
- { this->fSecondRunFile = file ; } // Sets second run file
-
-
+ //input
+ void SetSecondRunList(TList* list) { this->fSecondRunList = list; }
+ TList* GetSecondRunList() { return this->fSecondRunList; }
+
private:
AliFlowLYZEventPlane(const AliFlowLYZEventPlane& aAnalysis); // copy constructor
AliFlowLYZEventPlane& operator=(const AliFlowLYZEventPlane& aAnalysis); // assignment operator
- TFile* fSecondRunFile ; //! pointer to file from second run
- TString fSecondRunFileName; //!
-
- Double_t fWR; // event weight
- Double_t fPsi; // reaction plane
+ TList* fSecondRunList; // list from Second LYZ run output
+ Double_t fWR; // event weight
+ Double_t fPsi; // reaction plane
TProfile* fSecondReDtheta; // holds Re of Dtheta
TProfile* fSecondImDtheta; // holds Im of Dtheta
#include "TList.h"
#include "AliFlowLYZConstants.h"
#include "AliFlowLYZHist1.h"
+#include "TBrowser.h"
class TH1D;
//-----------------------------------------------------------------------
AliFlowLYZHist1::AliFlowLYZHist1():
- TObject(),
+ TNamed(),
fHistGtheta(0),
fHistProReGtheta(0),
fHistProImGtheta(0),
//-----------------------------------------------------------------------
- AliFlowLYZHist1::AliFlowLYZHist1(Int_t theta):
- TObject(),
+ AliFlowLYZHist1::AliFlowLYZHist1(Int_t theta, const char *anInput,const char *aTitle):
+ TNamed(anInput,aTitle),
fHistGtheta(0),
fHistProReGtheta(0),
fHistProImGtheta(0),
}
+//-----------------------------------------------------------------------
+void AliFlowLYZHist1::Print(Option_t *option) const
+{
+ // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
+ // ===============================================
+ // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
+ printf( "Class.Print Name = %s, Histogram list:\n",GetName());
+
+ if (fHistList) {
+ fHistList->Print(option);
+ }
+ else
+ {
+ printf( "Empty histogram list \n");
+ }
+}
+
+//-----------------------------------------------------------------------
+ void AliFlowLYZHist1::Browse(TBrowser *b)
+{
+
+ if (!b) return;
+ if (fHistList) b->Add(fHistList,"AliFlowLYZHist1List");
+}
+
+
+
class TH1D;
class TCollection;
class TList;
-
+class TBrowser;
// Description: Class to organise histograms for Flow
// by the LeeYangZeros method in the first run.
// in the generating function.
-class AliFlowLYZHist1: public TObject {
+class AliFlowLYZHist1: public TNamed {
public:
AliFlowLYZHist1(); //default constructor
- AliFlowLYZHist1(Int_t theta); //constructor
+ AliFlowLYZHist1(Int_t theta, const char *name ,const char *title = "AliFlowLYZHist1"); //constructor
virtual ~AliFlowLYZHist1(); //destructor
+ Bool_t IsFolder() const {return kTRUE;};
+ void Browse(TBrowser *b);
+
void Fill(Double_t f, TComplex c); //fill the histograms
TH1D* FillGtheta(); //fills fHistGtheta
Double_t GetR0(); //get R0
Double_t GetBinCenter(Int_t i); //Get a bincentre of fHistGtheta
Int_t GetNBins(); //Gets Nbins
- TH1D* GetHistGtheta() {return this->fHistGtheta;} ;
+ TH1D* GetHistGtheta() {return this->fHistGtheta;} ;
TProfile* GetHistProReGtheta() {return this->fHistProReGtheta;} ;
TProfile* GetHistProImGtheta() {return this->fHistProImGtheta;} ;
- TList* GetHistList() {return this->fHistList;} ;
+ TList* GetHistList() {return this->fHistList;} ;
virtual Double_t Merge(TCollection *aList); //merge function
-
+ void Print(Option_t* option = "") const; //method to print stats
+
+
private:
AliFlowLYZHist1(const AliFlowLYZHist1& aAnalysis); //copy constructor
TList* fHistList; //list to hold all histograms
- ClassDef(AliFlowLYZHist1,0) // macro for rootcint
+ ClassDef(AliFlowLYZHist1,1) // macro for rootcint
};
#include "TString.h"
#include "TComplex.h"
#include "TList.h"
+#include "TBrowser.h"
class TH1D;
-
// Class to organize the histograms in the second run
// in the Lee Yang Zeros Flow analysis.
// Also contains methods to get values from the histograms
// which are called in AliFlowLeeYandZerosMaker::Finish().
// author: N. van der Kolk (kolk@nikhef.nl)
-
-
ClassImp(AliFlowLYZHist2)
//-----------------------------------------------------------------------
AliFlowLYZHist2::AliFlowLYZHist2():
- TObject(),
+ TNamed(),
fHistProReNumer(0),
fHistProImNumer(0),
fHistProReNumerPt(0),
//-----------------------------------------------------------------------
- AliFlowLYZHist2::AliFlowLYZHist2(Int_t theta):
- TObject(),
+ AliFlowLYZHist2::AliFlowLYZHist2(Int_t theta, const char *anInput,const char *aTitle):
+ TNamed(anInput,aTitle),
fHistProReNumer(0),
fHistProImNumer(0),
fHistProReNumerPt(0),
return (double)iCount;
}
+
+//-----------------------------------------------------------------------
+void AliFlowLYZHist2::Print(Option_t *option) const
+{
+ // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
+ // ===============================================
+ // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
+ printf( "Class.Print Name = %s, Histogram list:\n",GetName());
+
+ if (fHistList) {
+ fHistList->Print(option);
+ }
+ else
+ {
+ printf( "Empty histogram list \n");
+ }
+}
+
+//-----------------------------------------------------------------------
+ void AliFlowLYZHist2::Browse(TBrowser *b)
+{
+
+ if (!b) return;
+ if (fHistList) b->Add(fHistList,"AliFlowLYZHist2List");
+}
class TProfile2D;
class TCollection;
class TList;
+class TBrowser;
// Description: Class to organise histograms for Flow
// by the LeeYangZeros method in the second run.
// which are called in AliFlowLeeYandZerosMaker::Finish().
-class AliFlowLYZHist2: public TObject {
+class AliFlowLYZHist2: public TNamed {
public:
AliFlowLYZHist2(); //default constructor
- AliFlowLYZHist2(Int_t theta); //constructor
+ AliFlowLYZHist2(Int_t theta, const char *name = "AliFlowLYZHist2" ,const char *title = "AliFlowLYZHist2"); //constructor
virtual ~AliFlowLYZHist2(); //destructor
+ Bool_t IsFolder() const {return kTRUE;};
+ void Browse(TBrowser *b);
+
void Fill(Double_t d1,Double_t d2, TComplex c); //fill the histograms
Int_t GetNbinsX()
{Int_t iMaxEtaBins = fHistProReNumer->GetNbinsX(); return iMaxEtaBins;}
TList* GetHistList() {return this->fHistList;}
virtual Double_t Merge(TCollection *aList); //merge function
-
+ void Print(Option_t* option = "") const; //method to print stats
+
private:
AliFlowLYZHist2(const AliFlowLYZHist2& aAnalysis); //copy constructor
TProfile2D* fHistProImNumer2D; //holds Im of Numerator
TList* fHistList; //list to hold all histograms
- ClassDef(AliFlowLYZHist2,0) // macro for rootcint
+ ClassDef(AliFlowLYZHist2,1) // macro for rootcint
};
+++ /dev/null
-/**************************************************************************
- * 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. *
- **************************************************************************/
-
-/*
-$Log$
-*/
-
-#include "Riostream.h"
-#include "AliFlowLeeYangZerosMaker.h"
-#include "AliFlowEvent.h"
-#include "AliFlowSelection.h"
-#include "AliFlowConstants.h" //??
-#include "AliFlowLYZHist1.h"
-#include "AliFlowLYZHist2.h"
-#include "AliFlowLYZConstants.h" //??
-
-#include "TMath.h"
-#include "TComplex.h"
-#include "TProfile.h"
-#include "TObjArray.h"
-#include "TFile.h"
-#include "TVector.h"
-#include "TVector2.h"
-#include "TGraphErrors.h"
-#include "TCanvas.h"
-
-class TTree;
-class TH1F;
-class TH1D;
-class TVector3;
-class TProfile2D;
-class TObject;
-
-//class Riostream; //does not compile
-//class TMath; //does not compile
-//class TVector; //does not compile
-
-//Description: Maker to analyze Flow using the LeeYangZeros method
-// Equation numbers are from Big Paper (BP): Nucl. Phys. A 727, 373 (2003)
-// Practical Guide (PG): J. Phys. G: Nucl. Part. Phys. 30, S1213 (2004)
-// Adapted from StFlowLeeYangZerosMaker.cxx
-// by Markus Oldenberg and Art Poskanzer, LBNL
-// with advice from Jean-Yves Ollitrault and Nicolas Borghini
-//
-//Author: Naomi van der Kolk (kolk@nikhef.nl)
-
-
-
-ClassImp(AliFlowLeeYangZerosMaker)
-
- //-----------------------------------------------------------------------
-
- AliFlowLeeYangZerosMaker::AliFlowLeeYangZerosMaker():
- fFirstRun(kTRUE),
- fUseSum(kTRUE),
- fDebug(kFALSE),
- fHistFileName(0),
- fHistFile(0),
- fSummaryFile(0),
- firstRunFile(0)
-
-{
- //default constructor
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::AliFlowLeeYangZerosMaker default constructor****"<<endl;
-
- fFlowSelect = new AliFlowSelection();
- if (fDebug) { cerr<<"****fFlowSelect in constructor AliFlowLeeYangZerosMaker ("<<fFlowSelect<<")****"<<endl;}
- // output file (histograms)
- TString fHistFileName = "flowLYZAnalysPlot.root" ;
-}
-
-
-AliFlowLeeYangZerosMaker::AliFlowLeeYangZerosMaker(const AliFlowSelection* flowSelect):
- fFirstRun(kTRUE),
- fUseSum(kTRUE),
- fDebug(kFALSE),
- fHistFileName(0),
- fHistFile(0),
- fSummaryFile(0),
- firstRunFile(0)
-{
- //custum constructor
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::AliFlowLeeYangZerosMaker custum constructor****"<<endl;
-
- if(flowSelect) { fFlowSelect = new AliFlowSelection(*flowSelect); }
- else {
- fFlowSelect = new AliFlowSelection() ;
- if (fDebug) cerr<<"****fFlowSelect in constructor AliFlowLeeYangZerosMaker ("<<fFlowSelect<<")****"<<endl;
- }
- // output file (histograms)
- TString fHistFileName = "flowLYZAnalysPlot.root" ;
-}
-
- //-----------------------------------------------------------------------
-
-
- AliFlowLeeYangZerosMaker::~AliFlowLeeYangZerosMaker()
- {
- //default destructor
- if (fDebug) cout<<"****~AliFlowLeeYangZerosMaker****"<<endl;
- delete fHistFile;
- }
-
- //-----------------------------------------------------------------------
-
-Bool_t AliFlowLeeYangZerosMaker::Init()
-{
- //init method
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::Init()****"<<endl;
-
- // Open output files (->plots)
- fHistFile = new TFile(fHistFileName.Data(), "RECREATE");
- //fHistFile->cd() ; //all histograms will be saved in this file
-
- //for each harmonic ???
- fQsum.Set(0.,0.);
- fQ2sum = 0.;
-
- // Book histograms
- Int_t fNtheta = AliFlowLYZConstants::kTheta;
-
-
- //for control histograms
- fHistOrigMult = new TH1F("Control_FlowLYZ_OrigMult", "Control_FlowLYZ_OrigMult",1000, 0., 10000.);
- fHistOrigMult->SetXTitle("Original Multiplicity");
- fHistOrigMult->SetYTitle("Counts");
-
- fHistMult = new TH1F("Control_FlowLYZ_Mult", "Control_FlowLYZ_Mult",1000, 0., 10000.);
- fHistMult->SetXTitle("Multiplicity from selection");
- fHistMult->SetYTitle("Counts");
-
- fHistQ = new TH1F("Control_FlowLYZ_Q","Control_FlowLYZ_Q",500, 0., 10.);
- fHistQ->SetXTitle("Qvector");
- fHistQ->SetYTitle("Counts");
-
- fHistPt = new TH1F("Control_FlowLYZ_Pt","Control_FlowLYZ_Pt",200, 0., 10.);
- fHistPt->SetXTitle("Pt (GeV/c)");
- fHistPt->SetYTitle("Counts");
-
- fHistPhi = new TH1F("Control_FlowLYZ_Phi","Control_FlowLYZ_Phi",70, 0., 7.);
- fHistPhi->SetXTitle("Phi");
- fHistPhi->SetYTitle("Counts");
-
- fHistEta = new TH1F("Control_FlowLYZ_Eta","Control_FlowLYZ_Eta",40, 0., 2.);
- fHistEta->SetXTitle("Eta");
- fHistEta->SetYTitle("Counts");
-
- fHistQtheta = new TH1F("Control_FlowLYZ_Qtheta","Control_FlowLYZ_Qtheta",50,-1000.,1000.);
- fHistQtheta->SetXTitle("Qtheta");
- fHistQtheta->SetYTitle("Counts");
-
- if (fFirstRun){
- //for first loop over events
- fHistProR0thetaHar1 = new TProfile("First_FlowProLYZ_r0theta_Har1","First_FlowProLYZ_r0theta_Har1",fNtheta,-0.5,fNtheta-0.5);
- fHistProR0thetaHar1->SetXTitle("#theta");
- fHistProR0thetaHar1->SetYTitle("r_{0}^{#theta}");
-
- fHistProR0thetaHar2 = new TProfile("First_FlowProLYZ_r0theta_Har2","First_FlowProLYZ_r0theta_Har2",fNtheta,-0.5,fNtheta-0.5);
- fHistProR0thetaHar2->SetXTitle("#theta");
- fHistProR0thetaHar2->SetYTitle("r_{0}^{#theta}");
-
- fHistProVthetaHar1 = new TProfile("First_FlowProLYZ_Vtheta_Har1","First_FlowProLYZ_Vtheta_Har1",fNtheta,-0.5,fNtheta-0.5);
- fHistProVthetaHar1->SetXTitle("#theta");
- fHistProVthetaHar1->SetYTitle("V_{n}^{#theta}");
-
- fHistProVthetaHar2 = new TProfile("First_FlowProLYZ_Vtheta_Har2","First_FlowProLYZ_Vtheta_Har2",fNtheta,-0.5,fNtheta-0.5);
- fHistProVthetaHar2->SetXTitle("#theta");
- fHistProVthetaHar2->SetYTitle("V_{n}^{#theta}");
-
- fHistProVR0 = new TProfile("First_FlowProLYZ_vR0","First_FlowProLYZ_vR0",2,0.5,2.5,-100.,100.);
- fHistProVR0->SetXTitle("Harmonic");
- fHistProVR0->SetYTitle("v integrated from r0 (%)");
-
- fHistVR0 = new TH1D("First_FlowLYZ_vR0","First_FlowLYZ_vR0",2,0.5,2.5);
- fHistVR0->SetXTitle("Harmonic");
- fHistVR0->SetYTitle("v integrated from r0 (%)");
-
- fHistProV = new TProfile("First_FlowProLYZ_V","First_FlowProLYZ_V",2,0.5,2.5,-1000.,1000.);
- fHistProV->SetXTitle("Harmonic");
- fHistProV->SetYTitle("v integrated");
-
- //class AliFlowLYZHist1 defines the histograms: fHistProGtheta, fHistProReGtheta, fHistProImGtheta, fHistProR0theta
- for (Int_t j=0;j<AliFlowConstants::kHars;j++)
- {
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- fHist1[j][theta]=new AliFlowLYZHist1(theta,j+1);
- }
- }
- }
- else {
- //for second loop over events
- fHistProReDenomHar1 = new TProfile("Second_FlowProLYZ_ReDenom_Har1","Second_FlowProLYZ_ReDenom_Har1" , fNtheta, -0.5, fNtheta-0.5);
- fHistProReDenomHar1->SetXTitle("#theta");
- fHistProReDenomHar1->SetYTitle("Re(Q^{#theta}e^{ir_{0}^{#theta}Q^{#theta}})");
-
- fHistProReDenomHar2 = new TProfile("Second_FlowProLYZ_ReDenom_Har2","Second_FlowProLYZ_ReDenom_Har2" , fNtheta, -0.5, fNtheta-0.5);
- fHistProReDenomHar2->SetXTitle("#theta");
- fHistProReDenomHar2->SetYTitle("Re(Q^{#theta}e^{ir_{0}^{#theta}Q^{#theta}})");
-
- fHistProImDenomHar1 = new TProfile("Second_FlowProLYZ_ImDenom_Har1","Second_FlowProLYZ_ImDenom_Har1" , fNtheta, -0.5, fNtheta-0.5);
- fHistProImDenomHar1->SetXTitle("#theta");
- fHistProImDenomHar1->SetYTitle("Im(Q^{#theta}e^{ir_{0}^{#theta}Q^{#theta}})");
-
- fHistProImDenomHar2 = new TProfile("Second_FlowProLYZ_ImDenom_Har2","Second_FlowProLYZ_ImDenom_Har2" , fNtheta, -0.5, fNtheta-0.5);
- fHistProImDenomHar2->SetXTitle("#theta");
- fHistProImDenomHar2->SetYTitle("Im(Q^{#theta}e^{ir_{0}^{#theta}Q^{#theta}})");
-
- fHistProVetaHar1 = new TProfile("Second_FlowProLYZ_Veta_Har1","Second_FlowProLYZ_Veta_Har1",40,-2.,2.);
- fHistProVetaHar1->SetXTitle("rapidity");
- fHistProVetaHar1->SetYTitle("v (%)");
-
- fHistProVetaHar2 = new TProfile("Second_FlowProLYZ_Veta_Har2","Second_FlowProLYZ_Veta_Har2",40,-2.,2.);
- fHistProVetaHar2->SetXTitle("rapidity");
- fHistProVetaHar2->SetYTitle("v (%)");
-
- fHistProVPtHar1 = new TProfile("Second_FlowProLYZ_VPt_Har1","Second_FlowProLYZ_VPt_Har1",100,0.,10.);
- fHistProVPtHar1->SetXTitle("Pt");
- fHistProVPtHar1->SetYTitle("v (%)");
-
- fHistProVPtHar2 = new TProfile("Second_FlowProLYZ_VPt_Har2","Second_FlowProLYZ_VPt_Har2",100,0.,10.);
- fHistProVPtHar2->SetXTitle("Pt");
- fHistProVPtHar2->SetYTitle("v (%)");
-
- fHistVPtHar2 = new TH1D("Second_FlowLYZ_VPt_Har2","Second_FlowLYZ_VPt_Har2",100,0.,10.);
- fHistVPtHar2->SetXTitle("Pt");
- fHistVPtHar2->SetYTitle("v (%)");
-
- fHistProReDtheta = new TProfile("Second_FlowProLYZ_ReDtheta_Har2","Second_FlowProLYZ_ReDtheta_Har2",fNtheta, -0.5, fNtheta-0.5);
- fHistProReDtheta->SetXTitle("#theta");
- fHistProReDtheta->SetYTitle("Re(D^{#theta})");
-
- fHistProImDtheta = new TProfile("Second_FlowProLYZ_ImDtheta_Har2","Second_FlowProLYZ_ImDtheta_Har2",fNtheta, -0.5, fNtheta-0.5);
- fHistProImDtheta->SetXTitle("#theta");
- fHistProImDtheta->SetYTitle("Im(D^{#theta})");
-
-
- //class AliFlowLYZHist2 defines the histograms:
- for (Int_t j=0;j<AliFlowConstants::kHars;j++)
- {
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- fHist2[j][theta]=new AliFlowLYZHist2(theta,j+1);
- }
- }
-
- //read hists from first run file
- //firstRunFile = new TFile("fof_flowLYZAnal_firstrun.root","READ"); //default is read
- if (firstRunFile->IsZombie()){ //check if file exists
- cout << "Error opening file, run first with fFirstrun = kTRUE" << endl;
- exit(-1);
- } else if (firstRunFile->IsOpen()){
- cout<<"----firstRunFile is open----"<<endl<<endl;
- fHistProVthetaHar1 = (TProfile*)firstRunFile->Get("First_FlowProLYZ_Vtheta_Har1");
- fHistProVthetaHar2 = (TProfile*)firstRunFile->Get("First_FlowProLYZ_Vtheta_Har2");
- fHistProR0thetaHar2 = (TProfile*)firstRunFile->Get("First_FlowProLYZ_r0theta_Har2");
- fHistProV = (TProfile*)firstRunFile->Get("First_FlowProLYZ_V");
- }
- }
-
-
- if (fDebug) cout<<"****Histograms initialised****"<<endl;
- if (fDebug) cout<<"****fFlowSelect in Init() "<<fFlowSelect<<"****"<<endl;
-
- fEventNumber = 0; //set event counter to zero
- /*
- if (fUseSum)
- {
- //initialize LYZ summary class
- fLYZSummary = new AliFlowLYZSummary();
- fSummaryFile = new TFile("fFlowSummary.root","RECREATE","Flow LYZ summary file");
- fSummaryFile->SetCompressionLevel(1);
- fFlowTree = new TTree("FlowTree", "Flow Summary Tree");
- fFlowTree->SetAutoSave(1000000); // autosave when 1 Mbyte written
- fFlowTree->Branch("fLYZSummary","AliFlowLYZSummary",&fLYZSummary,25000,99);
- }
- */
- return kTRUE;
-}
-
- //-----------------------------------------------------------------------
-
-Bool_t AliFlowLeeYangZerosMaker::Make(AliFlowEvent* fFlowEvent)
-{
- //make method
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::Make()****"<<endl;
-
- //get tracks from event
- if (fFlowEvent) {
- fFlowTracks = fFlowEvent->TrackCollection();
- if (fDebug) cout<<"****fFlowSelect in Make() "<<fFlowSelect<<"****"<<endl;
- if (fDebug) fFlowSelect->PrintSelectionList() ;
- if (fDebug) fFlowSelect->PrintList() ;
-
- if (fFlowSelect && fFlowSelect->Select(fFlowEvent)) // check if event is selected
- {
- fFlowTracks = fFlowEvent->TrackCollection(); //get tracks from event
- fFlowEvent->SetSelections(fFlowSelect) ; // does the selection of tracks for r.p. calculation (sets flags in AliFlowTrack)
- if (fFirstRun){
- MakeControlHistograms(fFlowEvent);
- FillFromFlowEvent(fFlowEvent);
- }
- else {
- MakeControlHistograms(fFlowEvent);
- SecondFillFromFlowEvent(fFlowEvent);
- }
- }
- }
- else {
- cout<<"##### FlowLeeYangZero: FlowEvent pointer null"<<endl;
- return kFALSE;
- }
-
- fEventNumber++;
-
- return kTRUE;
-}
-
-
- //-----------------------------------------------------------------------
- Bool_t AliFlowLeeYangZerosMaker::Finish()
-{
- //finish method
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::Finish()****"<<endl;
-
- //define variables for both runs
- Float_t fR0 = 0;
- Float_t fv = 0;
- Float_t fVtheta = 0;
- Float_t fSigma2 = 0;
- Float_t fChi= 0;
- Float_t fJ01 = 2.405;
- Int_t fNtheta = AliFlowLYZConstants::kTheta;
- Float_t fV = 0;
-
- if (fFirstRun){
- for (Int_t j=0;j<AliFlowConstants::kHars;j++)
- //loop over harmonics j=0 ->first harmonic, j=1 ->second harmonic
- {
- Float_t fMeanMult = fHistMult->GetMean();
- //cerr<<"fMeanMult = "<<fMeanMult<<endl;
-
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- //get the first minimum r0
- fHist1[j][theta]->FillGtheta();
- fR0 = fHist1[j][theta]->GetR0();
- //cerr<<"fR0 = "<<fR0<<endl;
-
- //calculate integrated flow
- if (fR0!=0) fVtheta = fJ01/fR0;
- if (fMeanMult!=0.)
- {
- fv = fVtheta/fMeanMult;
- cerr<<"fv = "<<fv<<endl;
- }
-
-
- //fill the histograms
- if (j==0)
- {
- fHistProR0thetaHar1->Fill(theta,fR0);
- fHistProVthetaHar1->Fill(theta,fVtheta);
- fHistProV->Fill(j+1,fVtheta); //profile takes care of the averaging over theta.
- }
- if (j==1)
- {
- fHistProR0thetaHar2->Fill(theta,fR0);
- fHistProVthetaHar2->Fill(theta,fVtheta);
- fHistProV->Fill(j+1,fVtheta); //profile takes care of the averaging over theta.
- }
-
- fHistProVR0->Fill(j+1,fv*100); //*100 to get % //profile takes care of the averaging over theta.
-
- } //end of loop over theta
-
- //sigma2 and chi
- if (j==1) //second harmonic only temporarily
- {
- if (fEventNumber!=0) {
- fQsum /= fEventNumber;
- //cerr<<"fQsum.X() = "<<fQsum.X()<<endl;
- //cerr<<"fQsum.Y() = "<<fQsum.Y()<<endl;
- fQ2sum /= fEventNumber;
- //cerr<<"fQ2sum = "<<fQ2sum<<endl;
- fV = fHistProV->GetBinContent(j+1);
- fSigma2 = fQ2sum - TMath::Power(fQsum.X(),2.) - TMath::Power(fQsum.Y(),2.) - TMath::Power(fV,2.); //BP eq. 62
- if (fSigma2>0) fChi = fV/TMath::Sqrt(fSigma2);
- else fChi = -1.;
- cerr<<"fV = "<<fV<<" and chi = "<<fChi<<endl;
- }
-
- } //j==1
-
- } //end of loop over harmonics
-
- // recalculate statistical errors on integrated flow
- //combining 5 theta angles to 1 relative error BP eq. 89
- Double_t fRelErr2comb = 0.;
- Int_t nEvts = fEventNumber;
- for (Int_t theta=0;theta<fNtheta;theta++){
- Double_t fTheta = ((double)theta/fNtheta)*TMath::Pi();
- Double_t fApluscomb = TMath::Exp((fJ01*fJ01)/(2*fChi*fChi)*
- TMath::Cos(fTheta));
- Double_t fAmincomb = TMath::Exp(-(fJ01*fJ01)/(2*fChi*fChi)*
- TMath::Cos(fTheta));
- fRelErr2comb += (1/(2*nEvts*(fJ01*fJ01)*TMath::BesselJ1(fJ01)*
- TMath::BesselJ1(fJ01)))*
- (fApluscomb*TMath::BesselJ0(2*fJ01*TMath::Sin(fTheta/2)) +
- fAmincomb*TMath::BesselJ0(2*fJ01)*TMath::Cos(fTheta/2));
- }
- fRelErr2comb /= fNtheta;
- Double_t fRelErrcomb = TMath::Sqrt(fRelErr2comb);
- cerr<<"fRelErrcomb = "<<fRelErrcomb<<endl;
-
- //copy content of profile into TH1D and add error
- for(Int_t b=0;b<2;b++){
- Double_t fv2pro = fHistProVR0->GetBinContent(b+1);
- fHistVR0->SetBinContent(b+1,fv2pro);
- Double_t fv2Err = fv2pro*fRelErrcomb ;
- cerr<<"fv2pro +- fv2Err = "<<fv2pro<<" +- "<<fv2Err<<" for bin "<<b+1<<endl;
- fHistVR0->SetBinError(b+1,fv2Err);
- }
-
-
- if (fDebug) cout<<"****histograms filled****"<<endl;
- /*
- if (fUseSum)
- {
- if (fSummaryFile->IsOpen())
- {
- fSummaryFile->Write(0,TObject::kOverwrite);
- fSummaryFile->Close();
- }
- }
- */
-
- //save histograms in file //temp for testing selector
- fHistFile->cd();
- fHistFile->Write();
-
- return kTRUE;
- fEventNumber =0; //set to zero for second round over events
- } //firstrun
-
-
- else { //second run
-
- //calculate differential flow
- //declare variables
- Float_t fEta, fPt, fReRatio, fVeta, fVPt;
- Float_t fReDenom = 0;
- Float_t fImDenom = 0;
- Double_t fR0 = 0;
- TComplex fDenom, fNumer, fDtheta;
- Int_t m = 1;
- TComplex i = TComplex::I();
- Float_t fBesselRatio[3] = {1., 1.202, 2.69};
- Double_t fErrdifcomb = 0.; //set error to zero
- Double_t fErr2difcomb = 0.; //set error to zero
-
- /*
- //v1 integrated
- for (Int_t j=0;j<AliFlowConstants::kHars;j++) //loop over harmonics j=0 ->first harmonic, j=1 ->second harmonic
- {
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- //get the first minimum r0
- //fR0 = GetR0(fHist1[j][theta]->FillGtheta());
-
- fReDenom = fHistProReDenomHar2->GetBinContent(theta+1);
- fImDenom = fHistProImDenomHar2->GetBinContent(theta+1);
- TComplex fDenom(fReDenom,fImDenom);
-
- //complete!!
-
- }//end of loop over theta
- }//end of loop over harmonics
- */
-
- //differential flow
- //for (Int_t j=0;j<AliFlowConstants::kHars;j++) //loop over harmonics j=0 ->first harmonic, j=1 ->second harmonic
- //{
- Int_t j=1; //temp only harm 2
- //fFlowSelect->SetHarmonic(j); //not needed here?
-
- for (Int_t theta=0;theta<fNtheta;theta++) { //loop over theta
- if (j==0) {
- fR0 = fHistProR0thetaHar1->GetBinContent(theta+1);
- fVtheta = 2.405/fR0;
- //fVtheta = fHistProVthetaHar1->GetBinContent(theta+1); // BP Eq. 9 -> Vn^theta = j01/r0^theta
- fReDenom = fHistProReDenomHar1->GetBinContent(theta+1);
- fImDenom = fHistProImDenomHar1->GetBinContent(theta+1);
- }
- if (j==1) {
- fR0 = fHistProR0thetaHar2->GetBinContent(theta+1);
- if (fDebug) cerr<<"fR0 = "<<fR0<<endl;
- fVtheta = 2.405/fR0; // BP Eq. 9 -> Vn^theta = j01/r0^theta
- fReDenom = fHistProReDenomHar2->GetBinContent(theta+1);
- fImDenom = fHistProImDenomHar2->GetBinContent(theta+1);
- }
-
- fDenom(fReDenom,fImDenom);
-
-
- //for new method and use by others (only with the sum generating function):
- if (fUseSum) {
- fR0 = fHistProR0thetaHar2->GetBinContent(theta+1);
- fDtheta = fR0*fDenom;
- fHistProReDtheta->Fill(theta,fDtheta.Re());
- fHistProImDtheta->Fill(theta,fDtheta.Im());
- }
-
- fDenom *= TComplex::Power(i, m-1);
- //cerr<<"TComplex::Power(i, m-1) = "<<TComplex::Power(i, m-1).Rho()<<endl; //checked ok
-
- //v as a function of eta
- Int_t fEtaBins = AliFlowLYZConstants::kEtaBins;
- for (Int_t be=1;be<=fEtaBins;be++) {
- fEta = fHist2[j][theta]->GetBinCenter(be);
- fNumer = fHist2[j][theta]->GetfNumer(be);
- if (fNumer.Rho()==0) {
- if (fDebug) cerr<<"WARNING: modulus of fNumer is zero in Finish()"<<endl;
- fReRatio = 0;
- }
- else if (fDenom.Rho()==0) {
- if (fDebug) cerr<<"WARNING: modulus of fDenom is zero"<<endl;
- fReRatio = 0;
- }
- else {
- //if ( j==1 && theta==0) cerr<<"modulus of fNumer = "<<fNumer.Rho() <<endl; //always a number smaller than 1, or 0.
- fReRatio = (fNumer/fDenom).Re();
- }
-
- fVeta = fBesselRatio[m-1]*fReRatio*fVtheta*100.; //BP eq. 12
- //if ( j==1 && theta==0) cerr<<"eta = "<<fEta<<" cerr::fReRatio for eta = "<<fReRatio<<" cerr::fVeta for eta = "<<fVeta<<endl;
-
- if (j==0) fHistProVetaHar1->Fill(fEta,fVeta);
- if (j==1) fHistProVetaHar2->Fill(fEta,fVeta);
- } //loop over bins be
-
- //v as a function of Pt
- Int_t fPtBins = AliFlowLYZConstants::kPtBins;
- for (Int_t bp=1;bp<=fPtBins;bp++) {
- fPt = fHist2[j][theta]->GetBinCenterPt(bp);
- fNumer = fHist2[j][theta]->GetfNumerPt(bp);
- if (fNumer.Rho()==0) {
- if (fDebug) cerr<<"modulus of fNumer is zero"<<endl;
- fReRatio = 0;
- }
- else if (fDenom.Rho()==0) {
- if (fDebug) cerr<<"modulus of fDenom is zero"<<endl;
- fReRatio = 0;
- }
- else {
- //if ( j==1 && theta==0) cerr<<"modulus of fNumer = "<<fNumer.Rho() <<endl; //always a number smaller than 1, or 0.
- fReRatio = (fNumer/fDenom).Re();
- }
-
- fVPt = fBesselRatio[m-1]*fReRatio*fVtheta*100.; //BP eq. 12
- //cerr<<"fBesselRatio[m-1] = "<<fBesselRatio[m-1]<<endl; //checked ok
- //if ( j==1 && theta==0) cerr<<"pt = "<<fPt<<" cerr::fReRatio for pt = "<<fReRatio<<" cerr::fVPt for pt = "<<fVeta<<endl;
-
- if (j==0) fHistProVPtHar1->Fill(fPt,fVPt);
- if (j==1) fHistProVPtHar2->Fill(fPt,fVPt);
- } //loop over bins bp
-
- }//end of loop over theta
-
-
- //sigma2 and chi (for statistical error calculations)
- if (j==1) { //second harmonic only temporarily
-
- if (fEventNumber!=0) {
- fQsum /= fEventNumber;
- //cerr<<"fQsum.X() = "<<fQsum.X()<<endl;
- //cerr<<"fQsum.Y() = "<<fQsum.Y()<<endl;
- fQ2sum /= fEventNumber;
- cerr<<"fQ2sum = "<<fQ2sum<<endl;
- fV = fHistProV->GetBinContent(j+1);
- fSigma2 = fQ2sum - TMath::Power(fQsum.X(),2.) - TMath::Power(fQsum.Y(),2.) - TMath::Power(fV,2.); //BP eq. 62
- if (fSigma2>0) fChi = fV/TMath::Sqrt(fSigma2);
- else fChi = -1.;
- cerr<<"fV = "<<fV<<" and chi = "<<fChi<<endl;
- }
-
- } //j==1
-
-
- //copy content of profile into TH1D and add error
- for(Int_t b=0;b<100;b++){
- Double_t fv2pro = fHistProVPtHar2->GetBinContent(b);
- //calculate error
- for (Int_t theta=0;theta<fNtheta;theta++) {
- Double_t fTheta = ((double)theta/fNtheta)*TMath::Pi();
- Int_t Nprime = fHist2[j][theta]->GetNprimePt(b);
- //cerr<<"Nprime = "<<Nprime<<endl;
- if (Nprime!=0.) {
- Double_t fApluscomb = TMath::Exp((fJ01*fJ01)/(2*fChi*fChi)*
- TMath::Cos(fTheta));
- Double_t fAmincomb = TMath::Exp(-(fJ01*fJ01)/(2*fChi*fChi)*
- TMath::Cos(fTheta));
- fErr2difcomb += (TMath::Cos(fTheta)/(4*Nprime*TMath::BesselJ1(fJ01)*
- TMath::BesselJ1(fJ01)))*
- ((fApluscomb*TMath::BesselJ0(2*fJ01*TMath::Sin(fTheta/2))) -
- (fAmincomb*TMath::BesselJ0(2*fJ01*TMath::Cos(fTheta/2))));
- }
- } //loop over theta
-
- if (fErr2difcomb!=0.) {
- fErr2difcomb /= fNtheta;
- fErrdifcomb = TMath::Sqrt(fErr2difcomb)*100;
- //cerr<<"fErrdifcomb = "<<fErrdifcomb<<endl;
- }
- else {fErrdifcomb = 0.;}
-
- //fill TH1D
- if (j==1) {
- fHistVPtHar2->SetBinContent(b,fv2pro);
- fHistVPtHar2->SetBinError(b,fErrdifcomb);
- }
- //check that error is set
- //if (j==1) { cout<<"difference between calculated error and error in hostogram: "<< fErrdifcomb - fHistVPtHar2->GetBinError(b)<<endl; }
-
- } //loop over bins b
-
-
- //} //end of loop over harmonics //temporarily out
-
- //save histograms in file
- fHistFile->cd();
- fHistFile->Write();
- fHistFile->Close();
- //Note that when the file is closed, all histograms and Trees in memory associated with this file are deleted
- if (fDebug) cout<<"****Histograms saved and fHistFile closed, all histograms deleted****"<<endl;
-
- //close the first run file
- firstRunFile->Close();
-
-
- } //secondrun
-
- delete fFlowSelect;
- cout<<"----LYZ analysis finished....----"<<endl<<endl;
-
- return kTRUE;
-}
-
-
-//-----------------------------------------------------------------------
- Bool_t AliFlowLeeYangZerosMaker::MakeControlHistograms(AliFlowEvent* fFlowEvent)
-{
- //contol histograms of pt, eta, phi, Q, mult
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::MakeControlHistograms()****"<<endl;
-
- //define variables
- TVector2 fQ;
- Float_t fPt, fPhi, fEta, fQmult;
-
- if (!fFlowEvent){
- cout<<"##### FlowLeeYangZero: FlowEvent pointer null"<<endl;
- return kFALSE;
- }
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
- //set selection and harmonic
- fFlowSelect->SetSelection(1);
- fFlowSelect->SetHarmonic(1); //second harmonic
-
- //cerr<<"selection in MakeControlHistograms()"<<endl;
- //fFlowSelect->PrintSelectionList() ;
- //fFlowSelect->PrintList() ;
-
- fFlowTracks = fFlowEvent->TrackCollection();
- Int_t fNumberOfTracks = fFlowTracks->GetEntries();
- fHistOrigMult->Fill(fNumberOfTracks);
- //cerr<<"fNumberOfTracks = "<<fNumberOfTracks<<endl;
- Int_t fMult = fFlowEvent->Mult(fFlowSelect); // Multiplicity of tracks in the specified Selection
- fHistMult->Fill(fMult);
- //cerr<<"Mult = "<<fMult<<endl;
-
- fQ = fFlowEvent ->Q(fFlowSelect);
- fQmult = fQ.Mod()/TMath::Sqrt(fNumberOfTracks);
- fHistQ->Fill(fQmult);
-
- Int_t tempmult = 0; //for testing
- for (Int_t i=0;i<fNumberOfTracks;i++)
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- //if (fFlowSelect->SelectPart(fFlowTrack)) //if track is selected
- if (fFlowSelect->Select(fFlowTrack))
- {
- fPt = fFlowTrack->Pt();
- fHistPt->Fill(fPt);
- tempmult++;
- fPhi = fFlowTrack->Phi();
- if (fPhi<0.) fPhi+=2*TMath::Pi();
- fHistPhi->Fill(fPhi);
- fEta = fFlowTrack->Eta();
- fHistEta->Fill(fEta);
- }
- }
- if (fMult!=tempmult){cerr<<"ERROR: Mult() is not tempmult! "<<fMult<<" :: "<<tempmult<<endl<<endl;}
- //else {cerr<<"Mult()= tempmult "<<fMult<<" :: "<<tempmult<<endl<<endl;}
-
- return kTRUE;
-
-
-}
-
-
-
-//-----------------------------------------------------------------------
-
- Bool_t AliFlowLeeYangZerosMaker::FillFromFlowEvent(AliFlowEvent* fFlowEvent)
-{
- // Get event quantities from AliFlowEvent for all particles
-
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::FillFromFlowEvent()****"<<endl;
-
- if (!fFlowEvent){
- cout<<"##### FlowLeeYangZero: FlowEvent pointer null"<<endl;
- return kFALSE;
- }
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
-
- //define variables
- TComplex fExpo, fGtheta, fGthetaNew, fZ;
- //Int_t m;
- Int_t fNtheta = AliFlowLYZConstants::kTheta;
- Int_t fNbins = AliFlowLYZConstants::kNbins;
-
-
- //calculate flow
- fFlowSelect->SetSelection(1);
-
- for (Int_t j=0;j<AliFlowConstants::kHars;j++) //loop over harmonics j=0 ->first harmonic, j=1 ->second harmonic
- {
- Int_t m=1;
- fFlowSelect->SetHarmonic(j);
- Float_t fOrder = (double)((j+1)/m);
- //cerr<<"fOrder = "<<fOrder<<endl;
-
- //get the Q vector from the FlowEvent
- TVector2 fQ = fFlowEvent->Q(fFlowSelect);
- //for chi calculation:
- if (j==1) //second harmonic only temporarily
- {
- fQsum += fQ;
- fQ2sum += fQ.Mod2();
- //cerr<<"fQ2sum = "<<fQ2sum<<endl;
- }
-
- fFlowTracks = fFlowEvent->TrackCollection();
-
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- fTheta = ((float)theta/fNtheta)*TMath::Pi()/fOrder;
- //cerr<<"fTheta = "<<fTheta<<endl;
-
- //calculate fQtheta = cos(fOrder*(fPhi-fTheta);the projection of the Q vector on the reference direction fTheta
- fQtheta = GetQtheta(fFlowSelect,fFlowTracks,fTheta);
- //save fQtheta in AliFlowLYZSummary class
-
- //something
- //AliFlowLYZSummary::SetQtheta(theta,fQtheta);
-
- if (j==1 && theta==0) fHistQtheta->Fill(fQtheta);
- //cerr<<"fQtheta = "<<fQtheta<<endl;
-
- for (Int_t bin=1;bin<=fNbins;bin++)
- {
- Float_t fR = fHist1[j][theta]->GetBinCenter(bin); //bincentre of bins in histogram
- //if (theta == 0) cerr<<"cerr::fR = "<<fR<<endl;
- if (fUseSum)
- {
- //calculate the sum generating function
- fExpo(0.,fR*fQtheta); //Re=0 ; Im=fR*fQtheta
- fGtheta = TComplex::Exp(fExpo);
- }
- else
- {
- //calculate the product generating function
- fGtheta = GetGrtheta(fFlowSelect,fR,fTheta); //make this function
- if (fGtheta.Rho2() > 100.) break;
- }
-
- fHist1[j][theta]->Fill(fR,fGtheta); //fill real and imaginary part of fGtheta
-
- } //loop over bins
- } //loop over theta
- } //loop over harmonics
-
- return kTRUE;
-
-
-}
-
- //-----------------------------------------------------------------------
- Bool_t AliFlowLeeYangZerosMaker::SecondFillFromFlowEvent(AliFlowEvent* fFlowEvent)
-{
- //for differential flow
-
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::SecondFillFromFlowEvent()****"<<endl;
-
- if (!fFlowEvent){
- cout<<"##### FlowLeeYangZero: FlowEvent pointer null"<<endl;
- return kFALSE;
- }
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
- //define variables
- //TVector2 fQ;
- TComplex fExpo, fDenom, fNumer,fCosTermComplex;
- Float_t fOrder, fR0, fPhi, fEta, fPt;
- Double_t fCosTerm;
- Double_t m = 1.;
- Int_t fNtheta = AliFlowLYZConstants::kTheta;
-
- //calculate flow
- fFlowSelect->SetSelection(1);
-
- //cerr<<"selection in SecondFillFromFlowEvent()"<<endl;
- //fFlowSelect->PrintSelectionList() ;
- //fFlowSelect->PrintList() ;
-
-
- for (Int_t j=0;j<AliFlowConstants::kHars;j++) //loop over harmonics j=0 ->first harmonic, j=1 ->second harmonic
- {
- m=1;
- fFlowSelect->SetHarmonic(j);
- fOrder = (double)((j+1)/m);
-
- //get the Q vector from the FlowEvent
- TVector2 fQ = fFlowEvent->Q(fFlowSelect);
- //for chi calculation:
- if (j==1) //second harmonic only temporarily
- {
- fQsum += fQ;
- fQ2sum += fQ.Mod2();
- //cerr<<"fQ2sum = "<<fQ2sum<<endl;
- }
-
- fFlowTracks = fFlowEvent->TrackCollection();
-
- for (Int_t theta=0;theta<fNtheta;theta++)
- {
- fTheta = ((float)theta/fNtheta)*TMath::Pi()/fOrder;
-
- //calculate fQtheta = cos(fOrder*(fPhi-fTheta);the projection of the Q vector on the reference direction fTheta
- fQtheta = GetQtheta(fFlowSelect,fFlowTracks,fTheta);
- //cerr<<"fQtheta for fdenom = "<<fQtheta<<endl;
-
- /*
- if (j==0) //first harmonic
- {
- //denominator for differential v
- fR0 = fHistProR0thetaHar1->GetBinContent(theta+1);
- fExpo(0.,fR0*fQtheta);
- fDenom = fQtheta*(TComplex::Exp(fExpo)); //BP eq 12
- //cerr<<"fDenom.Re() = "<<fDenom.Re()<<endl;
- //cerr<<"fDenom.Im() = "<<fDenom.Im()<<endl;
-
- //denominator for differential v
- // ****put in product generating function!!
-
-
- fHistProReDenomHar1->Fill(theta,fDenom.Re()); //fill the real part of fDenom
- fHistProImDenomHar1->Fill(theta,fDenom.Im()); //fill the imaginary part of fDenom
- }
- */
-
- if (j==1) //second harmonic
- {
- //denominator for differential v
- fR0 = fHistProR0thetaHar2->GetBinContent(theta+1);
- //cerr<<"fR0 = "<<fR0 <<endl;
-
- if (fUseSum) //sum generating function
- {
- fExpo(0.,fR0*fQtheta);
- fDenom = fQtheta*(TComplex::Exp(fExpo)); //BP eq 12
- //loop over tracks in event
- fFlowTracks = fFlowEvent->TrackCollection();
- Int_t fNumberOfTracks = fFlowTracks->GetEntries();
- for (Int_t i=0;i<fNumberOfTracks;i++)
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- if (fFlowSelect->SelectPart(fFlowTrack)) //if track is selected
- {
- fPhi = fFlowTrack->Phi();
- fCosTerm = cos(m*fOrder*(fPhi-fTheta));
- //cerr<<"fCosTerm = "<<fCosTerm <<endl;
- fNumer = fCosTerm*(TComplex::Exp(fExpo));
- if (fNumer.Rho()==0) {cerr<<"WARNING: modulus of fNumer is zero in SecondFillFromFlowEvent"<<endl;}
- if (fDebug) cerr<<"modulus of fNumer is "<<fNumer.Rho()<<endl;
- fEta = fFlowTrack->Eta(); //rapidity
- fPt = fFlowTrack->Pt();
- fHist2[j][theta]->Fill(fEta,fPt,fNumer);
- } //if
- } //loop over tracks
-
- } //sum
- else //product generating function
- {
- fDenom = GetDiffFlow(fFlowSelect,fR0,theta);
-
- }//product
-
- fHistProReDenomHar2->Fill(theta,fDenom.Re()); //fill the real part of fDenom
- fHistProImDenomHar2->Fill(theta,fDenom.Im()); //fill the imaginary part of fDenom
- }//j==1
-
-
- }//end of loop over theta
-
- }//loop over harmonics
-
-
- return kTRUE;
-
-
-
-}
- //-----------------------------------------------------------------------
- Double_t AliFlowLeeYangZerosMaker::GetQtheta(AliFlowSelection* fFlowSelect, TObjArray* fFlowTracks, Float_t fTheta)
-{
- //calculate Qtheta. Qtheta is the sum over all particles of cos(fOrder*(fPhi-fTheta)) BP eq. 3
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::GetQtheta()****"<<endl;
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
-
- Double_t fQtheta = 0.;
- Int_t fHarmonic = fFlowSelect->Har();
- Double_t fOrder = (double)(fHarmonic+1);
-
- Int_t fNumberOfTracks = fFlowTracks->GetEntries();
- //cerr<<"GetQtheta::fNumberOfTracks = "<<fNumberOfTracks<<endl;
-
- for (Int_t i=0;i<fNumberOfTracks;i++) //loop over tracks in event
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- if(fFlowSelect->Select(fFlowTrack)) //if track is selected //gives the same number of particles as Mult(fFlowSelect) method
- {
- Float_t fPhi = fFlowTrack->Phi();
- fQtheta += cos(fOrder*(fPhi-fTheta));
- }
-
- }//loop over tracks
-
- return fQtheta;
-
-}
-
-//-----------------------------------------------------------------------
-TComplex AliFlowLeeYangZerosMaker::GetGrtheta(AliFlowSelection* fFlowSelect, Float_t fR, Float_t fTheta)
-{
- // Product Generating Function for LeeYangZeros method
- // PG Eq. 3 (J. Phys. G Nucl. Part. Phys 30 S1213 (2004))
-
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::GetGrtheta()****"<<endl;
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
- TComplex fG = TComplex::One();
- Int_t fHarmonic = fFlowSelect->Har();
- Double_t fOrder = (double)(fHarmonic+1);
- Double_t fWgt = 1.;
-
- Int_t fNumberOfTracks = fFlowTracks->GetEntries();
- //cerr<<"GetGrtheta::fNumberOfTracks = "<<fNumberOfTracks<<endl;
-
- for (Int_t i=0;i<fNumberOfTracks;i++) //loop over tracks in event
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- if (fFlowSelect->SelectPart(fFlowTrack)) //if track is selected
- {
- Float_t fPhi = fFlowTrack->Phi();
- Double_t fGIm = fR * fWgt*cos(fOrder*(fPhi - fTheta));
- TComplex fGi(1., fGIm);
- fG *= fGi; //product over all tracks
- }//if
- }//loop over tracks
-
- return fG;
-
- }
-
-
-//-----------------------------------------------------------------------
-TComplex AliFlowLeeYangZerosMaker::GetDiffFlow(AliFlowSelection* fFlowSelect, Float_t fR0, Int_t theta)
-{
- // Sum for the denominator for diff. flow for the Product Generating Function for LeeYangZeros method
- // PG Eq. 9 (J. Phys. G Nucl. Part. Phys 30 S1213 (2004))
- // Also for v1 mixed harmonics: DF Eq. 5
- // It is the deriverative of Grtheta at r0 divided by Grtheta at r0
-
- if (fDebug) cout<<"****AliFlowLeeYangZerosMaker::GetGrtheta()****"<<endl;
-
- if (!fFlowSelect){
- cout<<"##### FlowLeeYangZero: FlowSelect pointer null"<<endl;
- return kFALSE;
- }
-
-
- TComplex fG = TComplex::One();
- TComplex fdGr0(0.,0.);
- Int_t fHarmonic = fFlowSelect->Har();
- Double_t fOrder = (double)(fHarmonic+1);
- Double_t fWgt = 1.;
-
- Int_t fNumberOfTracks = fFlowTracks->GetEntries();
- //cerr<<"GetGrtheta::fNumberOfTracks = "<<fNumberOfTracks<<endl;
-
- Int_t fNtheta = AliFlowLYZConstants::kTheta;
- Float_t fTheta = ((float)theta/fNtheta)*TMath::Pi()/fOrder;
-
- for (Int_t i=0;i<fNumberOfTracks;i++) //loop over tracks in event
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- if (fFlowSelect->SelectPart(fFlowTrack)) //if track is selected
- {
- Float_t fPhi = fFlowTrack->Phi();
- Double_t fCosTerm = fWgt*cos(fOrder*(fPhi - fTheta));
- //GetGr0theta
- Double_t fGIm = fR0 * fCosTerm;
- TComplex fGi(1., fGIm);
- fG *= fGi; //product over all tracks
- //GetdGr0theta
- TComplex fCosTermComplex(1., fR0*fCosTerm);
- fdGr0 +=(fCosTerm / fCosTermComplex); //sum over all tracks
- }//if
- }//loop over tracks
-
- for (Int_t i=0;i<fNumberOfTracks;i++)
- {
- fFlowTrack = (AliFlowTrack*)fFlowTracks->At(i) ; //get object at array position i
- if (fFlowSelect->SelectPart(fFlowTrack)) //if track is selected
- {
- Float_t fPhi = fFlowTrack->Phi();
- Double_t fCosTerm = cos(fOrder*(fPhi-fTheta));
- TComplex fCosTermComplex(1.,fR0*fCosTerm);
-
- TComplex fNumer = fG*fCosTerm/fCosTermComplex; //PG Eq. 9
- Float_t fEta = fFlowTrack->Eta(); //rapidity
- Float_t fPt = fFlowTrack->Pt();
- fHist2[1][theta]->Fill(fEta,fPt,fNumer);
- }//if
- }//loop over tracks
-
- TComplex fDenom = fG*fdGr0;
- return fDenom;
-
- }
-
-//-------------------------------------------------
+++ /dev/null
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-
-/* $Id$ */
-
-
-#ifndef AliFlowLeeYangZerosMaker_H
-#define AliFlowLeeYangZerosMaker_H
-
-//#include <iostream>
-//using namespace std; //required for resolving the 'cout' symbol
-
-//#include "Riostream.h"
-#include "TVector2.h" //called explicitly
-#include "AliFlowSelection.h" //only with AliFlowSelection* called ???
-#include "AliFlowTrack.h" //only with AliFlowTrack* called ???
-
-//class AliFlowTrack; //why doesn't compile?
-//class AliFlowSelection; //why doesn't compile?
-
-class AliFlowEvent;
-class AliFlowLYZHist1;
-class AliFlowLYZHist2;
-//class AliFlowLYZSummary; //class not yet finished
-
-class TH1F;
-class TH1D;
-class TProfile;
-class TProfile2D;
-class TObjArray;
-class TFile;
-class TTree;
-class TComplex;
-class Riostream;
-
-
-// Description: Maker to analyze Flow by the LeeYangZeros method
-// One needs to do two runs over the data;
-// First to calculate the integrated flow
-// and in the second to calculate the differential flow
-
-class AliFlowLeeYangZerosMaker {
-
- public:
-
- AliFlowLeeYangZerosMaker(); //default constructor
-
- AliFlowLeeYangZerosMaker(const AliFlowSelection* fFlowSelect); //custum constructor
-
- AliFlowLeeYangZerosMaker(const AliFlowLeeYangZerosMaker&); //copy constuctor (dummy)
-
- AliFlowLeeYangZerosMaker& operator=(const AliFlowLeeYangZerosMaker&); //assignment operator (dummy)
-
- virtual ~AliFlowLeeYangZerosMaker(); //destructor
-
- Bool_t Init(); //defines variables and histograms
- Bool_t Make(AliFlowEvent* fFlowEvent); //calculates variables and fills histograms
- Bool_t Finish(); //saves histograms
-
- Double_t GetQtheta(AliFlowSelection* fFlowSelect, TObjArray* fFlowTracks, Float_t fTheta);
-
- void SetFirstRun(Bool_t kt) { this->fFirstRun = kt ; }
- Bool_t GetFirstRun() const { return this->fFirstRun ; }
- void SetUseSum(Bool_t kt) { this->fUseSum = kt ; }
- Bool_t GetUseSum() const { return this->fUseSum ; }
- void SetDebug(Bool_t kt) { this->fDebug = kt ; }
- Bool_t GetDebug() const { return this->fDebug ; }
-
-
- // Output
- void SetHistFileName(TString name) { this->fHistFileName = name ; } // Sets output file name
- TString GetHistFileName() const { return this->fHistFileName ; } // Gets output file name
- TFile* GetHistFile() const { return this->fHistFile ; } // Gets output file
-
- // input for second run
- void SetFirstRunFileName(TString name) { this->firstRunFileName = name ; } // Sets input file name
- TString GetFirstRunFileName() const { return this->firstRunFileName ; } // Gets output file name
- void SetFirstRunFile(TFile* file) { this->firstRunFile = file ; } // Sets first run file
-
- private:
- Bool_t MakeControlHistograms(AliFlowEvent* fFlowEvent);
- Bool_t FillFromFlowEvent(AliFlowEvent* fFlowEvent);
- Bool_t SecondFillFromFlowEvent(AliFlowEvent* fFlowEvent);
-
- //Double_t GetQtheta(AliFlowSelection* fFlowSelect, TObjArray* fFlowTracks, Float_t fTheta) const;
- TComplex GetGrtheta(AliFlowSelection* fFlowSelect, Float_t fR, Float_t fTheta);
- TComplex GetDiffFlow(AliFlowSelection* fFlowSelect, Float_t fR, Int_t theta);
- Float_t GetR0(TH1D* fHistGtheta);
-
-
-
-#ifndef __CINT__
- TVector2 fQ; // flow vector
- TVector2 fQsum; // flow vector sum
- Float_t fQ2sum; // flow vector sum squared
- Double_t fQtheta; // flow vector projected on ref. angle theta
-
-#endif /*__CINT__*/
-
- Int_t fEventNumber; // event counter
- Int_t fMult; // multiplicity
- Int_t fNbins; // number of bins
- Float_t fTheta; // ref. angle theta
-
- AliFlowEvent* fFlowEvent; //! pointer to AliFlowEvent
- AliFlowSelection* fFlowSelect; //! pointer to AliFlowSelection
- TObjArray* fFlowTracks ; //! pointer to the TrackCollection
- AliFlowTrack* fFlowTrack ; //! pointer to the AliFlowTrack
- //AliFlowLYZSummary* fLYZSummary; //!
- TTree* fFlowTree; //!
-
- Bool_t fFirstRun ; //! flag for lyz analysis: true=first run over data, false=second run
- Bool_t fUseSum ; //! flag for lyz analysis: true=use sum gen.function, false=use product gen.function
- Bool_t fDebug ; //! flag for lyz analysis: more print statements
-
- TString fHistFileName; //!
- TFile* fHistFile; //!
- TFile* fSummaryFile; //!
- TDirectory * fFistLoop; //!
- TString firstRunFileName; //!
- TFile* firstRunFile; //!
- // for single histograms
- TH1F* fHistOrigMult; //!
- TH1F* fHistMult; //!
- TH1F* fHistQ; //!
- TH1F* fHistPt; //!
- TH1F* fHistEta; //!
- TH1F* fHistPhi; //!
- TH1F* fHistQtheta; //!
-
-
- TProfile* fHistProVthetaHar1; //!
- TProfile* fHistProVthetaHar2; //!
- TProfile* fHistProVetaHar1; //!
- TProfile* fHistProVetaHar2; //!
- TProfile* fHistProVPtHar1; //!
- TProfile* fHistProVPtHar2; //!
- TH1D* fHistVPtHar2; //!
- TProfile* fHistProVR0; //!
- TH1D* fHistVR0; //!
- TProfile* fHistProV; //!
- TProfile* fHistProR0thetaHar1; //!
- TProfile* fHistProR0thetaHar2; //!
- TProfile* fHistProReDenomHar1; //!
- TProfile* fHistProReDenomHar2; //!
- TProfile* fHistProImDenomHar1; //!
- TProfile* fHistProImDenomHar2; //!
- TProfile* fHistProReDtheta; //!
- TProfile* fHistProImDtheta; //!
-
-
- //class AliFlowLYZHist1 defines the histograms: fHistProGtheta, fHistProReGtheta, fHistProImGtheta
- AliFlowLYZHist1* fHist1[2][5]; //!
-
- //class AliFlowLYZHist1 defines the histograms: fHistProReNumer, fHistProImNumer, fHistProReNumerPt,
- //fHistProImNumerPt, fHistProReNumer2D, fHistProImNumer2D.
- AliFlowLYZHist2* fHist2[2][5]; //!
-
- ClassDef(AliFlowLeeYangZerosMaker,0) // macro for rootcint
- };
-
-
-#endif
// from CreateESDChain.C - instead of #include "CreateESDChain.C"
-TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
+TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 10, Int_t offset = 0) ;
+TChain* CreateAODChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 10, Int_t offset = 0) ;
void LookupWrite(TChain* chain, const char* target) ;
+//////////////////////////////////////////////////////////
+//
+// WARNING: THIS MACRO IS NOT TESTED FOR THE OPTION "CUM"
+//
+//////////////////////////////////////////////////////////
+
+
//SETTING THE ANALYSIS
//Flow analysis method can be:
// LYZEP = Lee Yang Zeroes Event Plane
// CUM = Cumulants
// MCEP = Flow calculated from the real MC event plane (only for simulated data)
-const TString method = "MCEP";
+const TString method = "SP";
//Type of analysis can be:
// ESD, AOD, MC, ESDMC0, ESDMC1
-const TString type = "ESD";
+const TString type = "AOD";
+//Bolean to fill/not fill the QA histograms
+Bool_t QA = kTRUE;
//SETTING THE CUTS
const Int_t mintrackrefsTPC2 = 2;
const Int_t mintrackrefsITS2 = 3;
const Int_t charge2 = 1;
-const Int_t PDG2 = 321;
+const Int_t PDG2 = 211;
const Int_t minclustersTPC2 = 50;
const Int_t maxnsigmatovertex2 = 3;
-void runAliAnalysisTaskFlow(Int_t nRuns = 10, const Char_t* dataDir="/data/alice1/kolk/TherminatorFIX", Int_t offset = 0)
-
+//void runAliAnalysisTaskFlow(Int_t nRuns = 10, const Char_t* dataDir="/data/alice1/kolk/TherminatorFIX", Int_t offset = 0)
+void runAliAnalysisTaskFlow(Int_t nRuns = 3, const Char_t* dataDir="/data/alice1/kolk/AOD", Int_t offset = 0)
{
TStopwatch timer;
timer.Start();
cerr<<"libPWG2flow.so loaded..."<<endl;
// create the TChain. CreateESDChain() is defined in CreateESDChain.C
- TChain* chain = CreateESDChain(dataDir, nRuns, offset);
- cout<<"chain ("<<chain<<")"<<endl;
+ if (type!="AOD") { TChain* chain = CreateESDChain(dataDir, nRuns, offset);
+ cout<<"chain ("<<chain<<")"<<endl; }
+ else { TChain* chain = CreateAODChain(dataDir, nRuns, offset);
+ cout<<"chain ("<<chain<<")"<<endl; }
//____________________________________________//
//Create cuts using correction framework
+ //Set TList for the QA histograms
+ if (QA) {
+ TList* qaInt = new TList();
+ TList* qaDiff = new TList();
+ }
+
//############# cuts on MC
AliCFTrackKineCuts* mcKineCuts1 = new AliCFTrackKineCuts("mcKineCuts1","MC-level kinematic cuts");
mcKineCuts1->SetPtRange(ptmin1,ptmax1);
mcKineCuts1->SetRapidityRange(ymin1,ymax1);
mcKineCuts1->SetChargeMC(charge1);
+ if (QA) { mcKineCuts1->SetQAOn(qaInt); }
AliCFTrackKineCuts* mcKineCuts2 = new AliCFTrackKineCuts("mcKineCuts2","MC-level kinematic cuts");
mcKineCuts2->SetPtRange(ptmin2,ptmax2);
mcKineCuts2->SetRapidityRange(ymin2,ymax2);
mcKineCuts2->SetChargeMC(charge2);
+ if (QA) { mcKineCuts2->SetQAOn(qaDiff); }
AliCFParticleGenCuts* mcGenCuts1 = new AliCFParticleGenCuts("mcGenCuts1","MC particle generation cuts");
mcGenCuts1->SetRequireIsPrimary();
mcGenCuts1->SetRequirePdgCode(PDG1);
+ if (QA) { mcGenCuts1->SetQAOn(qaInt); }
AliCFParticleGenCuts* mcGenCuts2 = new AliCFParticleGenCuts("mcGenCuts2","MC particle generation cuts");
mcGenCuts2->SetRequireIsPrimary();
mcGenCuts2->SetRequirePdgCode(PDG2);
-
- //############# Acceptance Cuts ????????
- AliCFAcceptanceCuts *mcAccCuts = new AliCFAcceptanceCuts("mcAccCuts","MC acceptance cuts");
- mcAccCuts->SetMinNHitITS(mintrackrefsITS1);
- mcAccCuts->SetMinNHitTPC(mintrackrefsTPC1);
+ if (QA) { mcGenCuts2->SetQAOn(qaDiff); }
+
+ //############# Acceptance Cuts
+ AliCFAcceptanceCuts *mcAccCuts1 = new AliCFAcceptanceCuts("mcAccCuts1","MC acceptance cuts");
+ mcAccCuts1->SetMinNHitITS(mintrackrefsITS1);
+ mcAccCuts1->SetMinNHitTPC(mintrackrefsTPC1);
+ if (QA) { mcAccCuts1->SetQAOn(qaInt); }
+
+ AliCFAcceptanceCuts *mcAccCuts2 = new AliCFAcceptanceCuts("mcAccCuts2","MC acceptance cuts");
+ mcAccCuts2->SetMinNHitITS(mintrackrefsITS2);
+ mcAccCuts2->SetMinNHitTPC(mintrackrefsTPC2);
+ if (QA) { mcAccCuts2->SetQAOn(qaDiff); }
//############# Rec-Level kinematic cuts
AliCFTrackKineCuts *recKineCuts1 = new AliCFTrackKineCuts("recKineCuts1","rec-level kine cuts");
recKineCuts1->SetPtRange(ptmin1,ptmax1);
recKineCuts1->SetRapidityRange(ymin1,ymax1);
recKineCuts1->SetChargeRec(charge1);
+ if (QA) { recKineCuts1->SetQAOn(qaInt); }
AliCFTrackKineCuts *recKineCuts2 = new AliCFTrackKineCuts("recKineCuts2","rec-level kine cuts");
recKineCuts2->SetPtRange(ptmin2,ptmax2);
recKineCuts2->SetRapidityRange(ymin2,ymax2);
recKineCuts2->SetChargeRec(charge2);
+ if (QA) { recKineCuts2->SetQAOn(qaDiff); }
+
+ AliCFTrackQualityCuts *recQualityCuts1 = new AliCFTrackQualityCuts("recQualityCuts1","rec-level quality cuts");
+ recQualityCuts1->SetMinNClusterTPC(minclustersTPC1);
+ recQualityCuts1->SetRequireITSRefit(kTRUE);
+ if (QA) { recQualityCuts1->SetQAOn(qaInt); }
- AliCFTrackQualityCuts *recQualityCuts = new AliCFTrackQualityCuts("recQualityCuts","rec-level quality cuts");
- recQualityCuts->SetMinNClusterTPC(minclustersTPC1);
- recQualityCuts->SetRequireITSRefit(kTRUE);
+ AliCFTrackQualityCuts *recQualityCuts2 = new AliCFTrackQualityCuts("recQualityCuts2","rec-level quality cuts");
+ recQualityCuts2->SetMinNClusterTPC(minclustersTPC2);
+ recQualityCuts2->SetRequireITSRefit(kTRUE);
+ if (QA) { recQualityCuts2->SetQAOn(qaDiff); }
- AliCFTrackIsPrimaryCuts *recIsPrimaryCuts = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts");
- recIsPrimaryCuts->SetMaxNSigmaToVertex(maxnsigmatovertex1);
+ AliCFTrackIsPrimaryCuts *recIsPrimaryCuts1 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts1","rec-level isPrimary cuts");
+ recIsPrimaryCuts1->SetMaxNSigmaToVertex(maxnsigmatovertex1);
+ if (QA) { recIsPrimaryCuts1->SetQAOn(qaInt); }
+
+ AliCFTrackIsPrimaryCuts *recIsPrimaryCuts2 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts2","rec-level isPrimary cuts");
+ recIsPrimaryCuts2->SetMaxNSigmaToVertex(maxnsigmatovertex2);
+ if (QA) { recIsPrimaryCuts2->SetQAOn(qaDiff); }
AliCFTrackCutPid* cutPID1 = new AliCFTrackCutPid("cutPID1","ESD_PID") ;
AliCFTrackCutPid* cutPID2 = new AliCFTrackCutPid("cutPID2","ESD_PID") ;
case 2212 : cutPID2->SetParticleType(AliPID::kProton , kTRUE); break;
default : printf("UNDEFINED PID\n"); break;
}
-
+ if (QA) { cutPID1->SetQAOn(qaInt);
+ cutPID2->SetQAOn(qaDiff); }
printf("CREATE MC KINE CUTS\n");
TObjArray* mcList1 = new TObjArray(0);
mcList2->AddLast(mcGenCuts2);
printf("CREATE ACCEPTANCE CUTS\n");
- TObjArray* accList = new TObjArray(0) ;
- accList->AddLast(mcAccCuts);
+ TObjArray* accList1 = new TObjArray(0) ;
+ accList1->AddLast(mcAccCuts1);
+
+ TObjArray* accList2 = new TObjArray(0) ;
+ accList2->AddLast(mcAccCuts2);
printf("CREATE RECONSTRUCTION CUTS\n");
TObjArray* recList1 = new TObjArray(0) ;
recList1->AddLast(recKineCuts1);
- recList1->AddLast(recQualityCuts);
- recList1->AddLast(recIsPrimaryCuts);
+ recList1->AddLast(recQualityCuts1);
+ recList1->AddLast(recIsPrimaryCuts1);
TObjArray* recList2 = new TObjArray(0) ;
recList2->AddLast(recKineCuts2);
- recList2->AddLast(recQualityCuts);
- recList2->AddLast(recIsPrimaryCuts);
+ recList2->AddLast(recQualityCuts2);
+ recList2->AddLast(recIsPrimaryCuts2);
printf("CREATE PID CUTS\n");
TObjArray* fPIDCutList1 = new TObjArray(0) ;
printf("CREATE INTERFACE AND CUTS\n");
AliCFManager* cfmgr1 = new AliCFManager();
- cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1);
- //cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
- cfmgr1->SetParticleCutsList(AliCFManager::kPartRecCuts,recList1);
- cfmgr1->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList1);
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1); //on MC
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList1);
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartRecCuts,recList1); //on ESD
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList1); //on ESD
AliCFManager* cfmgr2 = new AliCFManager();
cfmgr2->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList2);
- //cfmgr2->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
+ cfmgr2->SetParticleCutsList(AliCFManager::kPartAccCuts,accList2);
cfmgr2->SetParticleCutsList(AliCFManager::kPartRecCuts,recList2);
cfmgr2->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList2);
inputFileName += "_firstrun.root";
cout<<"The input file is "<<inputFileName.Data()<<endl;
TFile* fInputFile = new TFile(inputFileName.Data(),"READ");
- if(!fInputFile || fInputFile->IsZombie()) { cerr << " ERROR: NO First Run file... " << endl ; }
- cout<<"input file read..."<<endl;
+ if(!fInputFile || fInputFile->IsZombie()) {
+ cerr << " ERROR: NO First Run file... " << endl ; }
+ else {
+ TList* fInputList = (TList*)fInputFile->Get("cobj1");
+ if (!fInputList) {cout<<"list is NULL pointer!"<<endl;}
+ }
+ cout<<"input file/list read..."<<endl;
}
if (method == "LYZEP") {
cout<<"The input file is "<<inputFileName.Data()<<endl;
TFile* fInputFile = new TFile(inputFileName.Data(),"READ");
if(!fInputFile || fInputFile->IsZombie()) { cerr << " ERROR: NO First Run file... " << endl ; }
- cout<<"input file read..."<<endl;
+ else {
+ TList* fInputList = (TList*)fInputFile->Get("cobj1");
+ if (!fInputList) {cout<<"list is NULL pointer!"<<endl;}
+ }
+ cout<<"input file/list read..."<<endl;
}
//____________________________________________//
// 1st MC EP task
if (method == "SP"){
- AliAnalysisTaskScalarProduct *task1 = new AliAnalysisTaskScalarProduct("TaskScalarProduct");
- task1->SetAnalysisType(type);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kTRUE); }
+ else { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kFALSE); }
+ taskSP->SetAnalysisType(type);
+ taskSP->SetCFManager1(cfmgr1);
+ taskSP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskSP->SetQAList1(qaInt);
+ taskSP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskSP);
}
else if (method == "LYZ1"){
- AliAnalysisTaskLeeYangZeros *task1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE);
- task1->SetAnalysisType(type);
- task1->SetFirstRunLYZ(kTRUE);
- task1->SetUseSumLYZ(kTRUE);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kTRUE);}
+ else { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kFALSE);}
+ taskLYZ1->SetAnalysisType(type);
+ taskLYZ1->SetFirstRunLYZ(kTRUE);
+ taskLYZ1->SetUseSumLYZ(kTRUE);
+ taskLYZ1->SetCFManager1(cfmgr1);
+ taskLYZ1->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZ1->SetQAList1(qaInt);
+ taskLYZ1->SetQAList2(qaDiff);}
+ mgr->AddTask(taskLYZ1);
}
else if (method == "LYZ2"){
- AliAnalysisTaskLeeYangZeros *task1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE);
- task1->SetAnalysisType(type);
- task1->SetFirstRunLYZ(kFALSE);
- task1->SetUseSumLYZ(kTRUE);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kTRUE);}
+ else { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kFALSE); }
+ taskLYZ2->SetAnalysisType(type);
+ taskLYZ2->SetFirstRunLYZ(kFALSE);
+ taskLYZ2->SetUseSumLYZ(kTRUE);
+ taskLYZ2->SetCFManager1(cfmgr1);
+ taskLYZ2->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZ2->SetQAList1(qaInt);
+ taskLYZ2->SetQAList2(qaDiff); }
+ mgr->AddTask(taskLYZ2);
}
else if (method == "LYZEP"){
- AliAnalysisTaskLYZEventPlane *task1 = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane");
- task1->SetAnalysisType(type);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kTRUE); }
+ else { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kFALSE); }
+ taskLYZEP->SetAnalysisType(type);
+ taskLYZEP->SetCFManager1(cfmgr1);
+ taskLYZEP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZEP->SetQAList1(qaInt);
+ taskLYZEP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskLYZEP);
}
else if (method == "CUM"){
- AliAnalysisTaskCumulants *task1 = new AliAnalysisTaskCumulants("TaskCumulants");
- task1->SetAnalysisType(type);
- //task1->SetCFManager1(cfmgr1);
- //task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskCumulants *taskCUM = new AliAnalysisTaskCumulants("TaskCumulants",kTRUE);}
+ else { AliAnalysisTaskCumulants *taskCUM = new AliAnalysisTaskCumulants("TaskCumulants",kFALSE);}
+ taskCUM->SetAnalysisType(type);
+ taskCUM->SetCFManager1(cfmgr1);
+ taskCUM->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskCUM->SetQAList1(qaInt);
+ taskCUM->SetQAList2(qaDiff); }
+ mgr->AddTask(taskCUM);
}
else if (method == "MCEP"){
- AliAnalysisTaskMCEventPlane *task1 = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane");
- task1->SetAnalysisType(type);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
+ if (QA) { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kTRUE);}
+ else { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kFALSE);}
+ taskMCEP->SetAnalysisType(type);
+ taskMCEP->SetCFManager1(cfmgr1);
+ taskMCEP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskMCEP->SetQAList1(qaInt);
+ taskMCEP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskMCEP);
}
AliAnalysisDataContainer *coutput1 =
mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer,outputName);
+ if (QA) {
+ TString qaNameInt = "QAforInt_";
+ qaNameInt += method;
+ qaNameInt += "_";
+ qaNameInt += type;
+ qaNameInt += ".root";
+ AliAnalysisDataContainer *coutput2 =
+ mgr->CreateContainer("QAint", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameInt);
+
+ TString qaNameDiff = "QAforDiff_";
+ qaNameDiff += method;
+ qaNameDiff += "_";
+ qaNameDiff += type;
+ qaNameDiff += ".root";
+ AliAnalysisDataContainer *coutput3 =
+ mgr->CreateContainer("QAdiff", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiff);
+ }
//____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- if (method == "LYZ2" || method == "LYZEP") {
- mgr->ConnectInput(task1,1,cinput2); }
- mgr->ConnectOutput(task1,0,coutput1);
-
- if (method == "LYZ2" || method == "LYZEP"){
- cinput2->SetData(fInputFile);}
-
+
+ if (method == "SP") {
+ mgr->ConnectInput(taskSP,0,cinput1);
+ mgr->ConnectOutput(taskSP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskSP,1,coutput2);
+ mgr->ConnectOutput(taskSP,2,coutput3); }
+ }
+ else if (method == "LYZ1") {
+ mgr->ConnectInput(taskLYZ1,0,cinput1);
+ mgr->ConnectOutput(taskLYZ1,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZ1,1,coutput2);
+ mgr->ConnectOutput(taskLYZ1,2,coutput3); }
+ }
+ else if (method == "LYZ2") {
+ mgr->ConnectInput(taskLYZ2,0,cinput1);
+ mgr->ConnectInput(taskLYZ2,1,cinput2);
+ mgr->ConnectOutput(taskLYZ2,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZ2,1,coutput2);
+ mgr->ConnectOutput(taskLYZ2,2,coutput3); }
+ cinput2->SetData(fInputList);
+ }
+ else if (method == "LYZEP") {
+ mgr->ConnectInput(taskLYZEP,0,cinput1);
+ mgr->ConnectInput(taskLYZEP,1,cinput2);
+ mgr->ConnectOutput(taskLYZEP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZEP,1,coutput2);
+ mgr->ConnectOutput(taskLYZEP,2,coutput3); }
+ cinput2->SetData(fInputList);
+ }
+ else if (method == "CUM") {
+ mgr->ConnectInput(taskCUM,0,cinput1);
+ mgr->ConnectOutput(taskCUM,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskCUM,1,coutput2);
+ mgr->ConnectOutput(taskCUM,2,coutput3); }
+ }
+ else if (method == "MCEP") {
+ mgr->ConnectInput(taskMCEP,0,cinput1);
+ mgr->ConnectOutput(taskMCEP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskMCEP,1,coutput2);
+ mgr->ConnectOutput(taskMCEP,2,coutput3); }
+ }
+
if (!mgr->InitAnalysis()) return;
mgr->PrintStatus();
mgr->StartAnalysis("local",chain);
return chain;
}
+// Helper macros for creating chains
+// from: CreateESDChain.C,v 1.10 jgrosseo Exp
+
+TChain* CreateAODChain(const char* aDataDir, Int_t aRuns, Int_t offset)
+{
+ // creates chain of files in a given directory or file containing a list.
+ // In case of directory the structure is expected as:
+ // <aDataDir>/<dir0>/AliAOD.root
+ // <aDataDir>/<dir1>/AliAOD.root
+ // ...
+
+ if (!aDataDir)
+ return 0;
+
+ Long_t id, size, flags, modtime;
+ if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
+ {
+ printf("%s not found.\n", aDataDir);
+ return 0;
+ }
+
+ TChain* chain = new TChain("aodTree");
+ TChain* chaingAlice = 0;
+
+ if (flags & 2)
+ {
+ TString execDir(gSystem->pwd());
+ TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
+ TList* dirList = baseDir->GetListOfFiles();
+ Int_t nDirs = dirList->GetEntries();
+ gSystem->cd(execDir);
+
+ Int_t count = 0;
+
+ 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)
+ continue;
+
+ if (offset > 0)
+ {
+ --offset;
+ continue;
+ }
+
+ if (count++ == aRuns)
+ break;
+
+ TString presentDirName(aDataDir);
+ presentDirName += "/";
+ presentDirName += presentDir->GetName();
+ chain->Add(presentDirName + "/AliAOD.root/aodTree");
+ cerr<<presentDirName<<endl;
+ }
+
+ }
+ else
+ {
+ // Open the input stream
+ ifstream in;
+ in.open(aDataDir);
+
+ Int_t count = 0;
+
+ // Read the input list of files and add them to the chain
+ TString aodfile;
+ while(in.good()) {
+ in >> aodfile;
+ if (!aodfile.Contains("root")) continue; // protection
+
+ if (offset > 0)
+ {
+ --offset;
+ continue;
+ }
+
+ if (count++ == aRuns)
+ break;
+
+ // add aod file
+ chain->Add(aodfile);
+ }
+
+ in.close();
+ }
+
+ return chain;
+}
+
+
+
void LookupWrite(TChain* chain, const char* target)
{
// looks up the chain and writes the remaining files to the text file target
+++ /dev/null
-// from CreateESDChain.C - instead of #include "CreateESDChain.C"
-TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
-void LookupWrite(TChain* chain, const char* target) ;
-
-
-
-void runAliAnalysisTaskLYZEventPlane(Int_t nRuns = 2, TString type = "ESD", const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFIX", Int_t offset = 0)
-{
- TStopwatch timer;
- timer.Start();
-
- // include path (to find the .h files when compiling)
- gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
- gSystem->AddIncludePath("-I$ROOTSYS/include") ;
-
- // load needed libraries
- gSystem->Load("libTree.so");
- gSystem->Load("libESD.so");
- cerr<<"libESD loaded..."<<endl;
- gSystem->Load("libANALYSIS.so");
- cerr<<"libANALYSIS.so loaded..."<<endl;
- gSystem->Load("libPWG2flow.so");
-
- // gROOT->LoadMacro("AliFlowLYZConstants.cxx+");
- // gROOT->LoadMacro("AliFlowCommonConstants.cxx+");
- // gROOT->LoadMacro("AliFlowVector.cxx+");
- // gROOT->LoadMacro("AliFlowTrackSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimpleMaker.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHist.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHistResults.cxx+");
- // gROOT->LoadMacro("AliFlowLYZEventPlane.cxx+");
- // gROOT->LoadMacro("AliFlowAnalysisWithLYZEventPlane.cxx+");
- // gROOT->LoadMacro("AliAnalysisTaskLYZEventPlane.cxx+");
-
- // create the TChain. CreateESDChain() is defined in CreateESDChain.C
- TChain* chain = CreateESDChain(dataDir, nRuns, offset);
- cout<<"chain ("<<chain<<")"<<endl;
-
-
- // read the input files
- TString firstRunFileName = "outputFromLeeYangZerosAnalysis" ;
- firstRunFileName += type;
- firstRunFileName += "_firstrun.root";
- TFile* fFirstRunFile = new TFile(firstRunFileName.Data(),"READ");
- if(!fFirstRunFile || fFirstRunFile->IsZombie()) { cerr << " ERROR: NO First Run file... " << endl ; }
- TString secondRunFileName = "outputFromLeeYangZerosAnalysis" ;
- secondRunFileName += type;
- secondRunFileName += "_secondrun.root";
- TFile* fSecondRunFile = new TFile(secondRunFileName.Data(),"READ");
- if(!fSecondRunFile || fSecondRunFile->IsZombie()) { cerr << " ERROR: NO Second Run file... " << endl ; }
- cout<<"input files read..."<<endl;
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
- AliVEventHandler* esdH = new AliESDInputHandler;
- mgr->SetInputEventHandler(esdH);
- AliMCEventHandler *mc = new AliMCEventHandler();
- mgr->SetMCtruthEventHandler(mc);
- //____________________________________________//
- // 1st LYZ EP task
- AliAnalysisTaskLYZEventPlane *task1 = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane");
- task1->SetAnalysisType(type);
- mgr->AddTask(task1);
-
- // Create containers for input/output
- AliAnalysisDataContainer *cinput1 =
- mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
- AliAnalysisDataContainer *cinput2 =
- mgr->CreateContainer("cobj2",TList::Class(),AliAnalysisManager::kInputContainer);
- AliAnalysisDataContainer *cinput3 =
- mgr->CreateContainer("cobj3",TList::Class(),AliAnalysisManager::kInputContainer);
- AliAnalysisDataContainer *coutput1 =
- mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer);
-
- //____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- mgr->ConnectInput(task1,1,cinput2);
- mgr->ConnectInput(task1,2,cinput3);
- mgr->ConnectOutput(task1,0,coutput1);
-
- cinput2->SetData(fFirstRunFile);
- cinput3->SetData(fSecondRunFile);
-
- if (!mgr->InitAnalysis()) return;
- mgr->PrintStatus();
- mgr->StartAnalysis("local",chain);
-
- timer.Stop();
- timer.Print();
-}
-
-
-// Helper macros for creating chains
-// from: CreateESDChain.C,v 1.10 jgrosseo Exp
-
-TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
-{
- // creates chain of files in a given directory or file containing a list.
- // In case of directory the structure is expected as:
- // <aDataDir>/<dir0>/AliESDs.root
- // <aDataDir>/<dir1>/AliESDs.root
- // ...
-
- if (!aDataDir)
- return 0;
-
- Long_t id, size, flags, modtime;
- if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
- {
- printf("%s not found.\n", aDataDir);
- return 0;
- }
-
- TChain* chain = new TChain("esdTree");
- TChain* chaingAlice = 0;
-
- if (flags & 2)
- {
- TString execDir(gSystem->pwd());
- TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
- TList* dirList = baseDir->GetListOfFiles();
- Int_t nDirs = dirList->GetEntries();
- gSystem->cd(execDir);
-
- Int_t count = 0;
-
- 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)
- continue;
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- TString presentDirName(aDataDir);
- presentDirName += "/";
- presentDirName += presentDir->GetName();
- chain->Add(presentDirName + "/AliESDs.root/esdTree");
- cerr<<presentDirName<<endl;
- }
-
- }
- else
- {
- // Open the input stream
- ifstream in;
- in.open(aDataDir);
-
- Int_t count = 0;
-
- // Read the input list of files and add them to the chain
- TString esdfile;
- while(in.good()) {
- in >> esdfile;
- if (!esdfile.Contains("root")) continue; // protection
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- // add esd file
- chain->Add(esdfile);
- }
-
- in.close();
- }
-
- return chain;
-}
-
-void LookupWrite(TChain* chain, const char* target)
-{
- // looks up the chain and writes the remaining files to the text file target
-
- chain->Lookup();
-
- TObjArray* list = chain->GetListOfFiles();
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- ofstream outfile;
- outfile.open(target);
-
- while ((obj = iter->Next()))
- outfile << obj->GetTitle() << "#AliESDs.root" << endl;
-
- outfile.close();
-
- delete iter;
-}
-
+++ /dev/null
-// from CreateESDChain.C - instead of #include "CreateESDChain.C"
-TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
-void LookupWrite(TChain* chain, const char* target) ;
-
-
-
-void runAliAnalysisTaskLeeYangZeros(Int_t nRuns = 2, Bool_t firstrun = kTRUE, TString type = "ESD", Bool_t usesum = kTRUE, const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFIX", Int_t offset = 0)
-{
- TStopwatch timer;
- timer.Start();
-
- // include path (to find the .h files when compiling)
- gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
- gSystem->AddIncludePath("-I$ROOTSYS/include") ;
-
- // load needed libraries
- gSystem->Load("libTree.so");
- gSystem->Load("libESD.so");
- cerr<<"libESD loaded..."<<endl;
- gSystem->Load("libANALYSIS.so");
- cerr<<"libANALYSIS.so loaded..."<<endl;
- gSystem->Load("libPWG2flow.so");
-
- // gROOT->LoadMacro("AliFlowLYZConstants.cxx+");
- // gROOT->LoadMacro("AliFlowCommonConstants.cxx+");
- // gROOT->LoadMacro("AliFlowVector.cxx+");
- // gROOT->LoadMacro("AliFlowTrackSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimpleMaker.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHist.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHistResults.cxx+");
- // gROOT->LoadMacro("AliFlowLYZHist1.cxx+");
- // gROOT->LoadMacro("AliFlowLYZHist2.cxx+");
- // gROOT->LoadMacro("AliFlowAnalysisWithLeeYangZeros.cxx+");
- // gROOT->LoadMacro("AliAnalysisTaskLeeYangZeros.cxx+");
-
- // create the TChain. CreateESDChain() is defined in CreateESDChain.C
- TChain* chain = CreateESDChain(dataDir, nRuns, offset);
- cout<<"chain ("<<chain<<")"<<endl;
-
- if (!firstrun){
- // read the input files
- TString firstRunFileName = "outputFromLeeYangZerosAnalysis" ;
- firstRunFileName += type;
- firstRunFileName += "_firstrun.root";
- TFile* fFirstRunFile = new TFile(firstRunFileName.Data(),"READ");
- if(!fFirstRunFile || fFirstRunFile->IsZombie()) { cerr << " ERROR: NO First Run file... " << endl ; }
- cout<<"input files read..."<<endl;
- }
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
- AliVEventHandler* esdH = new AliESDInputHandler;
- mgr->SetInputEventHandler(esdH);
- AliMCEventHandler *mc = new AliMCEventHandler();
- mgr->SetMCtruthEventHandler(mc);
- //____________________________________________//
- // 1st LYZ task
- AliAnalysisTaskLeeYangZeros *task1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros", firstrun);
- task1->SetAnalysisType(type);
- task1->SetFirstRunLYZ(firstrun); //set to first or second run
- task1->SetUseSumLYZ(usesum); //set to sum gen.function or product gen.function
- mgr->AddTask(task1);
-
- // Create containers for input/output
- AliAnalysisDataContainer *cinput1 =
- mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
- if (!firstrun){ AliAnalysisDataContainer *cinput2 =
- mgr->CreateContainer("cobj2",TList::Class(),AliAnalysisManager::kInputContainer); }
- AliAnalysisDataContainer *coutput1 =
- mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer);
-
- //____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- if (!firstrun) { mgr->ConnectInput(task1,1,cinput2); }
- mgr->ConnectOutput(task1,0,coutput1);
-
- if (!firstrun){ cinput2->SetData(fFirstRunFile);}
-
- if (!mgr->InitAnalysis()) return;
- mgr->PrintStatus();
- mgr->StartAnalysis("local",chain);
-
- timer.Stop();
- timer.Print();
-}
-
-
-// Helper macros for creating chains
-// from: CreateESDChain.C,v 1.10 jgrosseo Exp
-
-TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
-{
- // creates chain of files in a given directory or file containing a list.
- // In case of directory the structure is expected as:
- // <aDataDir>/<dir0>/AliESDs.root
- // <aDataDir>/<dir1>/AliESDs.root
- // ...
-
- if (!aDataDir)
- return 0;
-
- Long_t id, size, flags, modtime;
- if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
- {
- printf("%s not found.\n", aDataDir);
- return 0;
- }
-
- TChain* chain = new TChain("esdTree");
- TChain* chaingAlice = 0;
-
- if (flags & 2)
- {
- TString execDir(gSystem->pwd());
- TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
- TList* dirList = baseDir->GetListOfFiles();
- Int_t nDirs = dirList->GetEntries();
- gSystem->cd(execDir);
-
- Int_t count = 0;
-
- 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)
- continue;
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- TString presentDirName(aDataDir);
- presentDirName += "/";
- presentDirName += presentDir->GetName();
- chain->Add(presentDirName + "/AliESDs.root/esdTree");
- cerr<<presentDirName<<endl;
- }
-
- }
- else
- {
- // Open the input stream
- ifstream in;
- in.open(aDataDir);
-
- Int_t count = 0;
-
- // Read the input list of files and add them to the chain
- TString esdfile;
- while(in.good()) {
- in >> esdfile;
- if (!esdfile.Contains("root")) continue; // protection
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- // add esd file
- chain->Add(esdfile);
- }
-
- in.close();
- }
-
- return chain;
-}
-
-void LookupWrite(TChain* chain, const char* target)
-{
- // looks up the chain and writes the remaining files to the text file target
-
- chain->Lookup();
-
- TObjArray* list = chain->GetListOfFiles();
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- ofstream outfile;
- outfile.open(target);
-
- while ((obj = iter->Next()))
- outfile << obj->GetTitle() << "#AliESDs.root" << endl;
-
- outfile.close();
-
- delete iter;
-}
-
+++ /dev/null
-// from CreateESDChain.C - instead of #include "CreateESDChain.C"
-TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
-void LookupWrite(TChain* chain, const char* target) ;
-
-
-
-void runAliAnalysisTaskMCEventPlane(Int_t nRuns = 2, TString type = "MC", const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFIX", Int_t offset = 0)
-
-{
- TStopwatch timer;
- timer.Start();
-
- // include path (to find the .h files when compiling)
- gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
- gSystem->AddIncludePath("-I$ROOTSYS/include") ;
-
- // load needed libraries
- gSystem->Load("libTree.so");
- gSystem->Load("libESD.so");
- cerr<<"libESD loaded..."<<endl;
- gSystem->Load("libANALYSIS.so");
- cerr<<"libANALYSIS.so loaded..."<<endl;
- gSystem->Load("libPWG2flow.so");
-
- // gROOT->LoadMacro("AliFlowCommonConstants.cxx+");
- // gROOT->LoadMacro("AliFlowVector.cxx+");
- // gROOT->LoadMacro("AliFlowTrackSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimple.cxx+");
- // gROOT->LoadMacro("AliFlowEventSimpleMaker.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHist.cxx+");
- // gROOT->LoadMacro("AliFlowCommonHistResults.cxx+");
- // gROOT->LoadMacro("AliFlowAnalysisWithMCEventPlane.cxx+");
- // gROOT->LoadMacro("AliAnalysisTaskMCEventPlane.cxx+");
-
- // create the TChain. CreateESDChain() is defined in CreateESDChain.C
- TChain* chain = CreateESDChain(dataDir, nRuns, offset);
- cout<<"chain ("<<chain<<")"<<endl;
-
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
- if (type == "ESD" || type == "ESDMC0" || type == "ESDMC1" || type == "MC" ) {
- AliVEventHandler* esdH = new AliESDInputHandler;
- mgr->SetInputEventHandler(esdH);
- }
- if (type == "AOD") {
- AliVEventHandler* aodH = new AliAODInputHandler;
- mgr->SetInputEventHandler(aodH);
- }
- AliMCEventHandler *mc = new AliMCEventHandler();
- mgr->SetMCtruthEventHandler(mc);
- //____________________________________________//
- // 1st MC EP task
- AliAnalysisTaskMCEventPlane *task1 = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane");
- task1->SetAnalysisType(type);
- mgr->AddTask(task1);
-
- // Create containers for input/output
- AliAnalysisDataContainer *cinput1 =
- mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
- AliAnalysisDataContainer *coutput1 =
- mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer);
-
- //____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- mgr->ConnectOutput(task1,0,coutput1);
-
- if (!mgr->InitAnalysis()) return;
- mgr->PrintStatus();
- mgr->StartAnalysis("local",chain);
-
- timer.Stop();
- timer.Print();
-}
-
-
-// Helper macros for creating chains
-// from: CreateESDChain.C,v 1.10 jgrosseo Exp
-
-TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
-{
- // creates chain of files in a given directory or file containing a list.
- // In case of directory the structure is expected as:
- // <aDataDir>/<dir0>/AliESDs.root
- // <aDataDir>/<dir1>/AliESDs.root
- // ...
-
- if (!aDataDir)
- return 0;
-
- Long_t id, size, flags, modtime;
- if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
- {
- printf("%s not found.\n", aDataDir);
- return 0;
- }
-
- TChain* chain = new TChain("esdTree");
- TChain* chaingAlice = 0;
-
- if (flags & 2)
- {
- TString execDir(gSystem->pwd());
- TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
- TList* dirList = baseDir->GetListOfFiles();
- Int_t nDirs = dirList->GetEntries();
- gSystem->cd(execDir);
-
- Int_t count = 0;
-
- 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)
- continue;
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- TString presentDirName(aDataDir);
- presentDirName += "/";
- presentDirName += presentDir->GetName();
- chain->Add(presentDirName + "/AliESDs.root/esdTree");
- cerr<<presentDirName<<endl;
- }
-
- }
- else
- {
- // Open the input stream
- ifstream in;
- in.open(aDataDir);
-
- Int_t count = 0;
-
- // Read the input list of files and add them to the chain
- TString esdfile;
- while(in.good()) {
- in >> esdfile;
- if (!esdfile.Contains("root")) continue; // protection
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- // add esd file
- chain->Add(esdfile);
- }
-
- in.close();
- }
-
- return chain;
-}
-
-void LookupWrite(TChain* chain, const char* target)
-{
- // looks up the chain and writes the remaining files to the text file target
-
- chain->Lookup();
-
- TObjArray* list = chain->GetListOfFiles();
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- ofstream outfile;
- outfile.open(target);
-
- while ((obj = iter->Next()))
- outfile << obj->GetTitle() << "#AliESDs.root" << endl;
-
- outfile.close();
-
- delete iter;
-}
-
+++ /dev/null
-// from CreateESDChain.C - instead of #include "CreateESDChain.C"
-TChain* CreateESDChain(const char* aDataDir = "ESDfiles.txt", Int_t aRuns = 20, Int_t offset = 0) ;
-void LookupWrite(TChain* chain, const char* target) ;
-
-
-
-void runAliAnalysisTaskScalarProduct(Int_t nRuns = 50, TString type = "ESD", const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFIX", Int_t offset = 0)
-
-{
- TStopwatch timer;
- timer.Start();
-
- // include path (to find the .h files when compiling)
- gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
- gSystem->AddIncludePath("-I$ROOTSYS/include") ;
-
- // load needed libraries
- gSystem->Load("libTree.so");
- gSystem->Load("libESD.so");
- cerr<<"libESD loaded..."<<endl;
- gSystem->Load("libANALYSIS.so");
- cerr<<"libANALYSIS.so loaded..."<<endl;
- gSystem->Load("libPWG2flow.so");
- cerr<<"libPWG2flow.so loaded..."<<endl;
-
- // create the TChain. CreateESDChain() is defined in CreateESDChain.C
- TChain* chain = CreateESDChain(dataDir, nRuns, offset);
- cout<<"chain ("<<chain<<")"<<endl;
-
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
-
- if (type == "ESD"){
- AliVEventHandler* esdH = new AliESDInputHandler;
- mgr->SetInputEventHandler(esdH); }
-
- if (type == "AOD"){
- AliVEventHandler* aodH = new AliAODInputHandler;
- mgr->SetInputEventHandler(aodH); }
-
- if (type == "MC"){
- AliMCEventHandler *mc = new AliMCEventHandler();
- mgr->SetMCtruthEventHandler(mc); }
- //____________________________________________//
- // 1st MC EP task
- AliAnalysisTaskScalarProduct *task1 = new AliAnalysisTaskScalarProduct("TaskScalarProduct");
- task1->SetAnalysisType(type);
- mgr->AddTask(task1);
-
- // Create containers for input/output
- AliAnalysisDataContainer *cinput1 =
- mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
- AliAnalysisDataContainer *coutput1 =
- // mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer);
- mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer,"outputFromScalarProductAnalysisESD.root");
-
- //____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- mgr->ConnectOutput(task1,0,coutput1);
-
- if (!mgr->InitAnalysis()) return;
- mgr->PrintStatus();
- mgr->StartAnalysis("local",chain);
-
- timer.Stop();
- timer.Print();
-}
-
-
-// Helper macros for creating chains
-// from: CreateESDChain.C,v 1.10 jgrosseo Exp
-
-TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
-{
- // creates chain of files in a given directory or file containing a list.
- // In case of directory the structure is expected as:
- // <aDataDir>/<dir0>/AliESDs.root
- // <aDataDir>/<dir1>/AliESDs.root
- // ...
-
- if (!aDataDir)
- return 0;
-
- Long_t id, size, flags, modtime;
- if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
- {
- printf("%s not found.\n", aDataDir);
- return 0;
- }
-
- TChain* chain = new TChain("esdTree");
- TChain* chaingAlice = 0;
-
- if (flags & 2)
- {
- TString execDir(gSystem->pwd());
- TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
- TList* dirList = baseDir->GetListOfFiles();
- Int_t nDirs = dirList->GetEntries();
- gSystem->cd(execDir);
-
- Int_t count = 0;
-
- 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)
- continue;
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- TString presentDirName(aDataDir);
- presentDirName += "/";
- presentDirName += presentDir->GetName();
- chain->Add(presentDirName + "/AliESDs.root/esdTree");
- cerr<<presentDirName<<endl;
- }
-
- }
- else
- {
- // Open the input stream
- ifstream in;
- in.open(aDataDir);
-
- Int_t count = 0;
-
- // Read the input list of files and add them to the chain
- TString esdfile;
- while(in.good()) {
- in >> esdfile;
- if (!esdfile.Contains("root")) continue; // protection
-
- if (offset > 0)
- {
- --offset;
- continue;
- }
-
- if (count++ == aRuns)
- break;
-
- // add esd file
- chain->Add(esdfile);
- }
-
- in.close();
- }
-
- return chain;
-}
-
-void LookupWrite(TChain* chain, const char* target)
-{
- // looks up the chain and writes the remaining files to the text file target
-
- chain->Lookup();
-
- TObjArray* list = chain->GetListOfFiles();
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- ofstream outfile;
- outfile.open(target);
-
- while ((obj = iter->Next()))
- outfile << obj->GetTitle() << "#AliESDs.root" << endl;
-
- outfile.close();
-
- delete iter;
-}
--- /dev/null
+//RUN SETTINGS
+
+//Flow analysis method can be:
+// SP = Scalar Product
+// LYZ1 = Lee Yang Zeroes first run
+// LYZ2 = Lee Yang Zeroes second run
+// LYZEP = Lee Yang Zeroes Event Plane
+// CUM = Cumulants
+// MCEP = Flow calculated from the real MC event plane (only for simulated data)
+const TString method = "SP";
+
+
+//analysis type can be ESD, AOD, MC, ESDMC0, ESDMC1
+const TString type = "ESD";
+
+//Bolean to fill/not fill the QA histograms
+Bool_t QA = kFALSE;
+
+//SETTING THE CUTS
+
+//for integrated flow
+const Double_t ptmin1 = 0.0;
+const Double_t ptmax1 = 1000.0;
+const Double_t ymin1 = -2.;
+const Double_t ymax1 = 2.;
+const Int_t mintrackrefsTPC1 = 2;
+const Int_t mintrackrefsITS1 = 3;
+const Int_t charge1 = 1;
+const Int_t PDG1 = 211;
+const Int_t minclustersTPC1 = 50;
+const Int_t maxnsigmatovertex1 = 3;
+
+//for differential flow
+const Double_t ptmin2 = 0.0;
+const Double_t ptmax2 = 1000.0;
+const Double_t ymin2 = -2.;
+const Double_t ymax2 = 2.;
+const Int_t mintrackrefsTPC2 = 2;
+const Int_t mintrackrefsITS2 = 3;
+const Int_t charge2 = 1;
+const Int_t PDG2 = 211;
+const Int_t minclustersTPC2 = 50;
+const Int_t maxnsigmatovertex2 = 3;
+
+//ESD (pp)
+//void runProofFlow(const Char_t* data="/COMMON/COMMON/LHC08c11_10TeV_0.5T", Int_t nRuns=-1, Int_t offset=0) {
+
+//AOD
+//void runProofFlow(const Char_t* data="/PWG2/pulvir/PDC08_pythia10TeV_ESD", Int_t nRuns=100, Int_t offset=0) {
+//void runProofFlow(const Char_t* data="/PWG2/mvala/pp_09_run82xxT_ESD", Int_t nRuns=100, Int_t offset=0) {
+//void runProofFlow(const Char_t* data="/PWG2/pulvir/PDC08_pythia10TeV_MC", Int_t nRuns=100, Int_t offset=0) {
+void runProofFlow(const Char_t* data="/PWG2/akisiel/Therminator_midcentral", Int_t nRuns=-1, Int_t offset=0) {
+// void runProofFlow(const Char_t* data="/PWG2/nkolk/myDataSet", Int_t nRuns=-1, Int_t offset=0) {
+
+
+ TStopwatch timer;
+ timer.Start();
+
+ printf("*** Connect to PROOF ***\n");
+ // TProof::Open("snelling@alicecaf.cern.ch");
+ TProof::Open("snelling@localhost");
+
+
+ gProof->UploadPackage("AF-v4-15");
+ gProof->EnablePackage("AF-v4-15");
+ // gProof->UploadPackage("STEERBase.par");
+ // gProof->EnablePackage("STEERBase");
+ // gProof->UploadPackage("ESD.par");
+ // gProof->EnablePackage("ESD");
+ // gProof->UploadPackage("AOD.par");
+ // gProof->EnablePackage("AOD");
+ // gProof->UploadPackage("ANALYSIS.par");
+ // gProof->EnablePackage("ANALYSIS");
+ // gProof->UploadPackage("ANALYSISalice.par");
+ // gProof->EnablePackage("ANALYSISalice");
+ // gProof->UploadPackage("PWG2AOD.par");
+ // gProof->EnablePackage("PWG2AOD");
+ gProof->UploadPackage("CORRFW.par");
+ gProof->EnablePackage("CORRFW");
+ //gProof->ClearPackage("PWG2flow");
+ gProof->UploadPackage("PWG2flow.par");
+ gProof->EnablePackage("PWG2flow");
+
+
+//____________________________________________//
+ //Create cuts using correction framework
+
+ //Set TList for the QA histograms
+ if (QA) {
+ TList* qaInt = new TList();
+ TList* qaDiff = new TList();
+ }
+
+ //############# cuts on MC
+ AliCFTrackKineCuts* mcKineCuts1 = new AliCFTrackKineCuts("mcKineCuts1","MC-level kinematic cuts");
+ mcKineCuts1->SetPtRange(ptmin1,ptmax1);
+ mcKineCuts1->SetRapidityRange(ymin1,ymax1);
+ mcKineCuts1->SetChargeMC(charge1);
+ if (QA) { mcKineCuts1->SetQAOn(qaInt); }
+
+ AliCFTrackKineCuts* mcKineCuts2 = new AliCFTrackKineCuts("mcKineCuts2","MC-level kinematic cuts");
+ mcKineCuts2->SetPtRange(ptmin2,ptmax2);
+ mcKineCuts2->SetRapidityRange(ymin2,ymax2);
+ mcKineCuts2->SetChargeMC(charge2);
+ if (QA) { mcKineCuts2->SetQAOn(qaDiff); }
+
+ AliCFParticleGenCuts* mcGenCuts1 = new AliCFParticleGenCuts("mcGenCuts1","MC particle generation cuts");
+ mcGenCuts1->SetRequireIsPrimary();
+ mcGenCuts1->SetRequirePdgCode(PDG1);
+ if (QA) { mcGenCuts1->SetQAOn(qaInt); }
+
+ AliCFParticleGenCuts* mcGenCuts2 = new AliCFParticleGenCuts("mcGenCuts2","MC particle generation cuts");
+ mcGenCuts2->SetRequireIsPrimary();
+ mcGenCuts2->SetRequirePdgCode(PDG2);
+ if (QA) { mcGenCuts2->SetQAOn(qaDiff); }
+
+ //############# Acceptance Cuts
+ AliCFAcceptanceCuts *mcAccCuts1 = new AliCFAcceptanceCuts("mcAccCuts1","MC acceptance cuts");
+ mcAccCuts1->SetMinNHitITS(mintrackrefsITS1);
+ mcAccCuts1->SetMinNHitTPC(mintrackrefsTPC1);
+ if (QA) { mcAccCuts1->SetQAOn(qaInt); }
+
+ AliCFAcceptanceCuts *mcAccCuts2 = new AliCFAcceptanceCuts("mcAccCuts2","MC acceptance cuts");
+ mcAccCuts2->SetMinNHitITS(mintrackrefsITS2);
+ mcAccCuts2->SetMinNHitTPC(mintrackrefsTPC2);
+ if (QA) { mcAccCuts2->SetQAOn(qaDiff); }
+
+ //############# Rec-Level kinematic cuts
+ AliCFTrackKineCuts *recKineCuts1 = new AliCFTrackKineCuts("recKineCuts1","rec-level kine cuts");
+ recKineCuts1->SetPtRange(ptmin1,ptmax1);
+ recKineCuts1->SetRapidityRange(ymin1,ymax1);
+ recKineCuts1->SetChargeRec(charge1);
+ if (QA) { recKineCuts1->SetQAOn(qaInt); }
+
+ AliCFTrackKineCuts *recKineCuts2 = new AliCFTrackKineCuts("recKineCuts2","rec-level kine cuts");
+ recKineCuts2->SetPtRange(ptmin2,ptmax2);
+ recKineCuts2->SetRapidityRange(ymin2,ymax2);
+ recKineCuts2->SetChargeRec(charge2);
+ if (QA) { recKineCuts2->SetQAOn(qaDiff); }
+
+ AliCFTrackQualityCuts *recQualityCuts1 = new AliCFTrackQualityCuts("recQualityCuts1","rec-level quality cuts");
+ recQualityCuts1->SetMinNClusterTPC(minclustersTPC1);
+ // recQualityCuts1->SetRequireITSRefit(kTRUE);
+ if (QA) { recQualityCuts1->SetQAOn(qaInt); }
+
+ AliCFTrackQualityCuts *recQualityCuts2 = new AliCFTrackQualityCuts("recQualityCuts2","rec-level quality cuts");
+ recQualityCuts2->SetMinNClusterTPC(minclustersTPC2);
+ // recQualityCuts2->SetRequireITSRefit(kTRUE);
+ if (QA) { recQualityCuts2->SetQAOn(qaDiff); }
+
+ AliCFTrackIsPrimaryCuts *recIsPrimaryCuts1 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts1","rec-level isPrimary cuts");
+ recIsPrimaryCuts1->SetMaxNSigmaToVertex(maxnsigmatovertex1);
+ if (QA) { recIsPrimaryCuts1->SetQAOn(qaInt); }
+
+ AliCFTrackIsPrimaryCuts *recIsPrimaryCuts2 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts2","rec-level isPrimary cuts");
+ recIsPrimaryCuts2->SetMaxNSigmaToVertex(maxnsigmatovertex2);
+ if (QA) { recIsPrimaryCuts2->SetQAOn(qaDiff); }
+
+ AliCFTrackCutPid* cutPID1 = new AliCFTrackCutPid("cutPID1","ESD_PID") ;
+ AliCFTrackCutPid* cutPID2 = new AliCFTrackCutPid("cutPID2","ESD_PID") ;
+ int n_species = AliPID::kSPECIES ;
+ Double_t* prior = new Double_t[n_species];
+
+ prior[0] = 0.0244519 ;
+ prior[1] = 0.0143988 ;
+ prior[2] = 0.805747 ;
+ prior[3] = 0.0928785 ;
+ prior[4] = 0.0625243 ;
+
+ cutPID1->SetPriors(prior);
+ cutPID1->SetProbabilityCut(0.0);
+ cutPID1->SetDetectors("TPC TOF");
+ switch(TMath::Abs(PDG1)) {
+ case 11 : cutPID1->SetParticleType(AliPID::kElectron, kTRUE); break;
+ case 13 : cutPID1->SetParticleType(AliPID::kMuon , kTRUE); break;
+ case 211 : cutPID1->SetParticleType(AliPID::kPion , kTRUE); break;
+ case 321 : cutPID1->SetParticleType(AliPID::kKaon , kTRUE); break;
+ case 2212 : cutPID1->SetParticleType(AliPID::kProton , kTRUE); break;
+ default : printf("UNDEFINED PID\n"); break;
+ }
+
+ cutPID2->SetPriors(prior);
+ cutPID2->SetProbabilityCut(0.0);
+ cutPID2->SetDetectors("TPC TOF");
+ switch(TMath::Abs(PDG2)) {
+ case 11 : cutPID2->SetParticleType(AliPID::kElectron, kTRUE); break;
+ case 13 : cutPID2->SetParticleType(AliPID::kMuon , kTRUE); break;
+ case 211 : cutPID2->SetParticleType(AliPID::kPion , kTRUE); break;
+ case 321 : cutPID2->SetParticleType(AliPID::kKaon , kTRUE); break;
+ case 2212 : cutPID2->SetParticleType(AliPID::kProton , kTRUE); break;
+ default : printf("UNDEFINED PID\n"); break;
+ }
+ if (QA) { cutPID1->SetQAOn(qaInt);
+ cutPID2->SetQAOn(qaDiff); }
+
+ printf("CREATE MC KINE CUTS\n");
+ TObjArray* mcList1 = new TObjArray(0);
+ mcList1->AddLast(mcKineCuts1);
+ mcList1->AddLast(mcGenCuts1);
+
+ TObjArray* mcList2 = new TObjArray(0);
+ mcList2->AddLast(mcKineCuts2);
+ mcList2->AddLast(mcGenCuts2);
+
+ printf("CREATE ACCEPTANCE CUTS\n");
+ TObjArray* accList1 = new TObjArray(0) ;
+ accList1->AddLast(mcAccCuts1);
+
+ TObjArray* accList2 = new TObjArray(0) ;
+ accList2->AddLast(mcAccCuts2);
+
+ printf("CREATE RECONSTRUCTION CUTS\n");
+ TObjArray* recList1 = new TObjArray(0) ;
+ recList1->AddLast(recKineCuts1);
+ recList1->AddLast(recQualityCuts1);
+ recList1->AddLast(recIsPrimaryCuts1);
+
+ TObjArray* recList2 = new TObjArray(0) ;
+ recList2->AddLast(recKineCuts2);
+ recList2->AddLast(recQualityCuts2);
+ recList2->AddLast(recIsPrimaryCuts2);
+
+ printf("CREATE PID CUTS\n");
+ TObjArray* fPIDCutList1 = new TObjArray(0) ;
+ fPIDCutList1->AddLast(cutPID1);
+
+ TObjArray* fPIDCutList2 = new TObjArray(0) ;
+ fPIDCutList2->AddLast(cutPID2);
+
+ printf("CREATE INTERFACE AND CUTS\n");
+ AliCFManager* cfmgr1 = new AliCFManager();
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1);
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList1);
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartRecCuts,recList1);
+ cfmgr1->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList1);
+
+ AliCFManager* cfmgr2 = new AliCFManager();
+ cfmgr2->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList2);
+ cfmgr2->SetParticleCutsList(AliCFManager::kPartAccCuts,accList2);
+ cfmgr2->SetParticleCutsList(AliCFManager::kPartRecCuts,recList2);
+ cfmgr2->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList2);
+
+ if (method == "LYZ2"){
+ // read the input file from the first run
+ TString inputFileName = "outputLYZ1analysis" ;
+ inputFileName += type;
+ inputFileName += "_firstrun.root";
+ cout<<"The input file is "<<inputFileName.Data()<<endl;
+ TFile* fInputFile = new TFile(inputFileName.Data(),"READ");
+ if(!fInputFile || fInputFile->IsZombie()) {
+ cerr << " ERROR: NO First Run file... " << endl ; }
+ else {
+ TList* fInputList = (TList*)fInputFile->Get("cobj1");
+ if (!fInputList) {cout<<"list is NULL pointer!"<<endl;}
+ }
+ cout<<"input file/list read..."<<endl;
+ }
+
+ if (method == "LYZEP") {
+ // read the input file from the second LYZ run
+ TString inputFileName = "outputLYZ2analysis" ;
+ inputFileName += type;
+ inputFileName += "_secondrun.root";
+ cout<<"The input file is "<<inputFileName.Data()<<endl;
+ TFile* fInputFile = new TFile(inputFileName.Data(),"READ");
+ if(!fInputFile || fInputFile->IsZombie()) {
+ cerr << " ERROR: NO First Run file... " << endl ; }
+ else {
+ TList* fInputList = (TList*)fInputFile->Get("cobj1");
+ if (!fInputList) {cout<<"list is NULL pointer!"<<endl;}
+ }
+ cout<<"input file/list read..."<<endl;
+ }
+
+ //____________________________________________//
+ // Make the analysis manager
+ AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
+
+ if (type == "ESD"){
+ AliVEventHandler* esdH = new AliESDInputHandler;
+ mgr->SetInputEventHandler(esdH);
+ if (method == "MCEP") {
+ AliMCEventHandler *mc = new AliMCEventHandler();
+ mgr->SetMCtruthEventHandler(mc);} }
+
+ if (type == "AOD"){
+ AliVEventHandler* aodH = new AliAODInputHandler;
+ mgr->SetInputEventHandler(aodH);
+ if (method == "MCEP") {
+ AliMCEventHandler *mc = new AliMCEventHandler();
+ mgr->SetMCtruthEventHandler(mc);} }
+
+ if (type == "MC" || type == "ESDMC0" || type == "ESDMC1"){
+ AliVEventHandler* esdH = new AliESDInputHandler;
+ mgr->SetInputEventHandler(esdH);
+
+ AliMCEventHandler *mc = new AliMCEventHandler();
+ mgr->SetMCtruthEventHandler(mc); }
+
+ //____________________________________________//
+ // tasks
+ if (method == "SP"){
+ if (QA) { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kTRUE); }
+ else { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kFALSE); }
+ taskSP->SetAnalysisType(type);
+ taskSP->SetCFManager1(cfmgr1);
+ taskSP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskSP->SetQAList1(qaInt);
+ taskSP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskSP);
+ }
+ else if (method == "LYZ1"){
+ if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kTRUE);}
+ else { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kFALSE);}
+ taskLYZ1->SetAnalysisType(type);
+ taskLYZ1->SetFirstRunLYZ(kTRUE);
+ taskLYZ1->SetUseSumLYZ(kTRUE);
+ taskLYZ1->SetCFManager1(cfmgr1);
+ taskLYZ1->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZ1->SetQAList1(qaInt);
+ taskLYZ1->SetQAList2(qaDiff);}
+ mgr->AddTask(taskLYZ1);
+ }
+ else if (method == "LYZ2"){
+ if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kTRUE);}
+ else { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kFALSE); }
+ taskLYZ2->SetAnalysisType(type);
+ taskLYZ2->SetFirstRunLYZ(kFALSE);
+ taskLYZ2->SetUseSumLYZ(kTRUE);
+ taskLYZ2->SetCFManager1(cfmgr1);
+ taskLYZ2->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZ2->SetQAList1(qaInt);
+ taskLYZ2->SetQAList2(qaDiff); }
+ mgr->AddTask(taskLYZ2);
+ }
+ else if (method == "LYZEP"){
+ if (QA) { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kTRUE); }
+ else { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kFALSE); }
+ taskLYZEP->SetAnalysisType(type);
+ taskLYZEP->SetCFManager1(cfmgr1);
+ taskLYZEP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskLYZEP->SetQAList1(qaInt);
+ taskLYZEP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskLYZEP);
+ }
+ else if (method == "CUM"){
+ if (QA) { AliAnalysisTaskCumulants *taskCUM = new AliAnalysisTaskCumulants("TaskCumulants",kTRUE);}
+ else { AliAnalysisTaskCumulants *taskCUM = new AliAnalysisTaskCumulants("TaskCumulants",kFALSE);}
+ taskCUM->SetAnalysisType(type);
+ taskCUM->SetCFManager1(cfmgr1);
+ taskCUM->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskCUM->SetQAList1(qaInt);
+ taskCUM->SetQAList2(qaDiff); }
+ mgr->AddTask(taskCUM);
+ }
+ else if (method == "MCEP"){
+ if (QA) { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kTRUE);}
+ else { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kFALSE);}
+ taskMCEP->SetAnalysisType(type);
+ taskMCEP->SetCFManager1(cfmgr1);
+ taskMCEP->SetCFManager2(cfmgr2);
+ if (QA) {
+ taskMCEP->SetQAList1(qaInt);
+ taskMCEP->SetQAList2(qaDiff); }
+ mgr->AddTask(taskMCEP);
+ }
+
+
+ // Create containers for input/output
+ AliAnalysisDataContainer *cinput1 =
+ mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
+
+ if (method == "LYZ2" || method == "LYZEP"){
+ AliAnalysisDataContainer *cinput2 =
+ mgr->CreateContainer("cobj2",TList::Class(),AliAnalysisManager::kInputContainer); }
+
+
+ TString outputName = "output";
+ outputName+= method;
+ outputName+= "analysis";
+ outputName+= type;
+ if (method == "LYZ1") outputName+= "_firstrun";
+ if (method == "LYZ2") outputName+= "_secondrun";
+ outputName+= ".root";
+ AliAnalysisDataContainer *coutput1 =
+ mgr->CreateContainer("cobj1", TList::Class(),AliAnalysisManager::kOutputContainer,outputName);
+
+ if (QA) {
+ TString qaNameInt = "QAforInt_";
+ qaNameInt += method;
+ qaNameInt += "_";
+ qaNameInt += type;
+ qaNameInt += ".root";
+ AliAnalysisDataContainer *coutput2 =
+ mgr->CreateContainer("QAint", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameInt);
+
+ TString qaNameDiff = "QAforDiff_";
+ qaNameDiff += method;
+ qaNameDiff += "_";
+ qaNameDiff += type;
+ qaNameDiff += ".root";
+ AliAnalysisDataContainer *coutput3 =
+ mgr->CreateContainer("QAdiff", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiff);
+ }
+
+ //____________________________________________//
+
+ if (method == "SP") {
+ mgr->ConnectInput(taskSP,0,cinput1);
+ mgr->ConnectOutput(taskSP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskSP,1,coutput2);
+ mgr->ConnectOutput(taskSP,2,coutput3); }
+ }
+ else if (method == "LYZ1") {
+ mgr->ConnectInput(taskLYZ1,0,cinput1);
+ mgr->ConnectOutput(taskLYZ1,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZ1,1,coutput2);
+ mgr->ConnectOutput(taskLYZ1,2,coutput3); }
+ }
+ else if (method == "LYZ2") {
+ mgr->ConnectInput(taskLYZ2,0,cinput1);
+ mgr->ConnectInput(taskLYZ2,1,cinput2);
+ mgr->ConnectOutput(taskLYZ2,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZ2,1,coutput2);
+ mgr->ConnectOutput(taskLYZ2,2,coutput3); }
+ cinput2->SetData(fInputList);
+ }
+ else if (method == "LYZEP") {
+ mgr->ConnectInput(taskLYZEP,0,cinput1);
+ mgr->ConnectInput(taskLYZEP,1,cinput2);
+ mgr->ConnectOutput(taskLYZEP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskLYZEP,1,coutput2);
+ mgr->ConnectOutput(taskLYZEP,2,coutput3); }
+ cinput2->SetData(fInputList);
+ }
+ else if (method == "CUM") {
+ mgr->ConnectInput(taskCUM,0,cinput1);
+ mgr->ConnectOutput(taskCUM,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskCUM,1,coutput2);
+ mgr->ConnectOutput(taskCUM,2,coutput3); }
+ }
+ else if (method == "MCEP") {
+ mgr->ConnectInput(taskMCEP,0,cinput1);
+ mgr->ConnectOutput(taskMCEP,0,coutput1);
+ if (QA) { mgr->ConnectOutput(taskMCEP,1,coutput2);
+ mgr->ConnectOutput(taskMCEP,2,coutput3); }
+ }
+
+ if (!mgr->InitAnalysis()) return;
+ mgr->PrintStatus();
+ // old way with a chain
+ // mgr->StartAnalysis("proof",chain);
+ mgr->StartAnalysis("proof",data,nRuns,offset);
+
+ timer.Stop();
+ timer.Print();
+}
+++ /dev/null
-//RUN SETTINGS
-//analysis type can be ESD, AOD, MC, ESDMC0, ESDMC1
-const TString type = "ESD";
-
-
-//SETTING THE CUTS
-
-//for integrated flow
-const Double_t ptmin1 = 0.0;
-const Double_t ptmax1 = 1000.0;
-const Double_t ymin1 = -2.;
-const Double_t ymax1 = 2.;
-const Int_t mintrackrefsTPC1 = 2;
-const Int_t mintrackrefsITS1 = 3;
-const Int_t charge1 = 1;
-const Int_t PDG1 = 211;
-const Int_t minclustersTPC1 = 50;
-const Int_t maxnsigmatovertex1 = 3;
-
-//for differential flow
-const Double_t ptmin2 = 0.0;
-const Double_t ptmax2 = 1000.0;
-const Double_t ymin2 = -2.;
-const Double_t ymax2 = 2.;
-const Int_t mintrackrefsTPC2 = 2;
-const Int_t mintrackrefsITS2 = 3;
-const Int_t charge2 = 1;
-const Int_t PDG2 = 321;
-const Int_t minclustersTPC2 = 50;
-const Int_t maxnsigmatovertex2 = 3;
-
-
-//for AOD testing:
-//void runProofScalarProduct(const Char_t* data="/proofteam/mgheata/lhc08b1_30000x", Int_t nRuns=2000, Int_t offset=0) {
-
-//----------First physics data:
-void runProofScalarProduct(const Char_t* data="/PWG0/COMMON/run30000X_10TeV_0.5T", Int_t nRuns=2000, Int_t offset=0) {
- //void runProofScalarProduct(const Char_t* data="/PWG0/COMMON/run31000X_0.9TeV_0.5T", Int_t nRuns=200, Int_t offset=0) {
-
-//----------Test data:
- //void runProofScalarProduct(const Char_t* data="/COMMON/COMMON/sim_160000_esd", Int_t nRuns=20, Int_t offset=0) {
- //void runProofScalarProduct(const Char_t* data="/COMMON/COMMON/run82XX_part3", Int_t nRuns=20, Int_t offset=0) {
- //void runProofScalarProduct(const Char_t* data="/PWG2/COMMON/run82XX_test5", Int_t nRuns=20, Int_t offset=0) {
-
-//----------Therminator data:
- //void runProofScalarProduct(const Char_t* data="/PWG2/akisiel/LHC500C2030", Int_t nRuns=100, Int_t offset=0) {
-
- TStopwatch timer;
- timer.Start();
-
- printf("*** Connect to PROOF ***\n");
- // TProof::Open("nkolk@lxb6046.cern.ch");
- // TProof::Open("snelling@lxb6046.cern.ch");
- TProof::Open("snelling@localhost");
-
- gProof->UploadPackage("STEERBase.par");
- gProof->EnablePackage("STEERBase");
- gProof->UploadPackage("ESD.par");
- gProof->EnablePackage("ESD");
- gProof->UploadPackage("AOD.par");
- gProof->EnablePackage("AOD");
- gProof->UploadPackage("ANALYSIS.par");
- gProof->EnablePackage("ANALYSIS");
- gProof->UploadPackage("ANALYSISalice.par");
- gProof->EnablePackage("ANALYSISalice");
- gProof->UploadPackage("PWG2AOD.par");
- gProof->EnablePackage("PWG2AOD");
- gProof->UploadPackage("CORRFW.par");
- gProof->EnablePackage("CORRFW");
- gProof->ClearPackage("PWG2flow");
- gProof->UploadPackage("PWG2flow.par");
- gProof->EnablePackage("PWG2flow");
-
-
-//____________________________________________//
- //Create cuts using correction framework
-
- //############# cuts on MC
- AliCFTrackKineCuts* mcKineCuts1 = new AliCFTrackKineCuts("mcKineCuts1","MC-level kinematic cuts");
- mcKineCuts1->SetPtRange(ptmin1,ptmax1);
- mcKineCuts1->SetRapidityRange(ymin1,ymax1);
- mcKineCuts1->SetChargeMC(charge1);
-
- AliCFTrackKineCuts* mcKineCuts2 = new AliCFTrackKineCuts("mcKineCuts2","MC-level kinematic cuts");
- mcKineCuts2->SetPtRange(ptmin2,ptmax2);
- mcKineCuts2->SetRapidityRange(ymin2,ymax2);
- mcKineCuts2->SetChargeMC(charge2);
-
- AliCFParticleGenCuts* mcGenCuts1 = new AliCFParticleGenCuts("mcGenCuts1","MC particle generation cuts");
- mcGenCuts1->SetRequireIsPrimary();
- mcGenCuts1->SetRequirePdgCode(PDG1);
-
- AliCFParticleGenCuts* mcGenCuts2 = new AliCFParticleGenCuts("mcGenCuts2","MC particle generation cuts");
- mcGenCuts2->SetRequireIsPrimary();
- mcGenCuts2->SetRequirePdgCode(PDG2);
-
- //############# Acceptance Cuts ????????
- AliCFAcceptanceCuts *mcAccCuts = new AliCFAcceptanceCuts("mcAccCuts","MC acceptance cuts");
- mcAccCuts->SetMinNHitITS(mintrackrefsITS1);
- mcAccCuts->SetMinNHitTPC(mintrackrefsTPC1);
-
- //############# Rec-Level kinematic cuts
- AliCFTrackKineCuts *recKineCuts1 = new AliCFTrackKineCuts("recKineCuts1","rec-level kine cuts");
- recKineCuts1->SetPtRange(ptmin1,ptmax1);
- recKineCuts1->SetRapidityRange(ymin1,ymax1);
- recKineCuts1->SetChargeRec(charge1);
-
- AliCFTrackKineCuts *recKineCuts2 = new AliCFTrackKineCuts("recKineCuts2","rec-level kine cuts");
- recKineCuts2->SetPtRange(ptmin2,ptmax2);
- recKineCuts2->SetRapidityRange(ymin2,ymax2);
- recKineCuts2->SetChargeRec(charge2);
-
- AliCFTrackQualityCuts *recQualityCuts = new AliCFTrackQualityCuts("recQualityCuts","rec-level quality cuts");
- recQualityCuts->SetMinNClusterTPC(minclustersTPC1);
- recQualityCuts->SetRequireITSRefit(kTRUE);
-
- AliCFTrackIsPrimaryCuts *recIsPrimaryCuts = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts");
- recIsPrimaryCuts->SetMaxNSigmaToVertex(maxnsigmatovertex1);
-
- AliCFTrackCutPid* cutPID1 = new AliCFTrackCutPid("cutPID1","ESD_PID") ;
- AliCFTrackCutPid* cutPID2 = new AliCFTrackCutPid("cutPID2","ESD_PID") ;
- int n_species = AliPID::kSPECIES ;
- Double_t* prior = new Double_t[n_species];
-
- prior[0] = 0.0244519 ;
- prior[1] = 0.0143988 ;
- prior[2] = 0.805747 ;
- prior[3] = 0.0928785 ;
- prior[4] = 0.0625243 ;
-
- cutPID1->SetPriors(prior);
- cutPID1->SetProbabilityCut(0.0);
- cutPID1->SetDetectors("TPC TOF");
- switch(TMath::Abs(PDG1)) {
- case 11 : cutPID1->SetParticleType(AliPID::kElectron, kTRUE); break;
- case 13 : cutPID1->SetParticleType(AliPID::kMuon , kTRUE); break;
- case 211 : cutPID1->SetParticleType(AliPID::kPion , kTRUE); break;
- case 321 : cutPID1->SetParticleType(AliPID::kKaon , kTRUE); break;
- case 2212 : cutPID1->SetParticleType(AliPID::kProton , kTRUE); break;
- default : printf("UNDEFINED PID\n"); break;
- }
-
- cutPID2->SetPriors(prior);
- cutPID2->SetProbabilityCut(0.0);
- cutPID2->SetDetectors("TPC TOF");
- switch(TMath::Abs(PDG2)) {
- case 11 : cutPID2->SetParticleType(AliPID::kElectron, kTRUE); break;
- case 13 : cutPID2->SetParticleType(AliPID::kMuon , kTRUE); break;
- case 211 : cutPID2->SetParticleType(AliPID::kPion , kTRUE); break;
- case 321 : cutPID2->SetParticleType(AliPID::kKaon , kTRUE); break;
- case 2212 : cutPID2->SetParticleType(AliPID::kProton , kTRUE); break;
- default : printf("UNDEFINED PID\n"); break;
- }
-
- printf("CREATE MC KINE CUTS\n");
- TObjArray* mcList1 = new TObjArray(0);
- mcList1->AddLast(mcKineCuts1);
- mcList1->AddLast(mcGenCuts1);
-
- TObjArray* mcList2 = new TObjArray(0);
- mcList2->AddLast(mcKineCuts2);
- mcList2->AddLast(mcGenCuts2);
-
- printf("CREATE ACCEPTANCE CUTS\n");
- TObjArray* accList = new TObjArray(0) ;
- accList->AddLast(mcAccCuts);
-
- printf("CREATE RECONSTRUCTION CUTS\n");
- TObjArray* recList1 = new TObjArray(0) ;
- recList1->AddLast(recKineCuts1);
- recList1->AddLast(recQualityCuts);
- recList1->AddLast(recIsPrimaryCuts);
-
- TObjArray* recList2 = new TObjArray(0) ;
- recList2->AddLast(recKineCuts2);
- recList2->AddLast(recQualityCuts);
- recList2->AddLast(recIsPrimaryCuts);
-
- printf("CREATE PID CUTS\n");
- TObjArray* fPIDCutList1 = new TObjArray(0) ;
- fPIDCutList1->AddLast(cutPID1);
-
- TObjArray* fPIDCutList2 = new TObjArray(0) ;
- fPIDCutList2->AddLast(cutPID2);
-
- printf("CREATE INTERFACE AND CUTS\n");
- AliCFManager* cfmgr1 = new AliCFManager();
- cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1);
- //cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
- cfmgr1->SetParticleCutsList(AliCFManager::kPartRecCuts,recList1);
- cfmgr1->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList1);
-
- AliCFManager* cfmgr2 = new AliCFManager();
- cfmgr2->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList2);
- //cfmgr2->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
- cfmgr2->SetParticleCutsList(AliCFManager::kPartRecCuts,recList2);
- cfmgr2->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList2);
-
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
-
- if (type == "ESD"){
- AliVEventHandler* esdH = new AliESDInputHandler;
- //esdH->SetInactiveBranches("FMD CaloCluster"); //needed?
- mgr->SetInputEventHandler(esdH); }
-
- if (type == "AOD"){
- AliVEventHandler* aodH = new AliAODInputHandler;
- mgr->SetInputEventHandler(aodH); }
-
- if (type == "MC" || type == "ESDMC0" || type == "ESDMC1"){
- AliVEventHandler* esdH = new AliESDInputHandler;
- mgr->SetInputEventHandler(esdH);
-
- AliMCEventHandler *mc = new AliMCEventHandler();
- mgr->SetMCtruthEventHandler(mc); }
-
- //____________________________________________//
- // 1st Pt task
- AliAnalysisTaskScalarProduct *task1 = new AliAnalysisTaskScalarProduct("TaskScalarProduct");
- task1->SetAnalysisType(type);
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
- mgr->AddTask(task1);
-
-
- // Create containers for input/output
- AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
- TString outputName = "outputFromScalarProductAnalysis";
- outputName+=type;
- outputName+=".root";
- AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("clist1", TList::Class(),AliAnalysisManager::kOutputContainer,outputName);
-
- //____________________________________________//
- mgr->ConnectInput(task1,0,cinput1);
- mgr->ConnectOutput(task1,0,coutput1);
-
- if (!mgr->InitAnalysis()) return;
- mgr->PrintStatus();
- // old way with a chain
- // mgr->StartAnalysis("proof",chain);
- mgr->StartAnalysis("proof",data,nRuns,offset);
-
- timer.Stop();
- timer.Print();
-}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff