#include "TCanvas.h"
#include "TStyle.h"
#include "TRandom.h"
+#include "THashList.h"
#include "AliAnalysisManager.h"
#include "AliMCEventHandler.h"
fCentNMixed(10),
fNEMRPBins(9),
fPeriod(period),
+ fMaxAbsVertexZ(10.),
+ fManualV0EPCalc(false),
fOutputContainer(0x0),
fNonLinCorr(0),
fEvent(0x0),
{
const int nbins = 9;
Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
- fCentEdges = TArrayD(nbins+1, edges);
+ TArrayD centEdges(nbins+1, edges);
Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
- fCentNMixed = TArrayI(nbins, nMixed);
+ TArrayI centNMixed(nbins, nMixed);
+ SetCentralityBinning(centEdges, centNMixed);
for(Int_t i=0;i<kNCenBins;i++){
for(Int_t j=0;j<2; j++)
if(fOutputContainer != NULL){
delete fOutputContainer;
}
- fOutputContainer = new TList();
+ fOutputContainer = new THashList();
fOutputContainer->SetOwner(kTRUE);
//========QA histograms=======
//Single photon and pi0 spectrum
- const Int_t nPtPhot = 300 ;
- const Double_t ptPhotMax = 30 ;
+ const Int_t nPtPhot = 400 ;
+ const Double_t ptPhotMax = 40 ;
const Int_t nM = 500;
const Double_t mMin = 0.0;
const Double_t mMax = 1.0;
FillSelectedClusterHistograms();
LogProgress(8);
+ if( ! fCaloPhotonsPHOS->GetEntriesFast() )
+ return;
+ else
+ LogSelection(6, fInternalRunNumber);
+
+
// Step 9: Consider pi0 (photon/cluster) pairs.
ConsiderPi0s();
LogProgress(9);
// // hTotSelEvents->Draw();
// }
//________________________________________________________________________
-void AliAnalysisTaskPi0Flow::SetCentralityBinning(const TArrayD& edges)
+void AliAnalysisTaskPi0Flow::SetCentralityBinning(const TArrayD& edges, const TArrayI& nMixed)
{
// Define centrality bins by their edges
-
- int last = edges.GetSize()-1;
- if( edges.At(0) < 0.)
- AliError("lower edge less then 0");
- if( 90. < edges.At(last) )
- AliError("upper edge larger then 90.");
- for(int i=0; i<last-1; ++i)
- if(edges.At(i) > edges.At(i+1))
- AliError("edges are not sorted");
+ if( edges.At(0) < 0.) AliFatal("lower edge less then 0");
+ if( 90. < edges.At(edges.GetSize()-1) ) AliFatal("upper edge larger then 90.");
+ for(int i=0; i<edges.GetSize()-1; ++i)
+ if(edges.At(i) > edges.At(i+1)) AliFatal("edges are not sorted");
+ if( edges.GetSize() != nMixed.GetSize()+1) AliFatal("edges and nMixed don't have appropriate relative sizes");
fCentEdges = edges;
-}
-//________________________________________________________________________
-void AliAnalysisTaskPi0Flow::SetNMixedPerCentrality(const TArrayI& nMixed)
-{
- // Set number of mixed events for all centrality bins
-
fCentNMixed = nMixed;
}
+
+
//________________________________________________________________________
void AliAnalysisTaskPi0Flow::SetPHOSBadMap(Int_t mod, TH2I* badMapHist)
{
// Track Matching
Double_t dx=clu->GetTrackDx() ;
Double_t dz=clu->GetTrackDz() ;
- Bool_t cpvBit=kTRUE ; //No track matched by default
+ Bool_t cpvBit=kTRUE ; //No track matched by default. True means: not from charged, according to veto.
Bool_t cpvBit2=kTRUE ; //More Strict criterion
if( fEventESD ) {
snprintf(key,55,"hPi0Both_a07_cen%d",fCentBin) ;
FillHistogram(Form("hPi0Both_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
}
- if(fCentralityV0M>20.){
if(ph1->Module()==1 && ph2->Module()==1)
FillHistogram("hPi0M11",p12.M(),p12.Pt() );
else if(ph1->Module()==2 && ph2->Module()==2)
FillHistogram("hPi0M13",p12.M(),p12.Pt() );
else if(ph1->Module()==2 && ph2->Module()==3)
FillHistogram("hPi0M23",p12.M(),p12.Pt() );
- }
}
}
//_____________________________________________________________________________
void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x)const{
//FillHistogram
- TH1I * tmpI = dynamic_cast<TH1I*>(fOutputContainer->FindObject(key)) ;
- if(tmpI){
- tmpI->Fill(x) ;
- return ;
- }
- TH1F * tmpF = dynamic_cast<TH1F*>(fOutputContainer->FindObject(key)) ;
- if(tmpF){
- tmpF->Fill(x) ;
- return ;
- }
- TH1D * tmpD = dynamic_cast<TH1D*>(fOutputContainer->FindObject(key)) ;
- if(tmpD){
- tmpD->Fill(x) ;
- return ;
- }
- AliInfo(Form("can not find histogram <%s> ",key)) ;
+ TH1 * hist = dynamic_cast<TH1*>(fOutputContainer->FindObject(key)) ;
+ if(hist)
+ hist->Fill(x) ;
+ else
+ AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
}
//_____________________________________________________________________________
void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y)const{
//FillHistogram
- TObject * tmp = fOutputContainer->FindObject(key) ;
- if(!tmp){
- AliInfo(Form("can not find histogram <%s> ",key)) ;
- return ;
- }
- if(tmp->IsA() == TClass::GetClass("TH1F")){
- ((TH1F*)tmp)->Fill(x,y) ;
- return ;
- }
- if(tmp->IsA() == TClass::GetClass("TH2F")){
- ((TH2F*)tmp)->Fill(x,y) ;
- return ;
- }
- AliError(Form("Calling FillHistogram with 2 parameters for histo <%s> of type %s",key,tmp->IsA()->GetName())) ;
+ TH1 * th1 = dynamic_cast<TH1*> (fOutputContainer->FindObject(key));
+ if(th1)
+ th1->Fill(x, y) ;
+ else
+ AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
}
//_____________________________________________________________________________
void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
//Fills 1D histograms with key
- TObject * tmp = fOutputContainer->FindObject(key) ;
- if(!tmp){
- AliInfo(Form("can not find histogram <%s> ",key)) ;
- return ;
- }
- if(tmp->IsA() == TClass::GetClass("TH2F")){
- ((TH2F*)tmp)->Fill(x,y,z) ;
- return ;
+ TObject * obj = fOutputContainer->FindObject(key);
+
+ TH2 * th2 = dynamic_cast<TH2*> (obj);
+ if(th2) {
+ th2->Fill(x, y, z) ;
+ return;
}
- if(tmp->IsA() == TClass::GetClass("TH3F")){
- ((TH3F*)tmp)->Fill(x,y,z) ;
- return ;
+
+ TH3 * th3 = dynamic_cast<TH3*> (obj);
+ if(th3) {
+ th3->Fill(x, y, z) ;
+ return;
}
+
+ AliError(Form("can not find histogram (of instance TH2) <%s> ",key)) ;
}
//_____________________________________________________________________________
void AliAnalysisTaskPi0Flow::SetV0Calibration(){
// assigns: fMultV0, fV0Cpol, fV0Apol, fMeanQ, and fWidthQ
+ if ( ! fManualV0EPCalc ) {
+ if( fDebug >=2 )
+ AliInfo("Not setting V0Calibration, only needed for manual V0 EP Calculation");
+ return;
+ }
+
int runNumber = this->fRunNumber;
TString oadbfilename = "$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root";
return true; // reject
LogSelection(1, fInternalRunNumber);
- if ( TMath::Abs(fVertexVector.z()) > 10. )
+ if ( TMath::Abs(fVertexVector.z()) > fMaxAbsVertexZ )
return true; // reject
LogSelection(2, fInternalRunNumber);
void AliAnalysisTaskPi0Flow::EvalV0ReactionPlane(){
// set: fRPV0A and fRPV0C
- //VZERO data
- AliVVZERO* v0 = fEvent->GetVZEROData();
-
- //reset Q vector info
- Double_t Qxa2 = 0, Qya2 = 0;
- Double_t Qxc2 = 0, Qyc2 = 0;
-
- for (Int_t iv0 = 0; iv0 < 64; iv0++) {
- Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
- Float_t multv0 = v0->GetMultiplicity(iv0);
- if (iv0 < 32){ // V0C
- Qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
- Qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
- } else { // V0A
- Qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
- Qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
+ // Do Manual V0 EP Calculation
+ if ( fManualV0EPCalc )
+ {
+ //VZERO data
+ AliVVZERO* v0 = fEvent->GetVZEROData();
+
+ //reset Q vector info
+ Double_t Qxa2 = 0, Qya2 = 0;
+ Double_t Qxc2 = 0, Qyc2 = 0;
+
+ for (Int_t iv0 = 0; iv0 < 64; iv0++) {
+ Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
+ Float_t multv0 = v0->GetMultiplicity(iv0);
+ if (iv0 < 32){ // V0C
+ Qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
+ Qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
+ } else { // V0A
+ Qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
+ Qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
+ }
+ }
+
+ Int_t iC = -1;
+ // centrality bins
+ if(fCentralityV0M < 5) iC = 0;
+ else if(fCentralityV0M < 10) iC = 1;
+ else if(fCentralityV0M < 20) iC = 2;
+ else if(fCentralityV0M < 30) iC = 3;
+ else if(fCentralityV0M < 40) iC = 4;
+ else if(fCentralityV0M < 50) iC = 5;
+ else if(fCentralityV0M < 60) iC = 6;
+ else if(fCentralityV0M < 70) iC = 7;
+ else iC = 8;
+
+ //grab for each centrality the proper histo with the Qx and Qy to do the recentering
+ Double_t Qxamean2 = fMeanQ[iC][1][0];
+ Double_t Qxarms2 = fWidthQ[iC][1][0];
+ Double_t Qyamean2 = fMeanQ[iC][1][1];
+ Double_t Qyarms2 = fWidthQ[iC][1][1];
+
+ Double_t Qxcmean2 = fMeanQ[iC][0][0];
+ Double_t Qxcrms2 = fWidthQ[iC][0][0];
+ Double_t Qycmean2 = fMeanQ[iC][0][1];
+ Double_t Qycrms2 = fWidthQ[iC][0][1];
+
+ Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
+ Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
+ Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
+ Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
+
+ fRPV0A = TMath::ATan2(QyaCor2, QxaCor2)/2.;
+ fRPV0C = TMath::ATan2(QycCor2, QxcCor2)/2.;
}
- }
+ else // Use Official V0 EP Calculation.
+ {
+ AliEventplane *eventPlane = fEvent->GetEventplane();
+ if( ! eventPlane ) { AliError("Event has no event plane"); return; }
+ fRPV0A = eventPlane->GetEventplane("V0A", fEvent);
+ fRPV0C = eventPlane->GetEventplane("V0C", fEvent);
+ }
+
+ // Check that the A&C RP are within allowed range.
+ if( fDebug >= 3 && (fRPV0A<0 || fRPV0A>TMath::Pi() ) )
+ AliInfo(Form("RPV0A outside of permited range [0,pi]: %f, correcting", fRPV0A));
+ if( fDebug >= 3 && (fRPV0C<0 || fRPV0C>TMath::Pi() ) )
+ AliInfo(Form("RPV0C outside of permited range [0,pi]: %f, correcting", fRPV0C));
+ while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
+ while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
+ while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
+ while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
+
+ // Reaction plane histograms before flattening
+ if( fDebug >= 2 )
+ AliInfo(Form("V0 Reaction Plane before flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
- Int_t iC = -1;
- // centrality bins
- if(fCentralityV0M < 5) iC = 0;
- else if(fCentralityV0M < 10) iC = 1;
- else if(fCentralityV0M < 20) iC = 2;
- else if(fCentralityV0M < 30) iC = 3;
- else if(fCentralityV0M < 40) iC = 4;
- else if(fCentralityV0M < 50) iC = 5;
- else if(fCentralityV0M < 60) iC = 6;
- else if(fCentralityV0M < 70) iC = 7;
- else iC = 8;
-
- //grab for each centrality the proper histo with the Qx and Qy to do the recentering
- Double_t Qxamean2 = fMeanQ[iC][1][0];
- Double_t Qxarms2 = fWidthQ[iC][1][0];
- Double_t Qyamean2 = fMeanQ[iC][1][1];
- Double_t Qyarms2 = fWidthQ[iC][1][1];
-
- Double_t Qxcmean2 = fMeanQ[iC][0][0];
- Double_t Qxcrms2 = fWidthQ[iC][0][0];
- Double_t Qycmean2 = fMeanQ[iC][0][1];
- Double_t Qycrms2 = fWidthQ[iC][0][1];
-
- Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
- Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
- Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
- Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
-
- fRPV0A = TMath::ATan2(QyaCor2, QxaCor2)/2.;
- fRPV0C = TMath::ATan2(QycCor2, QxcCor2)/2.;
-
- while(fRPV0A<0)fRPV0A+=TMath::Pi() ;
- while(fRPV0A>TMath::Pi())fRPV0A-=TMath::Pi() ;
- while(fRPV0C<0)fRPV0C+=TMath::Pi() ;
- while(fRPV0C>TMath::Pi())fRPV0C-=TMath::Pi() ;
+ FillHistogram("phiRPV0A" ,fRPV0A,fCentralityV0M);
+ FillHistogram("phiRPV0C" ,fRPV0C,fCentralityV0M);
+ FillHistogram("phiRPV0AC",fRPV0A,fRPV0C,fCentralityV0M) ;
// Flattening
fRPV0A=ApplyFlatteningV0A(fRPV0A,fCentralityV0M) ;
- while(fRPV0A<0)fRPV0A+=TMath::Pi() ;
- while(fRPV0A>TMath::Pi())fRPV0A-=TMath::Pi() ;
+ while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
+ while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
+
fRPV0C=ApplyFlatteningV0C(fRPV0C,fCentralityV0M) ;
- while(fRPV0C<0)fRPV0C+=TMath::Pi() ;
- while(fRPV0C>TMath::Pi())fRPV0C-=TMath::Pi() ;
+ while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
+ while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
if( fDebug >= 2 )
- AliInfo(Form("V0 Reaction Plane is: A side: %f, C side: %f", fRPV0A, fRPV0C));
-
- FillHistogram("phiRPV0A",fRPV0A,fCentralityV0M);
- FillHistogram("phiRPV0C",fRPV0C,fCentralityV0M);
+ AliInfo(Form("V0 Reaction Plane after flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
FillHistogram("phiRPV0Aflat",fRPV0A,fCentralityV0M) ;
- FillHistogram("phiRPV0AC",fRPV0A,fRPV0C,fCentralityV0M) ;
FillHistogram("cos2V0AC",TMath::Cos(2.*(fRPV0A-fRPV0C)),fCentralityV0M) ;
if(fHaveTPCRP){
FillHistogram("phiRPV0ATPC",fRP,fRPV0A,fCentralityV0M) ;
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff