* provided "as is" without express or implied warranty. *
**************************************************************************/
-#include "Riostream.h" //needed as include
+/**************************************
+ * analysis task for cumulant method *
+ * *
+ * authors: Naomi van der Kolk *
+ * (kolk@nikhef.nl) *
+ * Raimond Snellings *
+ * (snelling@nikhef.nl) *
+ * Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ * ***********************************/
+
+#include "Riostream.h"
#include "TChain.h"
#include "TTree.h"
#include "TFile.h"
#include "TList.h"
#include "TH1.h"
-#include "TH3D.h"
#include "TProfile.h"
#include "TProfile2D.h"
#include "TProfile3D.h"
-#include "AliCumulantsFunctions.h"
-
#include "AliAnalysisTask.h"
#include "AliAnalysisDataSlot.h"
#include "AliAnalysisDataContainer.h"
-
-//class AliAnalysisTask;
#include "AliAnalysisManager.h"
#include "AliESDEvent.h"
#include "AliFlowCumuConstants.h"
#include "AliFlowCommonConstants.h"
#include "AliFlowCommonHistResults.h"
-
#include "AliCumulantsFunctions.h"
-// AliAnalysisTaskCumulants:
-// analysis task for
-// Cumulant method
-// with many authors (N.K. R.S. A.B.)
-// who do something
-
ClassImp(AliAnalysisTaskCumulants)
-//________________________________________________________________________
-AliAnalysisTaskCumulants::AliAnalysisTaskCumulants(const char *name) :
- AliAnalysisTask(name, ""),
- fESD(NULL),
- fAOD(NULL),
- fMyCumuAnalysis(NULL),
- fEventMaker(NULL),
- fAnalysisType("ESD"),
- fCFManager1(NULL),
- fCFManager2(NULL),
- fListHistos(NULL)
+//================================================================================================================
+
+AliAnalysisTaskCumulants::AliAnalysisTaskCumulants(const char *name):
+ AliAnalysisTask(name,""),
+ fESD(NULL),
+ fAOD(NULL),
+ fCA(NULL),//Cumulant Analysis (CA) object
+ fEventMaker(NULL),
+ fAnalysisType("ESD"),
+ fCFManager1(NULL),
+ fCFManager2(NULL),
+ fListHistos(NULL)
{
- // Constructor
- cout<<"AliAnalysisTaskCumulants::AliAnalysisTaskCumulants(const char *name)"<<endl;
-
- // Define input and output slots here
- // Input slot #0 works with a TChain
- DefineInput(0, TChain::Class());
-
- // Output slot #0 writes into a TList container
- DefineOutput(0, TList::Class());
+ //constructor
+ cout<<"AliAnalysisTaskCumulants::AliAnalysisTaskCumulants(const char *name)"<<endl;
+ //input and output slots:
+ //input slot #0 works with a TChain
+ DefineInput(0, TChain::Class());
+ //output slot #0 writes into a TList container
+ DefineOutput(0, TList::Class());
}
-//________________________________________________________________________
-AliAnalysisTaskCumulants::AliAnalysisTaskCumulants() :
- fESD(NULL),
- fAOD(NULL),
- fMyCumuAnalysis(NULL),
- fEventMaker(NULL),
- fAnalysisType("ESD"),
- fCFManager1(NULL),
- fCFManager2(NULL),
- fListHistos(NULL)
+
+AliAnalysisTaskCumulants::AliAnalysisTaskCumulants():
+ fESD(NULL),
+ fAOD(NULL),
+ fCA(NULL),//Cumulant Analysis (CA) object
+ fEventMaker(NULL),
+ fAnalysisType("ESD"),
+ fCFManager1(NULL),
+ fCFManager2(NULL),
+ fListHistos(NULL)
{
- // Constructor
- cout<<"AliAnalysisTaskCumulants::AliAnalysisTaskCumulants(const char *name)"<<endl;
- }
+ //dummy constructor
+ cout<<"AliAnalysisTaskCumulants::AliAnalysisTaskCumulants()"<<endl;
+}
+
+//================================================================================================================
-//________________________________________________________________________
void AliAnalysisTaskCumulants::ConnectInputData(Option_t *)
{
- // Connect ESD or AOD here
- // Called once
- cout<<"AliAnalysisTaskCumulants::ConnectInputData(Option_t *)"<<endl;
-
- TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
- if (!tree) {
- Printf("ERROR: Could not read chain from input slot 0");
- }
- else {
- // Disable all branches and enable only the needed ones
- if (fAnalysisType == "MC") {
- // we want to process only MC
+ //connect ESD or AOD (called once)
+ cout<<"AliAnalysisTaskCumulants::ConnectInputData(Option_t *)"<<endl;
+
+ TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
+ if (!tree)
+ {
+ Printf("ERROR: Could not read chain from input slot 0");
+ }
+ else
+ {
+ //disable all branches and enable only the needed ones
+ if (fAnalysisType == "MC") {
+ // we want to process only MC
tree->SetBranchStatus("*", kFALSE);
AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
}
}
-//________________________________________________________________________
+//================================================================================================================
+
void AliAnalysisTaskCumulants::CreateOutputObjects()
{
- // Called at every worker node to initialize
- cout<<"AliAnalysisTaskCumulants::CreateOutputObjects()"<<endl;
+ //called at every worker node to initialize
+ cout<<"AliAnalysisTaskCumulants::CreateOutputObjects()"<<endl;
- if (!(fAnalysisType == "AOD" || fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" || fAnalysisType == "MC")) {
- cout<<"WRONG ANALYSIS TYPE! only ESD, ESDMC0, ESDMC1, AOD and MC are allowed."<<endl;
- exit(1);
- }
- //event maker
- fEventMaker = new AliFlowEventSimpleMaker();
- //analyser
- fMyCumuAnalysis = new AliFlowAnalysisWithCumulants() ;
-
- fMyCumuAnalysis->CreateOutputObjects();
-
- if (fMyCumuAnalysis->GetHistList()) {
- // fSP->GetHistList()->Print();
- fListHistos = fMyCumuAnalysis->GetHistList();
- // fListHistos->Print();
- }
- else {Printf("ERROR: Could not retrieve histogram list"); }
+ //OpenFile(0);
+
+
+ if(!(fAnalysisType == "AOD" || fAnalysisType == "ESD" || fAnalysisType == "ESDMC0" || fAnalysisType == "ESDMC1" || fAnalysisType == "MC"))
+ {
+ cout<<"WRONG ANALYSIS TYPE! only ESD, ESDMC0, ESDMC1, AOD and MC are allowed."<<endl;
+ exit(1);
+ }
+
+ //event maker
+ fEventMaker = new AliFlowEventSimpleMaker();
+
+ //analyser
+ fCA = new AliFlowAnalysisWithCumulants();
+ fCA->CreateOutputObjects();
+
+ if(fCA->GetHistList())
+ {
+ fListHistos = fCA->GetHistList();
+ //fListHistos->Print();
}
+ else
+ {
+ Printf("ERROR: Could not retrieve histogram list");
+ }
+
+ //PostData(0,fListHistos);
+
+}
+
+//================================================================================================================
-//________________________________________________________________________
void AliAnalysisTaskCumulants::Exec(Option_t *)
{
- // Main loop
- // Called for each event
- if (fAnalysisType == "MC") {
+ //main loop (called for each event)
+ if (fAnalysisType == "MC") {
// Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
// This handler can return the current MC event
//cumulant analysis
AliFlowEventSimple* fEvent = fEventMaker->FillTracks(mcEvent,fCFManager1,fCFManager2);
- fMyCumuAnalysis->Make(fEvent);
+ fCA->Make(fEvent);
delete fEvent;
}
else if (fAnalysisType == "ESD") {
}
Printf("There are %d tracks in this event", fESD->GetNumberOfTracks());
- //Cumulant analysis
+ //cumulant analysis
AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fESD);
- fMyCumuAnalysis->Make(fEvent);
+ fCA->Make(fEvent);
delete fEvent;
}
else if (fAnalysisType == "ESDMC0") {
fCFManager1->SetEventInfo(mcEvent);
fCFManager2->SetEventInfo(mcEvent);
- //Cumulant analysis
+ //cumulant analysis
AliFlowEventSimple* fEvent=NULL;
if (fAnalysisType == "ESDMC0") {
fEvent = fEventMaker->FillTracks(fESD, mcEvent, fCFManager1, fCFManager2, 0); //0 = kine from ESD, 1 = kine from MC
} else if (fAnalysisType == "ESDMC1") {
fEvent = fEventMaker->FillTracks(fESD, mcEvent, fCFManager1, fCFManager2, 1); //0 = kine from ESD, 1 = kine from MC
}
- fMyCumuAnalysis->Make(fEvent);
+ fCA->Make(fEvent);
delete fEvent;
//delete mcEvent;
}
// analysis
//For the moment don't use CF //AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fAOD,fCFManager1,fCFManager2);
AliFlowEventSimple* fEvent = fEventMaker->FillTracks(fAOD);
- fMyCumuAnalysis->Make(fEvent);
+ fCA->Make(fEvent);
delete fEvent;
}
PostData(0,fListHistos);
}
-//________________________________________________________________________
+//================================================================================================================
+
void AliAnalysisTaskCumulants::Terminate(Option_t *)
{
- //=====================================================================================================
- // Accessing the output file which contains the merged results from all workers;
- fListHistos = (TList*)GetOutputData(0);
- //fListHistos->Print();//printing the list of stored histos
- //=====================================================================================================
-
- if(fListHistos)
- {
- // Profiles with values of generating functions
- TProfile2D *fIntFlowGenFun=dynamic_cast<TProfile2D*>(fListHistos->FindObject("fIntFlowGenFun"));
- TProfile3D *fDiffFlowGenFunRe=dynamic_cast<TProfile3D*>(fListHistos->FindObject("fDiffFlowGenFunRe"));
- TProfile3D *fDiffFlowGenFunIm=dynamic_cast<TProfile3D*>(fListHistos->FindObject("fDiffFlowGenFunIm"));
+ //accessing the output list which contains the merged 2D and 3D profiles from all worker nodes
+ fListHistos = (TList*)GetOutputData(0);
+ //fListHistos->Print();
+
+ if(fListHistos)
+ {
+ //profiles with avarage values of generating functions for int. and diff. flow
+ TProfile2D *intFlowGenFun = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fIntFlowGenFun"));
+ //TProfile3D *diffFlowGenFunRe=dynamic_cast<TProfile3D*>(fListHistos->FindObject("fDiffFlowGenFunRe"));
+ //TProfile3D *diffFlowGenFunIm=dynamic_cast<TProfile3D*>(fListHistos->FindObject("fDiffFlowGenFunIm"));
- // Histograms to store final results
- TH1D *fIntFlowResults=dynamic_cast<TH1D*>(fListHistos->FindObject("fIntFlowResults"));
- TH1D *fDiffFlowResults2=dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults2"));
- TH1D *fDiffFlowResults4=dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults4"));
- TH1D *fDiffFlowResults6=dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults6"));
- TH1D *fDiffFlowResults8=dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults8"));
-
- // Avarage multiplicity
- TProfile *AvMult=dynamic_cast<TProfile*>(fListHistos->FindObject("fHistProAvM"));
- Double_t BvM=AvMult->GetBinContent(1);//avarage multiplicity
-
- // Calling the function which from values of generating functions calculate integrated and differential flow and store the final results into histograms
- AliCumulantsFunctions FinalResults(fIntFlowGenFun,fDiffFlowGenFunRe,fDiffFlowGenFunIm,fIntFlowResults,fDiffFlowResults2,fDiffFlowResults4,fDiffFlowResults6,fDiffFlowResults8,BvM);
- FinalResults.Calculate();
- }
- else
- {
- cout<<"histogram list pointer is empty"<<endl;
- }
+ TProfile2D *diffFlowGenFunRe0 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe0"));
+ TProfile2D *diffFlowGenFunRe1 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe1"));
+ TProfile2D *diffFlowGenFunRe2 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe2"));
+ TProfile2D *diffFlowGenFunRe3 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe3"));
+ TProfile2D *diffFlowGenFunRe4 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe4"));
+ TProfile2D *diffFlowGenFunRe5 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe5"));
+ TProfile2D *diffFlowGenFunRe6 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe6"));
+ TProfile2D *diffFlowGenFunRe7 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunRe7"));
+ TProfile2D *diffFlowGenFunIm0 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm0"));
+ TProfile2D *diffFlowGenFunIm1 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm1"));
+ TProfile2D *diffFlowGenFunIm2 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm2"));
+ TProfile2D *diffFlowGenFunIm3 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm3"));
+ TProfile2D *diffFlowGenFunIm4 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm4"));
+ TProfile2D *diffFlowGenFunIm5 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm5"));
+ TProfile2D *diffFlowGenFunIm6 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm6"));
+ TProfile2D *diffFlowGenFunIm7 = dynamic_cast<TProfile2D*>(fListHistos->FindObject("fDiffFlowGenFunIm7"));
+
+ //histograms to store the final results
+ TH1D *intFlowResults = dynamic_cast<TH1D*>(fListHistos->FindObject("fIntFlowResults"));
+ TH1D *diffFlowResults2 = dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults2"));
+ TH1D *diffFlowResults4 = dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults4"));
+ TH1D *diffFlowResults6 = dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults6"));
+ TH1D *diffFlowResults8 = dynamic_cast<TH1D*>(fListHistos->FindObject("fDiffFlowResults8"));
+
+ //profile with avarage selected multiplicity for int. flow
+ TProfile *avMult = dynamic_cast<TProfile*>(fListHistos->FindObject("fAvMultIntFlow"));
+
+ //profile with avarage values of Q-vector components (1st bin: <Q_x>, 2nd bin: <Q_y>, 3rd bin: <(Q_x)^2>, 4th bin: <(Q_y)^2>)
+ TProfile *QVectorComponents = dynamic_cast<TProfile*>(fListHistos->FindObject("fQVectorComponents"));
+
+ //q-distribution
+ TH1D *qDist = dynamic_cast<TH1D*>(fListHistos->FindObject("fQDist"));
+
+ AliCumulantsFunctions finalResults(intFlowGenFun,NULL,NULL, intFlowResults,diffFlowResults2,diffFlowResults4,diffFlowResults6,diffFlowResults8,avMult,QVectorComponents,qDist,diffFlowGenFunRe0,diffFlowGenFunRe1,diffFlowGenFunRe2,
+diffFlowGenFunRe3,diffFlowGenFunRe4,diffFlowGenFunRe5,diffFlowGenFunRe6,diffFlowGenFunRe7,diffFlowGenFunIm0,diffFlowGenFunIm1,
+diffFlowGenFunIm2,diffFlowGenFunIm3,diffFlowGenFunIm4,diffFlowGenFunIm5,diffFlowGenFunIm6,diffFlowGenFunIm7);
+
+ finalResults.Calculate();
+
+ }
+ else
+ {
+ cout<<"histogram list pointer is empty"<<endl;
+ }
}
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
-* See cxx source for full Copyright notice */
-/* $Id$ */
-
+/*
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved.
+ * See cxx source for full Copyright notice
+ * $Id$
+ */
+
+/**************************************
+ * analysis task for cumulant method *
+ * *
+ * authors: Naomi van der Kolk *
+ * (kolk@nikhef.nl) *
+ * Raimond Snellings *
+ * (snelling@nikhef.nl) *
+ * Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ * ***********************************/
#ifndef AliAnalysisTaskCumulants_H
#define AliAnalysisTaskCumulants_H
-// AliAnalysisTaskCumulants:
-// analysis task for
-// Cumulant method
-// with many authors (N.K. R.S. A.B.)
-// who do something
+#include "AliAnalysisTask.h"
class AliESDEvent;
class AliAODEvent;
class AliFlowEventSimpleMaker;
class TFile;
-#include "AliAnalysisTask.h"
+//================================================================================================================
-class AliAnalysisTaskCumulants : public AliAnalysisTask {
+class AliAnalysisTaskCumulants : public AliAnalysisTask{
public:
- //AliAnalysisTaskCumulants(const char *name = "AliAnalysisTaskCumulants");
AliAnalysisTaskCumulants();
AliAnalysisTaskCumulants(const char *name);
- virtual ~AliAnalysisTaskCumulants() {}
+ virtual ~AliAnalysisTaskCumulants(){};
virtual void ConnectInputData(Option_t *);
virtual void CreateOutputObjects();
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 SetCFManager1(AliCFManager* cfmgr) {this->fCFManager1 = cfmgr; }
- AliCFManager* GetCFManager1() {return this->fCFManager1; }
- void SetCFManager2(AliCFManager* cfmgr) {this->fCFManager2 = cfmgr; }
- AliCFManager* GetCFManager2() {return this->fCFManager2; }
+ 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;}
private:
-
- AliAnalysisTaskCumulants(const AliAnalysisTaskCumulants& aAnalysis);
- AliAnalysisTaskCumulants& operator=(const AliAnalysisTaskCumulants& aAnalysis);
-
- AliESDEvent *fESD; //ESD object
- AliAODEvent* fAOD; //AOD object
- AliFlowAnalysisWithCumulants* fMyCumuAnalysis; //Cumulant analysis object
- AliFlowEventSimpleMaker* fEventMaker; //FlowEventSimple maker 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
- TList *fListHistos; //collection of output //NEW
-
- ClassDef(AliAnalysisTaskCumulants, 1); // example of analysis
+ AliAnalysisTaskCumulants(const AliAnalysisTaskCumulants& aatc);
+ AliAnalysisTaskCumulants& operator=(const AliAnalysisTaskCumulants& aatc);
+
+ AliESDEvent *fESD; //ESD object
+ AliAODEvent* fAOD; //AOD object
+ AliFlowAnalysisWithCumulants* fCA; //Cumulant Analysis (CA) object
+ AliFlowEventSimpleMaker* fEventMaker; //FlowEventSimple maker 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
+ TList *fListHistos; //collection of output
+
+ ClassDef(AliAnalysisTaskCumulants, 1);
};
+//================================================================================================================
+
#endif
-//*******************************
-// flow analysis with cumulants *
-// author: Ante Bilandzic *
-// email: anteb@nikhef.nl *
-//*******************************
+/*************************************************************************
+* Copyright(c) 1998-2008, 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. *
+**************************************************************************/
+
+/**********************************
+ * functions and equations needed *
+ * for calculation of cumulants *
+ * and final flow estimates *
+ * *
+ * author: Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ *********************************/
#define AliCumulantsFunctions_cxx
#include "AliFlowCommonConstants.h"
#include "TChain.h"
#include "TFile.h"
-#include "TList.h" //NEW
+#include "TList.h"
#include "TParticle.h"
#include "TProfile.h"
#include "TProfile2D.h"
#include "TProfile3D.h"
+#include "TF1.h"//VOLOSHIN
+#include "TAxis.h"//VOLOSHIN
+#include "TH1.h"
+#include "TH1D.h"
+
#include "AliFlowEventSimple.h"
#include "AliFlowTrackSimple.h"
#include "AliFlowCommonConstants.h"
#include "AliCumulantsFunctions.h"
-
ClassImp(AliCumulantsFunctions)
-//________________________________________________________________________
+//================================================================================================================_
AliCumulantsFunctions::AliCumulantsFunctions():
- fIG(NULL),
- fDGRe(NULL),
- fDGIm(NULL),
- fifr(NULL),
- fdfr2(NULL),
- fdfr4(NULL),
- fdfr6(NULL),
- fdfr8(NULL),
- fAvMult(0)
- {
- //default constructor
- }
+ fIntGenFun(NULL),
+ fDiffGenFunRe(NULL),
+ fDiffGenFunIm(NULL),
+ fifr(NULL),
+ fdfr2(NULL),
+ fdfr4(NULL),
+ fdfr6(NULL),
+ fdfr8(NULL),
+ fAvMult(NULL),
+ fQVector(NULL),
+ fQDistrib(NULL),//q-distribution
+ fdRe0(NULL),
+ fdRe1(NULL),
+ fdRe2(NULL),
+ fdRe3(NULL),
+ fdRe4(NULL),
+ fdRe5(NULL),
+ fdRe6(NULL),
+ fdRe7(NULL),
+ fdIm0(NULL),
+ fdIm1(NULL),
+ fdIm2(NULL),
+ fdIm3(NULL),
+ fdIm4(NULL),
+ fdIm5(NULL),
+ fdIm6(NULL),
+ fdIm7(NULL)
+{
+ //default constructor
+}
-AliCumulantsFunctions::~AliCumulantsFunctions(){
- //desctructor
+AliCumulantsFunctions::~AliCumulantsFunctions()
+{
+ //destructor
}
-AliCumulantsFunctions::AliCumulantsFunctions(TProfile2D *IntGen, TProfile3D *DiffGenRe, TProfile3D *DiffGenIm, TH1D *ifr, TH1D *dfr2, TH1D *dfr4, TH1D *dfr6, TH1D *dfr8, Double_t CvM):
- fIG(IntGen),
- fDGRe(DiffGenRe),
- fDGIm(DiffGenIm),
- fifr(ifr),
- fdfr2(dfr2),
- fdfr4(dfr4),
- fdfr6(dfr6),
- fdfr8(dfr8),
- fAvMult(CvM)
- {
- //custom constructor
- }
+AliCumulantsFunctions::AliCumulantsFunctions(TProfile2D *IntGenFun, TProfile3D *DiffGenFunRe, TProfile3D *DiffGenFunIm, TH1D *ifr, TH1D *dfr2, TH1D *dfr4, TH1D *dfr6, TH1D *dfr8, TProfile *AvMult, TProfile *QVector, TH1D *QDistrib, TProfile2D *dRe0, TProfile2D *dRe1, TProfile2D *dRe2, TProfile2D *dRe3, TProfile2D *dRe4, TProfile2D *dRe5, TProfile2D *dRe6, TProfile2D *dRe7, TProfile2D *dIm0, TProfile2D *dIm1, TProfile2D *dIm2, TProfile2D *dIm3, TProfile2D *dIm4, TProfile2D *dIm5, TProfile2D *dIm6, TProfile2D *dIm7):
+ fIntGenFun(IntGenFun),
+ fDiffGenFunRe(DiffGenFunRe),
+ fDiffGenFunIm(DiffGenFunIm),
+ fifr(ifr),
+ fdfr2(dfr2),
+ fdfr4(dfr4),
+ fdfr6(dfr6),
+ fdfr8(dfr8),
+ fAvMult(AvMult),
+ fQVector(QVector),
+ fQDistrib(QDistrib),//q-distribution
+ fdRe0(dRe0),
+ fdRe1(dRe1),
+ fdRe2(dRe2),
+ fdRe3(dRe3),
+ fdRe4(dRe4),
+ fdRe5(dRe5),
+ fdRe6(dRe6),
+ fdRe7(dRe7),
+ fdIm0(dIm0),
+ fdIm1(dIm1),
+ fdIm2(dIm2),
+ fdIm3(dIm3),
+ fdIm4(dIm4),
+ fdIm5(dIm5),
+ fdIm6(dIm6),
+ fdIm7(dIm7)
+{
+ //custom constructor
+}
-//___________________________________________________________________________
-void AliCumulantsFunctions::Calculate(){
- //calculate final flow estimates
-
- static const Int_t fgkQmax=AliFlowCumuConstants::kQmax;//needed for numerics
- static const Int_t fgkPmax=AliFlowCumuConstants::kPmax;//needed for numerics
- static const Int_t fgkFlow=AliFlowCumuConstants::kFlow;//integrated flow coefficient to be calculated
- static const Int_t fgkMltpl=AliFlowCumuConstants::kMltpl;//the multiple in p=m*n (diff. flow)
- static const Int_t fgknBins=100;//number of pt bins //to be improved
- Double_t fR0=AliFlowCumuConstants::fgR0;//needed for numerics
- Double_t fPtMax=AliFlowCommonConstants::GetPtMax();//maximum pt
- Double_t fPtMin=AliFlowCommonConstants::GetPtMin();//minimum pt
- Double_t fBinWidth=(fPtMax-fPtMin)/fgknBins;//width of pt bin (in GeV)
+//================================================================================================================
+
+void AliCumulantsFunctions::Calculate()
+{
+ //calculate cumulants and final integrated and differential flow estimates and store the results into output histograms
+ static const Int_t fgkQmax=AliFlowCumuConstants::kQmax; //needed for numerics
+ static const Int_t fgkPmax=AliFlowCumuConstants::kPmax; //needed for numerics
+ static const Int_t fgkFlow=AliFlowCumuConstants::kFlow; //integrated flow coefficient to be calculated
+ static const Int_t fgkMltpl=AliFlowCumuConstants::kMltpl; //the multiple in p=m*n (diff. flow)
+ static const Int_t fgknBins=100; //number of pt bins //to be improved
+ Double_t fR0=AliFlowCumuConstants::fgR0; //needed for numerics
+ //Double_t fPtMax=AliFlowCommonConstants::GetPtMax(); //maximum pt
+ //Double_t fPtMin=AliFlowCommonConstants::GetPtMin(); //minimum pt
+ //Double_t fBinWidth=(fPtMax-fPtMin)/fgknBins; //width of pt bin (in GeV)
- Double_t fBvG[fgkPmax][fgkQmax]={0.};
-
- for(Int_t p=0;p<fgkPmax;p++){
- for(Int_t q=0;q<fgkQmax;q++){
- fBvG[p][q]=fIG->GetBinContent(p+1,q+1);
- }
- }
-
-
- /////////////////////////////////////////////////////////////////////////////
- //////////////////gen. function for the cumulants////////////////////////////
- /////////////////////////////////////////////////////////////////////////////
-
- Double_t fC[fgkPmax][fgkQmax]={0.};
- for (Int_t p=0;p<fgkPmax;p++){
- for (Int_t q=0;q<fgkQmax;q++){
- fC[p][q]=1.*fAvMult*(pow(fBvG[p][q],(1./fAvMult))-1.); //to be improved
- }
+ //<G[p][q]>
+ Double_t AvG[fgkPmax][fgkQmax]={0.};
+ for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ AvG[p][q]=fIntGenFun->GetBinContent(p+1,q+1);
}
+ }
+
+ //avarage selected multiplicity
+ Double_t AvM = fAvMult->GetBinContent(1);
+
+ //mumber of events
+ Int_t nEvents = (Int_t)(fAvMult->GetBinEntries(1));
- /////////////////////////////////////////////////////////////////////////////
- ///////avaraging the gen. function for the cumulants over azimuth////////////
- /////////////////////////////////////////////////////////////////////////////
+ //<Q-vector stuff>
+ Double_t AvQx = fQVector->GetBinContent(1); //<Q_x>
+ Double_t AvQy = fQVector->GetBinContent(2); //<Q_y>
+ Double_t AvQ2x = fQVector->GetBinContent(3); //<(Q_x)^2>
+ Double_t AvQ2y = fQVector->GetBinContent(4); //<(Q_y)^2>
- Double_t fCAv[fgkPmax]={0.};
- for (Int_t p=0;p<fgkPmax;p++){
- Double_t fTempHere=0.;
- for (Int_t q=0;q<fgkQmax;q++){
- fTempHere+=1.*fC[p][q];
- }
- fCAv[p]=1.*fTempHere/fgkQmax;
- }
+ /////////////////////////////////////////////////////////////////////////////
+ //////////////////gen. function for the cumulants////////////////////////////
+ /////////////////////////////////////////////////////////////////////////////
- /////////////////////////////////////////////////////////////////////////////
- //////////////////////////////////final results//////////////////////////////
- /////////////////////////////////////////////////////////////////////////////
+ Double_t C[fgkPmax][fgkQmax]={0.};//C[p][q]
+ for (Int_t p=0;p<fgkPmax;p++){
+ for (Int_t q=0;q<fgkQmax;q++){
+ C[p][q]=1.*AvM*(pow(AvG[p][q],(1./AvM))-1.);
+ }
+ }
+
+ /////////////////////////////////////////////////////////////////////////////
+ ///////avaraging the gen. function for the cumulants over azimuth////////////
+ /////////////////////////////////////////////////////////////////////////////
- Double_t fCumulant[fgkPmax];//c_{iFlow}\{2(p+1)\}
+ Double_t AvC[fgkPmax]={0.};//<C[p][q]>
+ for (Int_t p=0;p<fgkPmax;p++){
+ Double_t tempHere=0.;
+ for (Int_t q=0;q<fgkQmax;q++){
+ tempHere+=1.*C[p][q];
+ }
+ AvC[p]=1.*tempHere/fgkQmax;
+ }
+
+ /////////////////////////////////////////////////////////////////////////////
+ //////////////////////////////////final results//////////////////////////////
+ /////////////////////////////////////////////////////////////////////////////
+
+ Double_t cumulant[fgkPmax];//array to store various order cumulants
- fCumulant[0]=(1./(fR0*fR0)) * (8.*fCAv[0] - 14.*fCAv[1] + (56./3.)*fCAv[2] - (35./2.)*fCAv[3] +
- (56./5.)*fCAv[4] - (14./3.)*fCAv[5] + (8./7.)*fCAv[6] - (1./8.)*fCAv[7]);
+ cumulant[0] = (1./(fR0*fR0)) * (8.*AvC[0] - 14.*AvC[1] + (56./3.)*AvC[2] - (35./2.)*AvC[3] +
+ (56./5.)*AvC[4] - (14./3.)*AvC[5] + (8./7.)*AvC[6] - (1./8.)*AvC[7]);
- fCumulant[1]=(1./pow(fR0,4.)) * ((-1924./35.)*fCAv[0] + (621./5.)*fCAv[1] - (8012./45.)*fCAv[2] +
- (691./4.)*fCAv[3] - (564./5.)*fCAv[4] + (2143./45.)*fCAv[5] -
- (412./35.)*fCAv[6] + (363./280.)*fCAv[7]);
+ cumulant[1] = (1./pow(fR0,4.)) * ((-1924./35.)*AvC[0] + (621./5.)*AvC[1] - (8012./45.)*AvC[2] +
+ (691./4.)*AvC[3] - (564./5.)*AvC[4] + (2143./45.)*AvC[5] -
+ (412./35.)*AvC[6] + (363./280.)*AvC[7]);
- fCumulant[2]=(1./pow(fR0,6.)) * (349.*fCAv[0] - (18353./20.)*fCAv[1] + (7173./5.)*fCAv[2] -
- 1457.*fCAv[3] + (4891./5.)*fCAv[4] - (1683./4.)*fCAv[5] +
- (527./5.)*fCAv[6] - (469./40.)*fCAv[7]);
+ cumulant[2] = (1./pow(fR0,6.)) * (349.*AvC[0] - (18353./20.)*AvC[1] + (7173./5.)*AvC[2] -
+ 1457.*AvC[3] + (4891./5.)*AvC[4] - (1683./4.)*AvC[5] +
+ (527./5.)*AvC[6] - (469./40.)*AvC[7]);
- fCumulant[3]=(1./pow(fR0,8.)) * ((-10528./5.)*fCAv[0] + (30578./5.)*fCAv[1] - (51456./5.)*fCAv[2] +
- 10993.*fCAv[3] - (38176./5.)*fCAv[4] + (16818./5.)*fCAv[5] -
- (4288./5.)*fCAv[6] + (967./10.)*fCAv[7]);
+ cumulant[3] = (1./pow(fR0,8.)) * ((-10528./5.)*AvC[0] + (30578./5.)*AvC[1] - (51456./5.)*AvC[2] +
+ 10993.*AvC[3] - (38176./5.)*AvC[4] + (16818./5.)*AvC[5] -
+ (4288./5.)*AvC[6] + (967./10.)*AvC[7]);
- fCumulant[4]=(1./pow(fR0,10.)) * (11500.*fCAv[0] - 35800.*fCAv[1] + 63900.*fCAv[2] - 71600.*fCAv[3] +
- 51620.*fCAv[4] - 23400.*fCAv[5] + 6100.*fCAv[6] - 700.*fCAv[7]);
+ cumulant[4] = (1./pow(fR0,10.)) * (11500.*AvC[0] - 35800.*AvC[1] + 63900.*AvC[2] - 71600.*AvC[3] +
+ 51620.*AvC[4] - 23400.*AvC[5] + 6100.*AvC[6] - 700.*AvC[7]);
- fCumulant[5]=(1./pow(fR0,12.)) * (-52560.*fCAv[0] + 172080.*fCAv[1] - 321840.*fCAv[2] + 376200.*fCAv[3] -
- 281520.*fCAv[4] + 131760.*fCAv[5] - 35280.*fCAv[6] + 4140.*fCAv[7]);
+ cumulant[5] = (1./pow(fR0,12.)) * (-52560.*AvC[0] + 172080.*AvC[1] - 321840.*AvC[2] + 376200.*AvC[3] -
+ 281520.*AvC[4] + 131760.*AvC[5] - 35280.*AvC[6] + 4140.*AvC[7]);
- fCumulant[6]=(1./pow(fR0,14.)) * (176400.*fCAv[0] - 599760.*fCAv[1] + 1164240.*fCAv[2] - 1411200.*fCAv[3] +
- 1093680.*fCAv[4] - 529200.*fCAv[5] + 146160.*fCAv[6] - 17640.*fCAv[7]);
+ cumulant[6] = (1./pow(fR0,14.)) * (176400.*AvC[0] - 599760.*AvC[1] + 1164240.*AvC[2] - 1411200.*AvC[3] +
+ 1093680.*AvC[4] - 529200.*AvC[5] + 146160.*AvC[6] - 17640.*AvC[7]);
- fCumulant[7]=(1./pow(fR0,16.)) * (-322560*fCAv[0] + 1128960.*fCAv[1] - 2257920.*fCAv[2] + 2822400.*fCAv[3] -
- 2257920.*fCAv[4] + 1128960.*fCAv[5] - 322560.*fCAv[6] + 40320.*fCAv[7]);
-
+ cumulant[7] = (1./pow(fR0,16.)) * (-322560*AvC[0] + 1128960.*AvC[1] - 2257920.*AvC[2] + 2822400.*AvC[3] -
+ 2257920.*AvC[4] + 1128960.*AvC[5] - 322560.*AvC[6] + 40320.*AvC[7]);
+
+ cout<<""<<endl;
+ cout<<"***************************"<<endl;
+ cout<<"cumulants:"<<endl;
- cout<<""<<endl;
- cout<<"***************************"<<endl;
- cout<<"cumulants:"<<endl;
+ cout<<" c_"<<fgkFlow<<"{2} = "<<cumulant[0]<<endl;
+ cout<<" c_"<<fgkFlow<<"{4} = "<<cumulant[1]<<endl;
+ cout<<" c_"<<fgkFlow<<"{6} = "<<cumulant[2]<<endl;
+ cout<<" c_"<<fgkFlow<<"{8} = "<<cumulant[3]<<endl;
+ cout<<"c_"<<fgkFlow<<"{10} = "<<cumulant[4]<<endl;
+ cout<<"c_"<<fgkFlow<<"{12} = "<<cumulant[5]<<endl;
+ cout<<"c_"<<fgkFlow<<"{14} = "<<cumulant[6]<<endl;
+ cout<<"c_"<<fgkFlow<<"{16} = "<<cumulant[7]<<endl;
- cout<<" c_"<<fgkFlow<<"{2} = "<<fCumulant[0]<<endl;
- cout<<" c_"<<fgkFlow<<"{4} = "<<fCumulant[1]<<endl;
- cout<<" c_"<<fgkFlow<<"{6} = "<<fCumulant[2]<<endl;
- cout<<" c_"<<fgkFlow<<"{8} = "<<fCumulant[3]<<endl;
- cout<<"c_"<<fgkFlow<<"{10} = "<<fCumulant[4]<<endl;
- cout<<"c_"<<fgkFlow<<"{12} = "<<fCumulant[5]<<endl;
- cout<<"c_"<<fgkFlow<<"{14} = "<<fCumulant[6]<<endl;
- cout<<"c_"<<fgkFlow<<"{16} = "<<fCumulant[7]<<endl;
-
- cout<<""<<endl;
- cout<<"integrated flow: "<<endl;
+ cout<<""<<endl;
+ cout<<"integrated flow: "<<endl;
+
+ Double_t V2=0.,V4=0.,V6=0.,V8=0.,V10=0.,V12=0.,V14=0.,V16=0.;
+
+ if(cumulant[0]>=0.){
+ V2 = pow(cumulant[0],(1./2.));
+ }
+ if(cumulant[1]<=0.){
+ V4 = pow(-cumulant[1],(1./4.));
+ }
+ if(cumulant[2]>=0.){
+ V6 = pow((1./4.)*cumulant[2],(1./6.));
+ }
+ if(cumulant[3]<=0.){
+ V8 = pow(-(1./33.)*cumulant[3],(1./8.));
+ }
+ if(cumulant[4]>=0.){
+ V10 = pow((1./456.)*cumulant[4],(1./10.));
+ }
+ if(cumulant[5]<=0.){
+ V12 = pow(-(1./9460.)*cumulant[5],(1./12.));
+ }
+ if(cumulant[6]>=0.){
+ V14 = pow((1./274800.)*cumulant[6],(1./14.));
+ }
+ if(cumulant[7]<=0.){
+ V16 = pow(-(1./10643745.)*cumulant[7],(1./16.));
+ }
+
+ Double_t SdQ[4]={0.};
+ Double_t ChiQ[4]={0.};
+
+ //v_2{2}
+ if((cumulant[0]>=0.) && (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(pow(cumulant[0],(1./2.))*AvM,2.)>0.))
+ {
+ ChiQ[0]=AvM*V2/pow(AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V2*AvM,2.),0.5);
+ SdQ[0]=pow(((1./(2.*AvM*nEvents))*((1.+1.*pow(ChiQ[0],2))/(1.*pow(ChiQ[0],2)))),0.5);
+ cout<<" v_"<<fgkFlow<<"{2} = "<<100.*V2<<"%, chi{2} = "<<ChiQ[0]<<", sd{2} = "<<100.*SdQ[0]<<"%"<<endl;
+ fifr->SetBinContent(1,100.*V2);
+ fifr->SetBinError(1,100.*SdQ[0]);
+ }
+ else
+ {
+ cout<<" v_"<<fgkFlow<<"{2} = Im"<<endl;
+ }
+
+ //v_2{4}
+ if((cumulant[1]<=0.) && (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(pow(-cumulant[1],(1./4.))*AvM,2.)>0.))
+ {
+ ChiQ[1]=AvM*V4/pow(AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V4*AvM,2.),0.5);
+ SdQ[1]=(1./(pow(2.*AvM*nEvents,.5)))*pow((1.+2.*pow(ChiQ[1],2)+(1./4.)*pow(ChiQ[1],4.)+(1./4.)*pow(ChiQ[1],6.))/((1./4.)*pow(ChiQ[1],6.)),.5);
+ cout<<" v_"<<fgkFlow<<"{4} = "<<100.*V4<<"%, chi{4} = "<<ChiQ[1]<<", sd{4} = "<<100.*SdQ[1]<<"%"<<endl;
+ fifr->SetBinContent(2,100.*V4);
+ fifr->SetBinError(2,100.*SdQ[1]);
+ }
+ else
+ {
+ cout<<" v_"<<fgkFlow<<"{4} = Im"<<endl;
+ }
+ //v_2{6}
+ if((cumulant[2]>=0.) && (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(pow((1./4.)*cumulant[2],(1./6.))*AvM,2.)>0.))
+ {
+ ChiQ[2]=AvM*V6/pow(AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V6*AvM,2.),0.5);
+ SdQ[2]=(1./(pow(2.*AvM*nEvents,.5)))*pow((3.+18.*pow(ChiQ[2],2)+9.*pow(ChiQ[2],4.)+28.*pow(ChiQ[2],6.)+12.*pow(ChiQ[2],8.)+24.*pow(ChiQ[2],10.))/(24.*pow(ChiQ[2],10.)),.5);
+ cout<<" v_"<<fgkFlow<<"{6} = "<<100.*V6<<"%, chi{6} = "<<ChiQ[2]<<", sd{6} = "<<100.*SdQ[2]<<"%"<<endl;
+ fifr->SetBinContent(3,100.*V6);
+ fifr->SetBinError(3,100.*SdQ[2]);
+ }
+ else
+ {
+ cout<<" v_"<<fgkFlow<<"{6} = Im"<<endl;
+ }
- Double_t fV2=0.,fV4=0.,fV6=0.,fV8=0.;
- Double_t fSdQ[4]={0.};
- Double_t fChiQ[4]={0.};
- if (fCumulant[0]>=0.){
- fV2=sqrt(fCumulant[0]);
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV2*fAvM,2.)>0.){
- //fChiQ[0]=fAvM*fV2/pow(fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV2*fAvM,2.),0.5);
- //fSdQ[0]=pow(((1./(2.*fAvM*nEvents))*((1.+1.*pow(fChiQ[0],2))/(1.*pow(fChiQ[0],2)))),0.5);
- cout<<" v_"<<fgkFlow<<"{2} = "<<100.*fV2<<"%, chi{2} = "<<fChiQ[0]<<", sd{2} = "<<100.*fSdQ[0]<<"%"<<endl;//to be improved (2->fgkFlow)
- //fCommonHistsRes2->FillIntegratedFlow(100.*fV2,100.*fSdQ[0]);
- //fCommonHistsRes2->FillChi(fChiQ[0]);
- fifr->SetBinContent(1,100.*fV2);
- //}
- //} else {
- //cout<<" v_"<<fgkFlow<<"{2} = Im"<<endl;
- } else {
- //cout<<" v_"<<fgkFlow<<"{8} = Im"<<endl;
- }
- if (fCumulant[1]<=0.){
- fV4=pow(-fCumulant[1],(1./4.));
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV4*fAvM,2.)>0.){
- //fChiQ[1]=fAvM*fV4/pow(fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV4*fAvM,2.),0.5);
- //fSdQ[1]=(1./(pow(2.*fAvM*nEvents,.5)))*pow((1.+2.*pow(fChiQ[1],2)+(1./4.)*pow(fChiQ[1],4.)+(1./4.)*pow(fChiQ[1],6.))/((1./4.)*pow(fChiQ[1],6.)),.5);
- cout<<" v_"<<fgkFlow<<"{4} = "<<100.*fV4<<"%, chi{4} = "<<fChiQ[1]<<", sd{4} = "<<100.*fSdQ[1]<<"%"<<endl;//to be improved (2->fgkFlow)
- //fCommonHistsRes4->FillInteg8ratedFlow(100.*fV4,100.*fSdQ[1]);
- //fCommonHistsRes4->FillChi(fChiQ[1]);
- fifr->SetBinContent(2,100.*fV4);
- //} else {
- //cout<<" v_"<<fgkFlow<<"{4} = Im"<<endl;
- //}
- // } else {
- // cout<<" v_"<<fgkFlow<<"{4} = Im"<<endl;
+ //v_2{8}
+ if((cumulant[3]<=0.) && (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(pow(-(1./33.)*cumulant[3],(1./8.))*AvM,2.)>0.))
+ {
+ ChiQ[3]=AvM*V8/pow(AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V8*AvM,2.),0.5);
+ SdQ[3]=(1./(pow(2.*AvM*nEvents,.5)))*pow((12.+96.*pow(ChiQ[3],2)+72.*pow(ChiQ[3],4.)+304.*pow(ChiQ[3],6.)+257.*pow(ChiQ[3],8.)+804.*pow(ChiQ[3],10.)+363.*pow(ChiQ[3],12.)+726.*pow(ChiQ[3],14.))/(726.*pow(ChiQ[3],14.)),.5);
+ cout<<" v_"<<fgkFlow<<"{8} = "<<100.*V8<<"%, chi{8} = "<<ChiQ[3]<<", sd{8} = "<<100.*SdQ[3]<<"%"<<endl;
+ fifr->SetBinContent(4,100.*V8);
+ fifr->SetBinError(4,100.*SdQ[3]);
}
- if (fCumulant[2]>=0.){
- fV6=pow((1./4.)*fCumulant[2],(1./6.));
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV6*fAvM,2.)>0.){
- //fChiQ[2]=fAvM*fV6/pow(fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV6*fAvM,2.),0.5);
- //fSdQ[2]=(1./(pow(2.*fAvM*nEvents,.5)))*pow((3.+18.*pow(fChiQ[2],2)+9.*pow(fChiQ[2],4.)+28.*pow(fChiQ[2],6.)+12.*pow(fChiQ[2],8.)+24.*pow(fChiQ[2],10.))/(24.*pow(fChiQ[2],10.)),.5);
- cout<<" v_"<<fgkFlow<<"{6} = "<<100.*fV6<<"%, chi{6} = "<<fChiQ[2]<<", sd{6} = "<<100.*fSdQ[2]<<"%"<<endl;//to be improved (2->fgkFlow)
- //fCommonHistsRes6->FillIntegratedFlow(100.*fV6,100.*fSdQ[2]);
- //fCommonHistsRes6->FillChi(fChiQ[2]);
- fifr->SetBinContent(3,100.*fV6);
- //} else {
- // cout<<" v_"<<fgkFlow<<"{6} = Im"<<endl;
- // }
- } else {
- //cout<<" v_"<<fgkFlow<<"{6} = Im"<<endl;
+ else
+ {
+ cout<<" v_"<<fgkFlow<<"{8} = Im"<<endl;
}
- if (fCumulant[3]<=0.){
- fV8=pow(-(1./33.)*fCumulant[3],(1./8.));
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV8*fAvM,2.)>0.){
- //fChiQ[3]=fAvM*fV8/pow(fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV8*fAvM,2.),0.5);
- //fSdQ[3]=(1./(pow(2.*fAvM*nEvents,.5)))*pow((12.+96.*pow(fChiQ[3],2)+72.*pow(fChiQ[3],4.)+304.*pow(fChiQ[3],6.)+257.*pow(fChiQ[3],8.)+804.*pow(fChiQ[3],10.)+363.*pow(fChiQ[3],12.)+726.*pow(fChiQ[3],14.))/(726.*pow(fChiQ[3],14.)),.5);
- cout<<" v_"<<fgkFlow<<"{8} = "<<100.*fV8<<"%, chi{8} = "<<fChiQ[3]<<", sd{8} = "<<100.*fSdQ[3]<<"%"<<endl;//to be improved (2->fgkFlow)
- //fCommonHistsRes8->FillIntegratedFlow(100.*fV8,100.*fSdQ[3]);
- //fCommonHistsRes8->FillChi(fChiQ[3]);
- fifr->SetBinContent(4,100.*fV8);
- //} else {
- //cout<<" v_"<<fgkFlow<<"{8} = Im"<<endl;a
- //}
- } else {
- //cout<<" v_"<<fgkFlow<<"{8} = Im"<<endl;
- }
- if (fCumulant[4]>=0.){//fHistProIntFlow
- cout<<"v_"<<fgkFlow<<"{10} = "<<100.*pow((1./456.)*fCumulant[4],(1./10.))<<"%"<<endl;//to be improved (2->fgkFlow)
- fifr->SetBinContent(5,100.*pow((1./456.)*fCumulant[4],(1./10.)));
+
+ //v_2{10}
+ if (cumulant[4]>=0.){
+ cout<<"v_"<<fgkFlow<<"{10} = "<<100.*V10<<"%"<<endl;
+ fifr->SetBinContent(5,100.*pow((1./456.)*cumulant[4],(1./10.)));
} else {
- cout<<"v_"<<fgkFlow<<"{10} = Im"<<endl; //to be improved (2->fgkFlow)
+ cout<<"v_"<<fgkFlow<<"{10} = Im"<<endl;
}
- if (fCumulant[5]<=0.){
- cout<<"v_"<<fgkFlow<<"{12} = "<<100.*pow(-(1./9460.)*fCumulant[5],(1./12.))<<"%"<<endl;//to be improved (2->fgkFlow)
- fifr->SetBinContent(6,100.*pow(-(1./9460.)*fCumulant[5],(1./12.)));
+
+ //v_2{12}
+ if (cumulant[5]<=0.){
+ cout<<"v_"<<fgkFlow<<"{12} = "<<100.*V12<<"%"<<endl;
+ fifr->SetBinContent(6,100.*pow(-(1./9460.)*cumulant[5],(1./12.)));
} else {
- cout<<"v_"<<fgkFlow<<"{12} = Im"<<endl; //to be improved (2->fgkFlow)
+ cout<<"v_"<<fgkFlow<<"{12} = Im"<<endl;
}
- if (fCumulant[6]>=0.){
- cout<<"v_"<<fgkFlow<<"{14} = "<<100.*pow((1./274800.)*fCumulant[6],(1./14.))<<"%"<<endl;//to be improved (2->fgkFlow)
- fifr->SetBinContent(7,100.*pow((1./274800.)*fCumulant[6],(1./14.)));
+
+ //v_2{14}
+ if (cumulant[6]>=0.){
+ cout<<"v_"<<fgkFlow<<"{14} = "<<100.*V14<<"%"<<endl;
+ fifr->SetBinContent(7,100.*pow((1./274800.)*cumulant[6],(1./14.)));
} else {
- cout<<"v_"<<fgkFlow<<"{14} = Im"<<endl; //to be improved (2->fgkFlow)
+ cout<<"v_"<<fgkFlow<<"{14} = Im"<<endl;
}
- if (fCumulant[7]<=0.){
- cout<<"v_"<<fgkFlow<<"{16} = "<<100.*pow(-(1./10643745.)*fCumulant[7],(1./16.))<<"%"<<endl;//to be improved (2->fgkFlow)
- fifr->SetBinContent(8,100.*pow(-(1./10643745.)*fCumulant[7],(1./16.)));
+
+ //v_2{16}
+ if (cumulant[7]<=0.){
+ cout<<"v_"<<fgkFlow<<"{16} = "<<100.*V16<<"%"<<endl;
+ fifr->SetBinContent(8,100.*pow(-(1./10643745.)*cumulant[7],(1./16.)));
} else {
- cout<<"v_"<<fgkFlow<<"{16} = Im"<<endl; //to be improved (2->fgkFlow)
+ cout<<"v_"<<fgkFlow<<"{16} = Im"<<endl;
}
- //cout<<"***************************"<<endl;
-
- //DIFFERENTIAL FLOW CALCULATIONS STARTS HERE!!!
+ cout<<"***************************"<<endl;
+
- Double_t fX[fgknBins][fgkPmax][fgkQmax]={0.};//see the text bellow relation (11) in PG
- Double_t fY[fgknBins][fgkPmax][fgkQmax]={0.};
+ /////////////////////////////////////////////////////////////////////////////
+ ///////////////////////DIFFERENTIAL FLOW CALCULATIONS////////////////////////
+ /////////////////////////////////////////////////////////////////////////////
+ Double_t X[fgknBins][fgkPmax][fgkQmax]={0.};
+ Double_t Y[fgknBins][fgkPmax][fgkQmax]={0.};
+
+ /*
+ //3D profiles
for(Int_t b=0;b<fgknBins;b++){
for(Int_t p=0;p<fgkPmax;p++){
for(Int_t q=0;q<fgkQmax;q++){
- fX[b][p][q]=fDGRe->GetBinContent(b+1,p+1,q+1)/fBvG[p][q];
- fY[b][p][q]=fDGIm->GetBinContent(b+1,p+1,q+1)/fBvG[p][q];
+ X[b][p][q]=fDiffGenFunRe->GetBinContent(b+1,p+1,q+1)/AvG[p][q];
+ Y[b][p][q]=fDiffGenFunIm->GetBinContent(b+1,p+1,q+1)/AvG[p][q];
}
}
}
+ */
+
+ for(Int_t b=0;b<fgknBins;b++){
+ //for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ X[b][0][q]=fdRe0->GetBinContent(b+1,q+1)/AvG[0][q];
+ Y[b][0][q]=fdIm0->GetBinContent(b+1,q+1)/AvG[0][q];
+ //--------------------------------------------------
+ X[b][1][q]=fdRe1->GetBinContent(b+1,q+1)/AvG[1][q];
+ Y[b][1][q]=fdIm1->GetBinContent(b+1,q+1)/AvG[1][q];
+ //--------------------------------------------------
+ X[b][2][q]=fdRe2->GetBinContent(b+1,q+1)/AvG[2][q];
+ Y[b][2][q]=fdIm2->GetBinContent(b+1,q+1)/AvG[2][q];
+ //--------------------------------------------------
+ X[b][3][q]=fdRe3->GetBinContent(b+1,q+1)/AvG[3][q];
+ Y[b][3][q]=fdIm3->GetBinContent(b+1,q+1)/AvG[3][q];
+ //--------------------------------------------------
+ X[b][4][q]=fdRe4->GetBinContent(b+1,q+1)/AvG[4][q];
+ Y[b][4][q]=fdIm4->GetBinContent(b+1,q+1)/AvG[4][q];
+ //--------------------------------------------------
+ X[b][5][q]=fdRe5->GetBinContent(b+1,q+1)/AvG[5][q];
+ Y[b][5][q]=fdIm5->GetBinContent(b+1,q+1)/AvG[5][q];
+ //--------------------------------------------------
+ X[b][6][q]=fdRe6->GetBinContent(b+1,q+1)/AvG[6][q];
+ Y[b][6][q]=fdIm6->GetBinContent(b+1,q+1)/AvG[6][q];
+ //--------------------------------------------------
+ X[b][7][q]=fdRe7->GetBinContent(b+1,q+1)/AvG[7][q];
+ Y[b][7][q]=fdIm7->GetBinContent(b+1,q+1)/AvG[7][q];
+ }
+ //}
+ }
-
- Double_t fD[fgknBins][fgkPmax]={0.};//implementing relation (11) from PG
+ Double_t D[fgknBins][fgkPmax]={0.};
for (Int_t b=0;b<fgknBins;b++){
for (Int_t p=0;p<fgkPmax;p++){
Double_t fTempHere3=0.;
for (Int_t q=0;q<fgkQmax;q++){
- fTempHere3+=cos(fgkMltpl*2.*q*TMath::Pi()/fgkQmax)*fX[b][p][q] + sin(fgkMltpl*2.*q*TMath::Pi()/fgkQmax)*fY[b][p][q];
+ fTempHere3+=cos(fgkMltpl*2.*q*TMath::Pi()/fgkQmax)*X[b][p][q] + sin(fgkMltpl*2.*q*TMath::Pi()/fgkQmax)*Y[b][p][q];
}
- fD[b][p]=1.*(pow(fR0*pow(p+1.,.5),fgkMltpl)/fgkQmax)*fTempHere3;
- if(fD[b][p]){
- //cout<<"And this "<<b<<" "<<fD[b][p]<<endl;
- }
+ D[b][p]=1.*(pow(fR0*pow(p+1.,.5),fgkMltpl)/fgkQmax)*fTempHere3;
}
}
- Double_t fDiffCumulant2[fgknBins]={0.};//implementing relation (12) from PG
- Double_t fDiffCumulant4[fgknBins]={0.};
- Double_t fDiffCumulant6[fgknBins]={0.};
- Double_t fDiffCumulant8[fgknBins]={0.};
+ Double_t DiffCumulant2[fgknBins]={0.};
+ Double_t DiffCumulant4[fgknBins]={0.};
+ Double_t DiffCumulant6[fgknBins]={0.};
+ Double_t DiffCumulant8[fgknBins]={0.};
for (Int_t b=0;b<fgknBins;b++){
- fDiffCumulant2[b]=(1./(fR0*fR0))*(4.*fD[b][0]-3.*fD[b][1]+(4./3.)*fD[b][2]-(1./4.)*fD[b][3]);
- fDiffCumulant4[b]=(1./pow(fR0,4.))*((-26./3.)*fD[b][0]+(19./2.)*fD[b][1]-(14./3.)*fD[b][2]+(11./12.)*fD[b][3]);
- fDiffCumulant6[b]=(1./pow(fR0,6.))*(18.*fD[b][0]-24.*fD[b][1]+14.*fD[b][2]-3.*fD[b][3]);
- fDiffCumulant8[b]=(1./pow(fR0,8.))*(-24.*fD[b][0]+36.*fD[b][1]-24.*fD[b][2]+6.*fD[b][3]);
+ DiffCumulant2[b]=(1./(fR0*fR0))*(4.*D[b][0]-3.*D[b][1]+(4./3.)*D[b][2]-(1./4.)*D[b][3]);
+ DiffCumulant4[b]=(1./pow(fR0,4.))*((-26./3.)*D[b][0]+(19./2.)*D[b][1]-(14./3.)*D[b][2]+(11./12.)*D[b][3]);
+ DiffCumulant6[b]=(1./pow(fR0,6.))*(18.*D[b][0]-24.*D[b][1]+14.*D[b][2]-3.*D[b][3]);
+ DiffCumulant8[b]=(1./pow(fR0,8.))*(-24.*D[b][0]+36.*D[b][1]-24.*D[b][2]+6.*D[b][3]);
}
- Double_t fv2[fgknBins],fv4[fgknBins],fv6[fgknBins],fv8[fgknBins];
- //Double_t fAvPt[fgknBins];
- Double_t fSddiff2[fgknBins],fSddiff4[fgknBins];
+ Double_t v2[fgknBins],v4[fgknBins],v6[fgknBins],v8[fgknBins];
+ //Double_t Sddiff2[fgknBins],Sddiff4[fgknBins];
- cout<<"number of pt bins: "<<fgknBins<<endl;
- cout<<"****************************************"<<endl;
+ //cout<<"number of pt bins: "<<fgknBins<<endl;
+ //cout<<"****************************************"<<endl;
for (Int_t b=0;b<fgknBins;b++){
- //if(fBinNoOfParticles[b]==0)continue;
- //fAvPt[b]=fBinMeanPt[b]/fBinNoOfParticles[b];
- cout<<"pt bin: "<<b*fBinWidth<<"-"<<(b+1)*fBinWidth<<" GeV"<<endl;
- //cout<<"number of particles in this pt bin: "<<tempNo->GetBinContent(b+1,1,1)<<endl;
- //cout<<"mean pt in this bin: "<<fAvPt[b]<<" GeV"<<endl;
- if(fCumulant[0]>=0){
- fv2[b]=100.*fDiffCumulant2[b]/pow(fCumulant[0],.5);
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV2*fAvM,2.)>0.){
- //fSddiff2[b]=pow((1./(2.*fBinNoOfParticles[b]))*((1.+pow(fChiQ[0],2.))/pow(fChiQ[0],2.)),0.5);
- cout<<"v'_2/2{2} = "<<fv2[b]<<"%, "<<" "<<"sd{2} = "<<100.*fSddiff2[b]<<"%"<<endl;
- //fDiffFlowResults2->SetBinContent(b+1,fv2[b],100.*fSddiff2[b]);
- fdfr2->SetBinContent(b+1,fv2[b]);
- //} else {
+ //cout<<"pt bin: "<<b*fBinWidth<<"-"<<(b+1)*fBinWidth<<" GeV"<<endl;
+
+ //v'_{2/2}{2}
+ if(cumulant[0]>=0)
+ {
+ v2[b]=100.*DiffCumulant2[b]/pow(cumulant[0],.5);
+ if (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V2*AvM,2.)>0.)
+ {
+ //Sddiff2[b]=pow((1./(2.*fBinNoOfParticles[b]))*((1.+pow(ChiQ[0],2.))/pow(ChiQ[0],2.)),0.5);
+ //cout<<"v'_2/2{2} = "<<v2[b]<<"%, "<<" "<<"sd{2} = "<<100.*Sddiff2[b]<<"%"<<endl;
+ //DiffFlowResults2->SetBinContent(b+1,v2[b],100.*Sddiff2[b]);
+ fdfr2->SetBinContent(b+1,v2[b]);
+ } else {
//cout<<"v'_2/2{2} = Im"<<endl;
- //}
+ }
}else{
- cout<<"v'_2/2{2} = Im"<<endl;
+ //cout<<"v'_2/2{2} = Im"<<endl;
}
- if(fCumulant[1]<=0){
- fv4[b]=-100.*fDiffCumulant4[b]/pow(-fCumulant[1],.75);
- //if (fAvQ2x+fAvQ2y-pow(fAvQx,2.)-pow(fAvQy,2.)-pow(fV4*fAvM,2.)>0.){
- //fSddiff4[b]=pow((1./(2.*fBinNoOfParticles[b]))*((2.+6.*pow(fChiQ[1],2.)+pow(fChiQ[1],4.)+pow(fChiQ[1],6.))/pow(fChiQ[1],6.)),0.5);
- cout<<"v'_2/2{4} = "<<fv4[b]<<"%, "<<" "<<"sd{4} = "<<100.*fSddiff4[b]<<"%"<<endl;
- //fCommonHistsRes4->FillDifferentialFlow(b+1,fv4[b],100.*fSddiff4[b]);
- fdfr4->SetBinContent(b+1,fv4[b]);
- //} else {
- // cout<<"v'_2/2{4} = Im"<<endl;
- //}
+
+ //v'_{2/2}{4}
+ if(cumulant[1]<=0)
+ {
+ v4[b]=-100.*DiffCumulant4[b]/pow(-cumulant[1],.75);
+ if (AvQ2x+AvQ2y-pow(AvQx,2.)-pow(AvQy,2.)-pow(V4*AvM,2.)>0.)
+ {
+ //Sddiff4[b]=pow((1./(2.*fBinNoOfParticles[b]))*((2.+6.*pow(ChiQ[1],2.)+pow(ChiQ[1],4.)+pow(ChiQ[1],6.))/pow(ChiQ[1],6.)),0.5);
+ //cout<<"v'_2/2{4} = "<<v4[b]<<"%, "<<" "<<"sd{4} = "<<100.*Sddiff4[b]<<"%"<<endl;
+ //fCommonHistsRes4->FillDifferentialFlow(b+1,v4[b],100.*Sddiff4[b]);
+ fdfr4->SetBinContent(b+1,v4[b]);
+ } else {
+ //cout<<"v'_2/2{4} = Im"<<endl;
+ }
}else{
- cout<<"v'_2/2{4} = Im"<<endl;
+ //cout<<"v'_2/2{4} = Im"<<endl;
}
- if(fCumulant[2]>=0){
- cout<<"v'_2/2{6} = "<<100.*fDiffCumulant6[b]/(4.*pow((1./4.)*fCumulant[2],(5./6.)))<<"%"<<endl;
- fv6[b]=100.*fDiffCumulant6[b]/(4.*pow((1./4.)*fCumulant[2],(5./6.)));
- //fCommonHistsRes6->FillDifferentialFlow(b+1,fv6[b],0.);
- fdfr6->SetBinContent(b+1,fv6[b]);
+
+ //v'_{2/2}{6}
+ if(cumulant[2]>=0){
+ //cout<<"v'_2/2{6} = "<<100.*DiffCumulant6[b]/(4.*pow((1./4.)*cumulant[2],(5./6.)))<<"%"<<endl;
+ v6[b]=100.*DiffCumulant6[b]/(4.*pow((1./4.)*cumulant[2],(5./6.)));
+ //fCommonHistsRes6->FillDifferentialFlow(b+1,v6[b],0.);
+ fdfr6->SetBinContent(b+1,v6[b]);
}else{
- cout<<"v'_2/2{6} = Im"<<endl;
+ //cout<<"v'_2/2{6} = Im"<<endl;
}
- if(fCumulant[3]<=0){
- cout<<"v'_2/2{8} = "<<-100.*fDiffCumulant8[b]/(33.*pow(-(1./33.)*fCumulant[3],(7./8.)))<<"%"<<endl;
- fv8[b]=-100.*fDiffCumulant8[b]/(33.*pow(-(1./33.)*fCumulant[3],(7./8.)));
- //fCommonHistsRes8->FillDifferentialFlow(b+1,fv8[b],0.);
- fdfr8->SetBinContent(b+1,fv8[b]);
+
+ //v'_{2/2}{8}
+ if(cumulant[3]<=0){
+ //cout<<"v'_2/2{8} = "<<-100.*DiffCumulant8[b]/(33.*pow(-(1./33.)*cumulant[3],(7./8.)))<<"%"<<endl;
+ v8[b]=-100.*DiffCumulant8[b]/(33.*pow(-(1./33.)*cumulant[3],(7./8.)));
+ //fCommonHistsRes8->FillDifferentialFlow(b+1,v8[b],0.);
+ fdfr8->SetBinContent(b+1,v8[b]);
}else{
- cout<<"v'_2/2{8} = Im"<<endl;
+ //cout<<"v'_2/2{8} = Im"<<endl;
}
- cout<<"****************************************"<<endl;
+ //cout<<"****************************************"<<endl;
}
+
}
+/*
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved.
+ * See cxx source for full Copyright notice
+ * $Id$
+ */
+
+/**********************************
+ * functions and equations needed *
+ * for calculation of cumulants *
+ * and final flow estimates *
+ * *
+ * author: Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ *********************************/
+
#ifndef AliCumulantsFunctions_H
#define AliCumulantsFunctions_H
-//*********************************
-// functions and equations needed *
-// for calcualation of cumulants *
-// and final flow estimates *
-// author: Ante Bilandzic *
-// email: anteb@nikhef.nl *
-//*********************************
-
#include "AliFlowCommonConstants.h"
#include "AliFlowCumuConstants.h"
class TList;
class TFile;
+//================================================================================================================
+
class AliCumulantsFunctions{
public:
AliCumulantsFunctions();
virtual ~AliCumulantsFunctions();
- AliCumulantsFunctions(TProfile2D *IntGen, TProfile3D *DiffGenRe, TProfile3D *DiffGenIm, TH1D *ifr, TH1D *dfr2, TH1D *dfr4, TH1D *dfr6, TH1D *dfr8, Double_t CvM);
+ AliCumulantsFunctions(TProfile2D *IntGenFun, TProfile3D *DiffGenFunRe, TProfile3D *DiffGenFunIm, TH1D *ifr, TH1D *dfr2, TH1D *dfr4, TH1D *dfr6, TH1D *dfr8, TProfile *AvMult, TProfile *fQVector, TH1D *fQDistrib, TProfile2D *fdRe0, TProfile2D *fdRe1, TProfile2D *fdRe2, TProfile2D *fdRe3, TProfile2D *fdRe4, TProfile2D *fdRe5, TProfile2D *fdRe6, TProfile2D *fdRe7, TProfile2D *fdIm0, TProfile2D *fdIm1, TProfile2D *fdIm2, TProfile2D *fdIm3, TProfile2D *fdIm4, TProfile2D *fdIm5, TProfile2D *fdIm6, TProfile2D *fdIm7);
+
void Calculate();
private:
AliCumulantsFunctions(const AliCumulantsFunctions& fun);
AliCumulantsFunctions& operator=(const AliCumulantsFunctions& fun);
- TProfile2D *fIG; //
- TProfile3D *fDGRe; //
- TProfile3D *fDGIm; //
- TH1D *fifr; //integrated flow final results
- TH1D *fdfr2; //differential flow final results //to be improved
- TH1D *fdfr4; //differential flow final results //to be improved
- TH1D *fdfr6; //differential flow final results //to be improved
- TH1D *fdfr8; //differential flow final results //to be improved
- Double_t fAvMult; //avarage multiplicity
+ TProfile2D *fIntGenFun; //avarage value of generating function for int. flow
+ TProfile3D *fDiffGenFunRe; //avarage value of generating function for diff. flow (real part)
+ TProfile3D *fDiffGenFunIm; //avarage value of generating function for diff. flow (imaginary part)
+ TH1D *fifr; //integrated flow final results
+ TH1D *fdfr2; //differential flow final results
+ TH1D *fdfr4; //differential flow final results
+ TH1D *fdfr6; //differential flow final results
+ TH1D *fdfr8; //differential flow final results
+
+ TProfile *fAvMult; //avarage selected multiplicity for int. flow
+ TProfile *fQVector; //avarage values of Q-vector components (1st bin: <Q_x>, 2nd bin: <Q_y>, 3rd bin: <(Q_x)^2>, 4th bin: <(Q_y)^2>)
+
+ TH1D *fQDistrib; //q-distribution
+
+ TProfile2D *fdRe0,*fdRe1,*fdRe2,*fdRe3,*fdRe4,*fdRe5,*fdRe6,*fdRe7;//differential flow
+ TProfile2D *fdIm0,*fdIm1,*fdIm2,*fdIm3,*fdIm4,*fdIm5,*fdIm6,*fdIm7;//differential flow
ClassDef(AliCumulantsFunctions, 0);
};
+
+//================================================================================================================
+
#endif
+/*************************************************************************
+* Copyright(c) 1998-2008, 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. *
+**************************************************************************/
+
+/********************************
+ * flow analysis with cumulants *
+ * *
+ * author: Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ *******************************/
+
#define AliFlowAnalysisWithCumulants_cxx
#include "Riostream.h"
#include "TFile.h"
#include "TList.h" //NEW
#include "TParticle.h"
-
-
-
-
-
-
-
+#include "TRandom3.h"
#include "TProfile.h"
#include "TProfile2D.h"
#include "TProfile3D.h"
-
-
-
-
-
-
#include "AliFlowEventSimple.h"
#include "AliFlowTrackSimple.h"
#include "AliFlowAnalysisWithCumulants.h"
#include "AliFlowCumuConstants.h"
class TH1;
-class TH3;
class TGraph;
class TPave;
class TLatex;
class TCanvas;
class TSystem;
class TROOT;
-
class AliFlowVector;
class TVector;
-//*******************************
-// flow analysis with cumulants *
-// author: Ante Bilandzic *
-// email: anteb@nikhef.nl *
-//*******************************
+//================================================================================================================
ClassImp(AliFlowAnalysisWithCumulants)
-//________________________________________________________________________
-
AliFlowAnalysisWithCumulants::AliFlowAnalysisWithCumulants():
- fTrack(NULL),
- //fHistFileName("cumulants.root"),//NEW: In the old version theis was not commented out
- //fHistFilefHistFile(NULL),//NEW: In the old version theis was not commented out
- fHistList(NULL), //NEW
- fAvM(0),
- fR0(0),
- fPtMax(0),
- fPtMin(0),
- fBinWidth(0),
- fAvQx(0),
- fAvQy(0),
- fAvQ2x(0),
- fAvQ2y(0),
- fNumberOfEvents(0),
-
-
-
- fHistProAvM(NULL),
- fIntFlowResults(NULL),
- fDiffFlowResults2(NULL),//to be improved
- fDiffFlowResults4(NULL),
- fDiffFlowResults6(NULL),
- fDiffFlowResults8(NULL),
- fIntFlowGenFun(NULL),
- fDiffFlowGenFunRe(NULL),
- fDiffFlowGenFunIm(NULL),
-
-
-
- fCommonHists(NULL)
- //fCommonHistsRes2(NULL),//to be improved
- //fCommonHistsRes4(NULL),//to be improved
- //fCommonHistsRes6(NULL),//to be improved
- //fCommonHistsRes8(NULL)//to be improved
- {
- //constructor
-
- fHistList = new TList(); //NEW
-
- fR0=AliFlowCumuConstants::fgR0;
- fPtMax=AliFlowCommonConstants::GetPtMax();
- fPtMin=AliFlowCommonConstants::GetPtMin();
- fBinWidth=(fPtMax-fPtMin)/fgknBins;
-
- }
+ fTrack(NULL),
+ fHistList(NULL),
+ fR0(0),
+ fPtMax(0),
+ fPtMin(0),
+ fBinWidth(0),
+ fAvQx(0),
+ fAvQy(0),
+ fAvQ2x(0),
+ fAvQ2y(0),
+ fAvMultIntFlow(NULL),
+ fQVectorComponents(NULL),
+ fIntFlowResults(NULL),
+ fDiffFlowResults2(NULL),
+ fDiffFlowResults4(NULL),
+ fDiffFlowResults6(NULL),
+ fDiffFlowResults8(NULL),
+ fIntFlowGenFun(NULL),
+ fDiffFlowGenFunRe(NULL),
+ fDiffFlowGenFunIm(NULL),
+ fDiffFlowGenFunRe0(NULL),
+ fDiffFlowGenFunRe1(NULL),
+ fDiffFlowGenFunRe2(NULL),
+ fDiffFlowGenFunRe3(NULL),
+ fDiffFlowGenFunRe4(NULL),
+ fDiffFlowGenFunRe5(NULL),
+ fDiffFlowGenFunRe6(NULL),
+ fDiffFlowGenFunRe7(NULL),
+ fDiffFlowGenFunIm0(NULL),
+ fDiffFlowGenFunIm1(NULL),
+ fDiffFlowGenFunIm2(NULL),
+ fDiffFlowGenFunIm3(NULL),
+ fDiffFlowGenFunIm4(NULL),
+ fDiffFlowGenFunIm5(NULL),
+ fDiffFlowGenFunIm6(NULL),
+ fDiffFlowGenFunIm7(NULL),
+ fCommonHists(NULL),
+ fQDist(NULL)//q-distribution
+{
+ //constructor
+ fHistList = new TList();
+ fR0=AliFlowCumuConstants::fgR0;
+ fPtMax=AliFlowCommonConstants::GetPtMax();
+ fPtMin=AliFlowCommonConstants::GetPtMin();
+ fBinWidth=(fPtMax-fPtMin)/fgknBins;
+}
-AliFlowAnalysisWithCumulants::~AliFlowAnalysisWithCumulants(){
- //desctructor
- delete fHistList; //NEW
+AliFlowAnalysisWithCumulants::~AliFlowAnalysisWithCumulants()
+{
+ //desctructor
+ delete fHistList;
}
-
-//___________________________________________________________________________
-void AliFlowAnalysisWithCumulants::CreateOutputObjects(){
- //output histograms
+//================================================================================================================
+
+void AliFlowAnalysisWithCumulants::CreateOutputObjects()
+{
+ //various output histograms
- fHistProAvM = new TProfile("fHistProAvM","Avarage Multiplicity",1,0,1,0,100000);
- fHistProAvM->SetXTitle("");
- fHistProAvM->SetYTitle("Avarage Multiplicity");
- fHistList->Add(fHistProAvM);
-
- fIntFlowResults = new TH1D("fIntFlowResults","Integrated Flow From Cumulants",8,0,8);
- fIntFlowResults->SetXTitle("");
- fIntFlowResults->SetYTitle("Integrated Flow [%]");
- fHistList->Add(fIntFlowResults);
+ //avarage multiplicity
+ fAvMultIntFlow = new TProfile("fAvMultIntFlow","Avarage multiplicity of selected particles used for int. flow",1,0,1,0,100000);
+ fAvMultIntFlow->SetXTitle("");
+ fAvMultIntFlow->SetYTitle("<multiplicity>");
+ fHistList->Add(fAvMultIntFlow);
- //to be improved, I should store the results as CommonHistResults
- fDiffFlowResults2 = new TH1D("fDiffFlowResults2","v'_2/2{2}",fgknBins,fPtMin,fPtMax);
- fDiffFlowResults2->SetXTitle("pt [GeV]");
- fDiffFlowResults2->SetYTitle("Differential Flow [%]");
- fHistList->Add(fDiffFlowResults2);
+ //Q-vector stuff
+ fQVectorComponents = new TProfile("fQVectorComponents","Avarage of Q-vector components",4,0.,4.);
+ fQVectorComponents->SetXTitle("");
+ fQVectorComponents->SetYTitle("< >");
+ fHistList->Add(fQVectorComponents);
- //to be improved, I should store the results as CommonHistResults
- fDiffFlowResults4 = new TH1D("fDiffFlowResults4","v'_2/2{4}",fgknBins,fPtMin,fPtMax);
- fDiffFlowResults4->SetXTitle("pt [GeV]");
- fDiffFlowResults4->SetYTitle("Differential Flow [%]");
- fHistList->Add(fDiffFlowResults4);
+ //final results for integrated flow (v_n{2}, v_n{4},..., v_n{16}) from cumulants (by default n=2)
+ fIntFlowResults = new TH1D("fIntFlowResults","Integrated Flow From Cumulants",8,0,8);
+ fIntFlowResults->SetXTitle("");
+ fIntFlowResults->SetYTitle("Integrated Flow [%]");
+ fHistList->Add(fIntFlowResults);
- //to be improved, I should store the results as CommonHistResults
- fDiffFlowResults6 = new TH1D("fDiffFlowResults6","v'_2/2{6}",fgknBins,fPtMin,fPtMax);
- fDiffFlowResults6->SetXTitle("pt [GeV]");
- fDiffFlowResults6->SetYTitle("Differential Flow [%]");
- fHistList->Add(fDiffFlowResults6);
-
- //to be improved, I should store the results as CommonHistResults
- fDiffFlowResults8 = new TH1D("fDiffFlowResults8","v'_2/2{8}",fgknBins,fPtMin,fPtMax);
- fDiffFlowResults8->SetXTitle("pt [GeV]");
- fDiffFlowResults8->SetYTitle("Differential Flow [%]");
- fHistList->Add(fDiffFlowResults8);
+ //final results for differential flow v_p/n{2} (by default p=n=2)
+ fDiffFlowResults2 = new TH1D("fDiffFlowResults2","v'_2/2{2}",fgknBins,fPtMin,fPtMax);
+ fDiffFlowResults2->SetXTitle("pt [GeV]");
+ fDiffFlowResults2->SetYTitle("Differential Flow [%]");
+ fHistList->Add(fDiffFlowResults2);
+
+ //final results for differential flow v_p/n{4} (by default p=n=2)
+ fDiffFlowResults4 = new TH1D("fDiffFlowResults4","v'_2/2{4}",fgknBins,fPtMin,fPtMax);
+ fDiffFlowResults4->SetXTitle("pt [GeV]");
+ fDiffFlowResults4->SetYTitle("Differential Flow [%]");
+ fHistList->Add(fDiffFlowResults4);
- fIntFlowGenFun = new TProfile2D("fIntFlowGenFun","G[p][q]",8,0.,8.,17,0.,17.);//to be improved (z_down and z_up)
- fIntFlowGenFun->SetXTitle("p");
- fIntFlowGenFun->SetYTitle("q");
- fHistList->Add(fIntFlowGenFun);
-
- fDiffFlowGenFunRe = new TProfile3D("fDiffFlowGenFunRe","Re(D[p][q])",fgknBins,fPtMin/fBinWidth,fPtMax/fBinWidth,8,0.,8.,17,0.,17.);
- fHistList->Add(fDiffFlowGenFunRe);
+ //final results for differential flow v_p/n{6} (by default p=n=2)
+ fDiffFlowResults6 = new TH1D("fDiffFlowResults6","v'_2/2{6}",fgknBins,fPtMin,fPtMax);
+ fDiffFlowResults6->SetXTitle("pt [GeV]");
+ fDiffFlowResults6->SetYTitle("Differential Flow [%]");
+ fHistList->Add(fDiffFlowResults6);
- fDiffFlowGenFunIm = new TProfile3D("fDiffFlowGenFunIm","Im(D[p][q])",fgknBins,fPtMin/fBinWidth,fPtMax/fBinWidth,8,0.,8.,17,0.,17.);
- fHistList->Add(fDiffFlowGenFunIm);
+ //final results for differential flow v_p/n{8} (by default p=n=2)
+ fDiffFlowResults8 = new TH1D("fDiffFlowResults8","v'_2/2{8}",fgknBins,fPtMin,fPtMax);
+ fDiffFlowResults8->SetXTitle("pt [GeV]");
+ fDiffFlowResults8->SetYTitle("Differential Flow [%]");
+ fHistList->Add(fDiffFlowResults8);
-
- /*
- fCommonHistsRes2 = new AliFlowCommonHistResults("Cumulants2");//to be improved
- fCommonHistsRes4 = new AliFlowCommonHistResults("Cumulants4");//to be improved
- fCommonHistsRes6 = new AliFlowCommonHistResults("Cumulants6");//to be improved
- fCommonHistsRes8 = new AliFlowCommonHistResults("Cumulants8");//to be improved
+ //avarage of the generating function for integrated flow <G[p][q]>
+ fIntFlowGenFun = new TProfile2D("fIntFlowGenFun","<G[p][q]>",fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
+ fIntFlowGenFun->SetXTitle("p");
+ fIntFlowGenFun->SetYTitle("q");
+ fHistList->Add(fIntFlowGenFun);
+ /*
+ //avarage of the real part of generating function for differential flow <Re(D[b][p][q])>
+ fDiffFlowGenFunRe = new TProfile3D("fDiffFlowGenFunRe","<Re(D[b][p][q])>",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe->SetXTitle("b");
+ fDiffFlowGenFunRe->SetYTitle("p");
+ fDiffFlowGenFunRe->SetZTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe);
- fHistList->Add(fCommonHistsRes2->GetHistList()); //NEW //to be improved
- fHistList->Add(fCommonHistsRes4->GetHistList()); //NEW //to be improved
- fHistList->Add(fCommonHistsRes6->GetHistList()); //NEW //to be improved
- fHistList->Add(fCommonHistsRes8->GetHistList()); //NEW //to be improved
+ //avarage of the imaginary part of generating function for differential flow <Im(D[b][p][q])>
+ fDiffFlowGenFunIm = new TProfile3D("fDiffFlowGenFunIm","<Im(D[b][p][q])>",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm->SetXTitle("b");
+ fDiffFlowGenFunIm->SetYTitle("p");
+ fDiffFlowGenFunIm->SetZTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm);
*/
-
- //control histograms
+
+ fDiffFlowGenFunRe0 = new TProfile2D("fDiffFlowGenFunRe0","Re(<D[b][0][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe0->SetXTitle("b");
+ fDiffFlowGenFunRe0->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe0);
+
+ fDiffFlowGenFunIm0 = new TProfile2D("fDiffFlowGenFunIm0","Im(<D[b][0][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm0->SetXTitle("b");
+ fDiffFlowGenFunIm0->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm0);
+
+ fDiffFlowGenFunRe1 = new TProfile2D("fDiffFlowGenFunRe1","Re(<D[b][1][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe1->SetXTitle("b");
+ fDiffFlowGenFunRe1->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe1);
+
+ fDiffFlowGenFunIm1 = new TProfile2D("fDiffFlowGenFunIm1","Im(<D[b][1][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm1->SetXTitle("b");
+ fDiffFlowGenFunIm1->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm1);
+
+ fDiffFlowGenFunRe2 = new TProfile2D("fDiffFlowGenFunRe2","Re(<D[b][2][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe2->SetXTitle("b");
+ fDiffFlowGenFunRe2->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe2);
+
+ fDiffFlowGenFunIm2 = new TProfile2D("fDiffFlowGenFunIm2","Im(<D[b][2][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm2->SetXTitle("b");
+ fDiffFlowGenFunIm2->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm2);
+
+ fDiffFlowGenFunRe3 = new TProfile2D("fDiffFlowGenFunRe3","Re(<D[b][3][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe3->SetXTitle("b");
+ fDiffFlowGenFunRe3->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe3);
+
+ fDiffFlowGenFunIm3 = new TProfile2D("fDiffFlowGenFunIm3","Im(<D[b][3][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm3->SetXTitle("b");
+ fDiffFlowGenFunIm3->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm3);
+
+ fDiffFlowGenFunRe4 = new TProfile2D("fDiffFlowGenFunRe4","Re(<D[b][4][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe4->SetXTitle("b");
+ fDiffFlowGenFunRe4->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe4);
+
+ fDiffFlowGenFunIm4 = new TProfile2D("fDiffFlowGenFunIm4","Im(<D[b][4][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm4->SetXTitle("b");
+ fDiffFlowGenFunIm4->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm4);
+
+ fDiffFlowGenFunRe5 = new TProfile2D("fDiffFlowGenFunRe5","Re(<D[b][5][q]>)",fgkQmax,0.,(Double_t)fgkQmax,fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth));
+ fDiffFlowGenFunRe5->SetXTitle("b");
+ fDiffFlowGenFunRe5->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe5);
+
+ fDiffFlowGenFunIm5 = new TProfile2D("fDiffFlowGenFunIm5","Im(<D[b][5][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm5->SetXTitle("b");
+ fDiffFlowGenFunIm5->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm5);
+
+ fDiffFlowGenFunRe6 = new TProfile2D("fDiffFlowGenFunRe6","Re(<D[b][6][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe6->SetXTitle("b");
+ fDiffFlowGenFunRe6->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe6);
+
+ fDiffFlowGenFunIm6 = new TProfile2D("fDiffFlowGenFunIm6","Im(<D[b][6][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm6->SetXTitle("b");
+ fDiffFlowGenFunIm6->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm6);
+
+ fDiffFlowGenFunRe7 = new TProfile2D("fDiffFlowGenFunRe7","Re(<D[b][7][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunRe7->SetXTitle("b");
+ fDiffFlowGenFunRe7->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunRe7);
+
+ fDiffFlowGenFunIm7 = new TProfile2D("fDiffFlowGenFunIm7","Im(<D[b][7][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
+ fDiffFlowGenFunIm7->SetXTitle("b");
+ fDiffFlowGenFunIm7->SetYTitle("q");
+ fHistList->Add(fDiffFlowGenFunIm7);
+
+ //common control histograms
fCommonHists = new AliFlowCommonHist("Cumulants");
- fHistList->Add(fCommonHists->GetHistList());
+ fHistList->Add(fCommonHists->GetHistList());
-}
-
-//________________________________________________________________________
-void AliFlowAnalysisWithCumulants::Make(AliFlowEventSimple* anEvent) {
- //running over data
+ //q-distribution
+ fQDist = new TH1D("fQDist","q-distribution",100,0,10);
+ fQDist->SetXTitle("q=Q/#sqrt{M}");
+ fQDist->SetYTitle("dN_{event}/dq");
+ fHistList->Add(fQDist);
- //---------------------------------------------------------
- //fill the common control histograms
- fCommonHists->FillControlHistograms(anEvent);
- //---------------------------------------------------------
+}//end of CreateOutputObjects()
+
+//================================================================================================================
+
+void AliFlowAnalysisWithCumulants::Make(AliFlowEventSimple* anEvent)
+{
+ //running over data
+ //fill the common control histograms
+ fCommonHists->FillControlHistograms(anEvent);
- //---------------------------------------------------------
- //initializing the generating function for integrated flow
- Double_t fG[fgkPmax][fgkQmax];
- for(Int_t p=0;p<fgkPmax;p++){
- for(Int_t q=0;q<fgkQmax;q++){
- fG[p][q]=1.;
- }
+ //initializing the generating function G[p][q] for integrated flow
+ Double_t G[fgkPmax][fgkQmax];
+ for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ G[p][q]=1.;
}
- //---------------------------------------------------------
-
-
- //---------------------------------------------------------
- //Q vector stuff
- AliFlowVector fQVector;
- fQVector.Set(0.,0.);
- fQVector.SetMult(0);
-
- fQVector=anEvent->GetQ();//get the Q vector for this event
+ }
- fAvQx+=fQVector.X();
- fAvQy+=fQVector.Y();
- fAvQ2x+=pow(fQVector.X(),2.);
- fAvQ2y+=pow(fQVector.Y(),2.);
- //----------------------------------------------------------
-
- Int_t nPrim = anEvent->NumberOfTracks();
- Int_t fEventNSelTracksIntFlow = anEvent->GetEventNSelTracksIntFlow();
- Int_t fSelTracksIntFlow = 0;
-
- cout<<"Number of input tracks for cumulant analysis: "<<nPrim<<endl;
- cout<<"Number of selected tracks for cumulant analysis: "<<fEventNSelTracksIntFlow<<endl;
+ Int_t nPrim = anEvent->NumberOfTracks();//total multiplicity
+ Int_t nEventNSelTracksIntFlow = anEvent->GetEventNSelTracksIntFlow();//selected multiplicity (parrticles used for int. flow)
+ Int_t nSelTracksIntFlow = 0;//cross-checking the selected multiplicity
- //------------------------------------------------------------------------------------
- //STARTING THE FIRST LOOP (CALCULATING THE GENERATING FUNCTION FOR INTEGRATED FLOW)
- for(Int_t i=0;i<nPrim;i++){
- fTrack=anEvent->GetTrack(i);
- if(fTrack&&fTrack->UseForIntegratedFlow()){
- fSelTracksIntFlow++;
- for(Int_t p=0;p<fgkPmax;p++){
- for(Int_t q=0;q<fgkQmax;q++){
- fG[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/fEventNSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax));
- }
- }
+ //first loop over data: evaluating the generating function G[p][q] for integrated flow
+ for(Int_t i=0;i<nPrim;i++){
+ fTrack=anEvent->GetTrack(i);
+ if(fTrack&&fTrack->UseForIntegratedFlow()){
+ nSelTracksIntFlow++;
+ for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ G[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax));
}
- }// ENDING THE FIRST LOOP OVER TRACKS
- //------------------------------------------------------------------------------------
-
-
- fHistProAvM->Fill(0.,fSelTracksIntFlow,1);
- fAvM=fHistProAvM->GetBinContent(1);
-
-
- //------------------------------------------------------------------------------------
- //STORING THE VALUE OF GENERATING FUNCTION FOR INTEGRATED FLOW INTO 2D PROFILE
- for(Int_t p=0;p<fgkPmax;p++){
- for(Int_t q=0;q<fgkQmax;q++){
- fIntFlowGenFun->Fill((Double_t)p,(Double_t)q,fG[p][q],1);
- }
+ }
}
- //------------------------------------------------------------------------------------
-
-
- //------------------------------------------------------------------------------------
- //STARTING THE SECOND LOOP (CALCULATING THE GENERATING FUNCTION FOR DIFFERENTIAL FLOW)
- // Remark 0: generating function for diff. flow is complex number, I need to calcuate separately real and imaginary part
- // Remark 1: note that I need here fG[p][q], the value of generating function for integrated flow for the CURRENT event
- // Remark 2: results are immediately stored in two 3D profiles, one for real and one for imaginary part
- for(Int_t i=0;i<nPrim;i++){
- fTrack=anEvent->GetTrack(i);
- if (fTrack && fTrack->UseForDifferentialFlow()){
- for(Int_t p=0;p<fgkPmax;p++){
- for(Int_t q=0;q<fgkQmax;q++){
- fDiffFlowGenFunRe->Fill(fTrack->Pt()/fBinWidth,(Double_t)p,(Double_t)q,fG[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/
- (1.+(2.*fR0*sqrt(p+1.)/fSelTracksIntFlow) *
- cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
- fDiffFlowGenFunIm->Fill(fTrack->Pt()/fBinWidth,(Double_t)p,(Double_t)q,fG[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/
- (1.+(2.*fR0*sqrt(p+1.)/fSelTracksIntFlow) *
- cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
- }
+ }//ending the first loop over data
+
+ //storing the value of G[p][q] in 2D profile in order to get automatically the avarage <G[p][q]>
+ for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ fIntFlowGenFun->Fill((Double_t)p,(Double_t)q,G[p][q],1);
+ }
+ }
+
+ //storing the selected multiplicity (if fAvMultIntFlow is not filled here then you had wrongly selected the particles used for integrated flow)
+ if(nSelTracksIntFlow==nEventNSelTracksIntFlow)
+ {
+ fAvMultIntFlow->Fill(0.,nSelTracksIntFlow,1);
+ }
+
+ //calculating Q-vector of event (needed for errors)
+ AliFlowVector fQVector;
+ fQVector.Set(0.,0.);
+ fQVector.SetMult(0);
+ fQVector=anEvent->GetQ(); //get the Q vector for this event
+ fQVectorComponents->Fill(0.,fQVector.X(),1); //in the 1st bin fill Q_x
+ fQVectorComponents->Fill(1.,fQVector.Y(),1); //in the 2nd bin fill Q_y
+ fQVectorComponents->Fill(2.,pow(fQVector.X(),2.),1); //in the 3rd bin fill (Q_x)^2
+ fQVectorComponents->Fill(3.,pow(fQVector.Y(),2.),1); //in the 4th bin fill (Q_y)^2
+
+ //q-distribution
+ Double_t qDist = fQVector.Mod()/sqrt(fQVector.GetMult());
+ fQDist->Fill(qDist,1);
+
+ /*
+ //two 3D profiles for differential flow
+ //second loop over data: evaluating the generating function D[b][p][q] for differential flow
+ //remark 0: D[b][p][q] is a complex number => real and imaginary part are calculated separately
+ //remark 1: note that below G[p][q] is needed, the value of generating function for integrated flow for the CURRENT event
+ //remark 2: results are stored in two 3D profiles in order to automatically get <Re(D[b][p][q])> and <Im(D[b][p][q])>
+ for(Int_t i=0;i<nPrim;i++){
+ fTrack=anEvent->GetTrack(i);
+ if (fTrack && fTrack->UseForDifferentialFlow()){
+ for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ //real part
+ fDiffFlowGenFunRe->Fill(fTrack->Pt()/fBinWidth,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm->Fill(fTrack->Pt()/fBinWidth,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
}
- }
- }// ENDING THE SECOND LOOP OVER TRACKS
- //------------------------------------------------------------------------------------
-
+ }
+ }
+ }//ending the second loop over data
+ */
+
+ //sixteen 2D profiles for differential flow
+ for(Int_t i=0;i<nPrim;i++){
+ fTrack=anEvent->GetTrack(i);
+ if (fTrack && fTrack->UseForDifferentialFlow()){
+ //for(Int_t p=0;p<fgkPmax;p++){
+ for(Int_t q=0;q<fgkQmax;q++){
+ //real part
+ fDiffFlowGenFunRe0->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[0][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(0.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm0->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[0][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(0.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe1->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[1][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(1.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm1->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[1][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(1.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe2->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[2][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(2.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm2->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[2][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(2.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe3->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[3][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(3.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm3->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[3][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(3.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe4->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[4][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(4.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm4->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[4][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(4.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe5->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[5][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(5.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm5->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[5][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(5.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe6->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[6][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(6.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm6->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[6][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(6.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //-----------------------------------------------------------------------
+ //real part
+ fDiffFlowGenFunRe7->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[7][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(7.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ //imaginary part
+ fDiffFlowGenFunIm7->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[7][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(7.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
+ }
+ //}
+ }
+ }//ending the second loop over data
+
}//end of Make()
-//________________________________________________________________________
-void AliFlowAnalysisWithCumulants::Finish(){
- //final results
-}
-
-
-
-
-
-
-
-
-
-
-
-
-
+//================================================================================================================
+void AliFlowAnalysisWithCumulants::Finish()
+{
+ //not needed for the time being...
+}
+/*
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved.
+ * See cxx source for full Copyright notice
+ * $Id$
+ */
+
+/********************************
+ * flow analysis with cumulants *
+ * *
+ * author: Ante Bilandzic *
+ * (anteb@nikhef.nl) *
+ *******************************/
+
#ifndef AliFlowAnalysisWithCumulants_H
#define AliFlowAnalysisWithCumulants_H
-//*******************************
-// flow analysis with cumulants *
-// author: Ante Bilandzic *
-// email: anteb@nikhef.nl *
-//*******************************
-
#include "AliFlowCommonConstants.h"
#include "AliFlowCumuConstants.h"
-class TH1;
-class TH3D;
-
class TObjArray;
+class TList;
+class TFile;
-
+class TH1;
class TProfile;
class TProfile2D;
class TProfile3D;
-
-class TList;
-class TFile;
-
class AliFlowEventSimple;
class AliFlowTrackSimple;
class AliFlowCommonHist;
class AliFlowCommonHistResults;
class AliFlowVector;
-class AliFlowAnalysisWithCumulants {
+//================================================================================================================
+
+class AliFlowAnalysisWithCumulants{
public:
AliFlowAnalysisWithCumulants();
- virtual ~AliFlowAnalysisWithCumulants();
+ virtual ~AliFlowAnalysisWithCumulants();
virtual void CreateOutputObjects();
virtual void Make(AliFlowEventSimple* anEvent);
virtual void Finish();
-
- // Output
- // Output
- TList* GetHistList() const { return this->fHistList ; } // Gets output histogram list //NEW
-
- private:
- AliFlowAnalysisWithCumulants(const AliFlowAnalysisWithCumulants& aAnalysis);
- AliFlowAnalysisWithCumulants& operator=(const AliFlowAnalysisWithCumulants& aAnalysis);
-
- AliFlowTrackSimple* fTrack;//track
- static const Int_t fgkQmax=AliFlowCumuConstants::kQmax;//needed for numerics
- static const Int_t fgkPmax=AliFlowCumuConstants::kPmax;//needed for numerics
- static const Int_t fgkFlow=AliFlowCumuConstants::kFlow;//integrated flow coefficient to be calculated
- static const Int_t fgkMltpl=AliFlowCumuConstants::kMltpl;//the multiple in p=m*n (diff. flow)
- static const Int_t fgknBins=100;//number of pt bins
+ TList* GetHistList() const {return this->fHistList;} //output histogram list
- //TString fHistFileName; //! The output file name // NEW: In the old version this was not commented out
- //TFile* fHistFile; // histogram file for Cumulants // NEW: In the old version this was not commented out
-
- TList* fHistList; //list to hold all output histograms //NEW
-
- Double_t fAvM;//avarage SELECTED multiplicity
-
- Double_t fR0;//needed for numerics
- Double_t fPtMax;//maximum pt
- Double_t fPtMin;//minimum pt
- Double_t fBinWidth;//width of pt bin (in GeV)
+ private:
+ AliFlowAnalysisWithCumulants(const AliFlowAnalysisWithCumulants& afawc);
+ AliFlowAnalysisWithCumulants& operator=(const AliFlowAnalysisWithCumulants& afawc);
+ AliFlowTrackSimple* fTrack; //track
+ static const Int_t fgkQmax=AliFlowCumuConstants::kQmax; //needed for numerics
+ static const Int_t fgkPmax=AliFlowCumuConstants::kPmax; //needed for numerics
+ static const Int_t fgkFlow=AliFlowCumuConstants::kFlow; //integrated flow coefficient to be calculated
+ static const Int_t fgkMltpl=AliFlowCumuConstants::kMltpl; //the multiple in p=m*n (diff. flow)
+ static const Int_t fgknBins=100; //number of pt bins
+ TList* fHistList; //list to hold all output histograms
+
+ Double_t fR0; //needed for numerics
+ Double_t fPtMax; //maximum pt
+ Double_t fPtMin; //minimum pt
+ Double_t fBinWidth; //width of pt bin (in GeV)
- Double_t fAvQx;//<Q_x>
- Double_t fAvQy;//<Q_y>
- Double_t fAvQ2x;//<(Q_x)^2>
- Double_t fAvQ2y;//<(Q_y)^2>
-
- Int_t fNumberOfEvents;//number of events
-
-
- TProfile* fHistProAvM; //avarage multiplicity
-
+ Double_t fAvQx; //<Q_x>
+ Double_t fAvQy; //<Q_y>
+ Double_t fAvQ2x; //<(Q_x)^2>
+ Double_t fAvQ2y; //<(Q_y)^2>
+
+ TProfile* fAvMultIntFlow; //avarage selected multiplicity
+
+ TProfile* fQVectorComponents; //avarages of Q-vector components (1st bin: <Q_x>, 2nd bin: <Q_y>, 3rd bin: <(Q_x)^2>, 4th bin: <(Q_y)^2>)
+
TH1D* fIntFlowResults; //integrated flow final results
- TH1D* fDiffFlowResults2; //differential flow final results //to be improved
- TH1D* fDiffFlowResults4; //differential flow final results //to be improved
- TH1D* fDiffFlowResults6; //differential flow final results //to be improved
- TH1D* fDiffFlowResults8; //differential flow final results //to be improved
+ TH1D* fDiffFlowResults2; //differential flow final results (2nd order estimate)
+ TH1D* fDiffFlowResults4; //differential flow final results (4th order estimate)
+ TH1D* fDiffFlowResults6; //differential flow final results (6th order estimate)
+ TH1D* fDiffFlowResults8; //differential flow final results (8th order estimate)
- TProfile2D* fIntFlowGenFun; // avarage of the generating function for integrated flow
- TProfile3D* fDiffFlowGenFunRe; // avarage of the generating function for differential flow (real part)
- TProfile3D* fDiffFlowGenFunIm; // avarage of the generating function for differential flow (imaginary part)
-
+ TProfile2D* fIntFlowGenFun; //avarage of the generating function for integrated flow
+ TProfile3D* fDiffFlowGenFunRe; //avarage of the generating function for differential flow (real part)
+ TProfile3D* fDiffFlowGenFunIm; //avarage of the generating function for differential flow (imaginary part)
- AliFlowCommonHist* fCommonHists;//control histograms
- //AliFlowCommonHistResults *fCommonHistsRes2, *fCommonHistsRes4, *fCommonHistsRes6, *fCommonHistsRes8;//final results//to be improved
+ TProfile2D *fDiffFlowGenFunRe0,*fDiffFlowGenFunRe1,*fDiffFlowGenFunRe2,*fDiffFlowGenFunRe3;//differential flow
+ TProfile2D *fDiffFlowGenFunRe4,*fDiffFlowGenFunRe5,*fDiffFlowGenFunRe6,*fDiffFlowGenFunRe7;//differential flow
+ TProfile2D *fDiffFlowGenFunIm0,*fDiffFlowGenFunIm1,*fDiffFlowGenFunIm2,*fDiffFlowGenFunIm3;//differential flow
+ TProfile2D *fDiffFlowGenFunIm4,*fDiffFlowGenFunIm5,*fDiffFlowGenFunIm6,*fDiffFlowGenFunIm7;//differential flow
+
+ AliFlowCommonHist* fCommonHists; //common control histograms
- Double_t fAvG[fgkPmax][fgkQmax];//avarage of the generating function used for integrated flow
- Int_t fBinEventEntries[fgknBins];//counts how many events have at least 1 particle in particular bin
- Int_t fBinNoOfParticles[fgknBins];//number of particles per bin
- Double_t fBinMeanPt[fgknBins];//mean pt per bin
- Double_t fBinEventDRe[fgknBins][fgkPmax][fgkQmax];//real part of the generating function used for differential flow
- Double_t fBinEventDIm[fgknBins][fgkPmax][fgkQmax];//imaginary part of the generating function used for differential flow
-
+ TH1D* fQDist; //q-distribution
+
ClassDef(AliFlowAnalysisWithCumulants, 0);
};
+
+//================================================================================================================
+
#endif
-//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;
-
-
// 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 runAliAnalysisTaskCumulants(Int_t nRuns = -1, TString type = "ESD", const Char_t* dataDir="/data/alice2/ante/ab", Int_t offset = 0)
+void runAliAnalysisTaskCumulants(Int_t nRuns =-1, TString type = "ESD", const Char_t* dataDir="/data/alice2/ante/ab2", Int_t offset = 0)
//void runAliAnalysisTaskCumulants(Int_t nRuns = -1, TString type = "ESD", const Char_t* dataDir="/data/alice2/ab2", Int_t offset = 0)
//void runAliAnalysisTaskCumulants(Int_t nRuns = -1, TString type = "ESD", const Char_t* dataDir="/data/alice2/LHyquid3_rot", Int_t offset = 0)
-void runAliAnalysisTaskCumulants(Int_t nRuns = -1, TString type = "ESD", const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFix", Int_t offset = 0)
+//void runAliAnalysisTaskCumulants(Int_t nRuns = 2, TString type = "MC", const Char_t* dataDir="/Users/snelling/alice_data/TherminatorFix", Int_t offset = 0)
{
TStopwatch timer;
timer.Start();
cerr<<"libESD loaded..."<<endl;
gSystem->Load("libANALYSIS.so");
cerr<<"libANALYSIS.so loaded..."<<endl;
- gSystem->Load("libCORRFW.so");
- cerr<<"libCORRFW.so loaded..."<<endl;
gSystem->Load("libPWG2flow.so");
cerr<<"libPWG2flow.so loaded..."<<endl;
-
-
-
-
-//____________________________________________//
- //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);
-
-
// create the TChain. CreateESDChain() is defined in CreateESDChain.C
TChain* chain = CreateESDChain(dataDir, nRuns, offset);
//____________________________________________//
// 1st cumulant task
AliAnalysisTaskCumulants *task1 = new AliAnalysisTaskCumulants("TaskCumulants");
- task1->SetAnalysisType("ESD");
- task1->SetCFManager1(cfmgr1);
- task1->SetCFManager2(cfmgr2);
+ task1->SetAnalysisType(type);
mgr->AddTask(task1);
// Create containers for input/output
delete iter;
}
+
const Int_t maxnsigmatovertex2 = 3;
- void runProofCumulants(const Char_t* data="/PWG2/akisiel/LHC500C2030", Int_t nRuns=100, Int_t offset=0) {
+ void runProofCumulants(const Char_t* data="/PWG2/akisiel/LHC500C2030", Int_t nRuns=4, Int_t offset=0) {
+ //void runProofCumulants(const Char_t* data="/PWG2/akisiel/LHC500C0005", Int_t nRuns=-1, Int_t offset=0){
+
+ //void runProofCumulants(const Char_t* data="/PWG2/pganoti/Pythia6At10TeV_05T_b", Int_t nRuns=1000, Int_t offset=0){
+ //void runProofCumulants(const Char_t* data="/PWG2/jgrosseo/sim_1600XX_esd", Int_t nRuns=-1, Int_t offset=0){
+
+ //void runProofCumulants(const Char_t* data="/COMMON/COMMON/run15035_PbPb", Int_t nRuns=-1, Int_t offset=0){
+
+ //void runProofCumulants(const Char_t* data="/PWG2/hricaud/LHC07f_160033DataSet", Int_t nRuns=-1, Int_t offset=0){
+
+ //void runProofCumulants(const Char_t* data="/PWG2/hricaud/LHC07f_160038_root_archiveDataSet", Int_t nRuns=-1, Int_t offset=0){
+
+ //void runProofCumulants(const Char_t* data="/PWG2/belikov/40825", Int_t nRuns=-1, Int_t offset=0){
+
+
+
TStopwatch timer;
timer.Start();
printf("*** Connect to PROOF ***\n");
- // 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");
+ TProof::Open("abilandz@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");
+
//____________________________________________//
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 ;
TObjArray* fPIDCutList2 = new TObjArray(0) ;
fPIDCutList2->AddLast(cutPID2);
- printf("CREATE INTERFACE AND CUTS\n");
+printf("CREATE INTERFACE AND CUTS\n");
AliCFManager* cfmgr1 = new AliCFManager();
cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1);
//cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList);
//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); }
-
- //____________________________________________//
- // Make the analysis manager
- AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
- AliESDInputHandler* esdH = new AliESDInputHandler;
- esdH->SetInactiveBranches("FMD CaloCluster");
- mgr->SetInputEventHandler(esdH);
+
+
+
+
+
+
+
+
//____________________________________________//
// 1st Pt task
AliAnalysisTaskCumulants *task1 = new AliAnalysisTaskCumulants("TaskCumulants");
- task1->SetAnalysisType("ESD");
+ 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);
- // AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("clist1", TList::Class(),AliAnalysisManager::kOutputContainer,"OutputFromCumulantAnlysisESD.root");
- AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("clist1", TList::Class(),AliAnalysisManager::kOutputContainer,"outputFromCumulantAnalysisESD.root");
+ AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
+ TString outputName = "outputFromCumulantAnalysis";
+ outputName+=type;
+ outputName+=".root";
+ AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("clist1", TList::Class(),AliAnalysisManager::kOutputContainer,outputName);
+
+ //added according to Andrei:
+ //coutput1->SetSpecialOutput();
+
//____________________________________________//
mgr->ConnectInput(task1,0,cinput1);
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff