fMinCentr(20),
fMaxCentr(80),
fEtaGap(kFALSE),
- fSubEvents(0)
+ fSubEvents(2)
{
// Default constructor
}
fMinCentr(20),
fMaxCentr(80),
fEtaGap(kFALSE),
- fSubEvents(0)
+ fSubEvents(2)
{
// standard constructor
Int_t pdg=421;
AliEventplane *pl=aod->GetEventplane();
if(!pl){
AliError("AliAnalysisTaskSEHFv2::UserExec:no eventplane! v2 analysis without eventplane not possible!\n");
+ fhEventsInfo->Fill(6);
return;
}
Double_t rpangleVZERO=0;
Double_t planereso=0;
Double_t deltaPsi=0;
+ Double_t rpangleeventC=0;
Double_t rpangleeventB=0;
Double_t rpangleeventA=0;
delete g;g=0x0;
if(fUseAfterBurner)fAfterBurner->SetEventPlane((Double_t)eventplane);
}else{
- if(fEventPlaneMeth!=kTPC){
-
- //VZERO EP and resolution
- rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0",aod,fV0EPorder));
- if(fEventPlaneMeth>kTPCVZERO){
- //Using V0A/C for VZERO resolution
- if(fEventPlaneMeth==kVZEROA){
- rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0A",aod,fV0EPorder));
- if(fSubEvents!=kSingleV0Side)rpangleeventB=GetPhi0Pi(pl->GetEventplane("V0C",aod,fV0EPorder));
- }
- else if(fEventPlaneMeth==kVZEROC){
- rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0C",aod,fV0EPorder));
- if(fSubEvents!=kSingleV0Side)rpangleeventB=GetPhi0Pi(pl->GetEventplane("V0A",aod,fV0EPorder));
- }
- // deltaPsi =rpangleeventA-rpangleeventB;
- eventplane=rpangleVZERO;
- }
- }
-
- // TPC event plane
+ // TPC event plane
rpangleTPC = pl->GetEventplane("Q");
if(rpangleTPC<0){
fhEventsInfo->Fill(6);
return;
}
rpangleeventA = rpangleTPC;
- if(fSubEvents>kFullTPC||fEtaGap||fEventPlaneMeth==kVZERO){
+ if(fSubEvents==2||fEventPlaneMeth==kVZERO){
qsub1 = pl->GetQsub1();
qsub2 = pl->GetQsub2();
if(!qsub1 || !qsub2){
fhEventsInfo->Fill(6);
return;
}
- if(fSubEvents==kPosTPC||fSubEvents==kSingleV0Side||fEventPlaneMeth==kVZERO){
- rpangleeventA = qsub1->Phi()/2.;
- if(fSubEvents==kSingleV0Side||fEventPlaneMeth==kVZERO)rpangleeventB = qsub2->Phi()/2.;
+ rpangleeventA = qsub1->Phi()/2.;
+ rpangleeventB = qsub2->Phi()/2.;
+ }
+ if(fEventPlaneMeth!=kTPC){
+ //VZERO EP and resolution
+ rpangleVZERO=GetPhi0Pi(pl->GetEventplane("V0",aod,fV0EPorder));
+ rpangleeventC=rpangleVZERO;
+ if(fEventPlaneMeth==kVZEROA||fEventPlaneMeth==kVZEROC||(fEventPlaneMeth==kTPCVZERO&&fSubEvents==3)){
+ rpangleeventB=GetPhi0Pi(pl->GetEventplane("V0A",aod,fV0EPorder));
+ rpangleeventC=GetPhi0Pi(pl->GetEventplane("V0C",aod,fV0EPorder));
+ if(fEventPlaneMeth==kVZEROA)rpangleVZERO=rpangleeventB;
+ else if(fEventPlaneMeth==kVZEROC)rpangleVZERO=rpangleeventC;
}
- if(fSubEvents==kNegTPC)rpangleeventA = qsub2->Phi()/2.;
}
-
- if(fEventPlaneMeth<=kTPCVZERO){
+
+
+ if(fEventPlaneMeth>kTPCVZERO)eventplane=rpangleVZERO;
+ else{
q = pl->GetQVector();
- deltaPsi = pl->GetQsubRes();
- }else{
- deltaPsi=eventplane-rpangleeventA;
+ eventplane=rpangleTPC;
}
+ deltaPsi =rpangleeventA-rpangleeventB;
}
-
+
if(TMath::Abs(deltaPsi)>TMath::Pi()/2.){
if(deltaPsi>0.) deltaPsi-=TMath::Pi();
else deltaPsi +=TMath::Pi();
((TH1F*)fOutput->FindObject(Form("hEvPlaneA%s",centrbinname.Data())))->Fill(rpangleeventA); //Angle of first subevent
((TH1F*)fOutput->FindObject(Form("hEvPlaneB%s",centrbinname.Data())))->Fill(rpangleeventB); //Angle of second subevent
}
- if(fEventPlaneMeth>kTPCVZERO){
- Double_t deltaSub=rpangleeventA-rpangleeventB;
+ if(fEventPlaneMeth>kTPCVZERO||fSubEvents==3){
+ Double_t deltaSub=rpangleeventA-rpangleeventC;
if(TMath::Abs(deltaSub)>TMath::Pi()/2.){// difference of subevents reaction plane angle cannot be bigger than phi/2
if(deltaSub>0.) deltaSub-=TMath::Pi();
else deltaSub +=TMath::Pi();
}
((TH1F*)fOutput->FindObject(Form("hEvPlaneReso2%s",centrbinname.Data())))->Fill(TMath::Cos(2.*deltaSub)); //RP resolution
- deltaSub =eventplane-rpangleeventB;
+ deltaSub =rpangleeventB-rpangleeventC;
if(TMath::Abs(deltaSub)>TMath::Pi()/2.){// difference of subevents reaction plane angle cannot be bigger than phi/2
if(deltaSub>0.) deltaSub-=TMath::Pi();
else deltaSub +=TMath::Pi();
}
fEventPlaneMeth=method;
}
-
+//________________________________________________________________________
+void AliAnalysisTaskSEHFv2::SetNTPCSubEvents(Int_t nsub){
+ if(nsub>3||nsub<2){
+ AliWarning("Only 2 or 3 subevents implemented. Setting 2 subevents\n");
+ nsub=2;
+ }
+ if(fEventPlaneMeth==kTPC&&nsub==3){
+ AliWarning("V0 must be enabled to run 3 TPC subevents. Enabling...\n");
+ fEventPlaneMeth=kTPCVZERO;
+ }
+ fSubEvents=nsub;
+}
//________________________________________________________________________
Float_t AliAnalysisTaskSEHFv2::GetPhi0Pi(Float_t phi){
// Sets the phi angle in the range 0-pi
enum DecChannel{kDplustoKpipi,kD0toKpi,kDstartoKpipi}; //more particles can be added
enum EventPlaneMeth{kTPC,kTPCVZERO,kVZERO,kVZEROA,kVZEROC}; //Event plane to be calculated in the task
- enum SubEvents{kFullTPC,kPosTPC,kNegTPC,kSingleV0Side}; //Sub-events for V0 EP
-
+ // enum SubEvents{kFullTPC,kPosTPC,kNegTPC,kSingleV0Side}; //Sub-events for V0 EP
+
AliAnalysisTaskSEHFv2();
AliAnalysisTaskSEHFv2(const char *name, AliRDHFCuts *rdCuts, Int_t decaychannel);
void SetUseAfterBurner(Bool_t ab){fUseAfterBurner=ab;}
void SetAfterBurner(AliHFAfterBurner *ab){fAfterBurner=ab;}
void SetEtaGapFeatureForEventplaneFromTracks (Bool_t etaGap) {fEtaGap = etaGap;}
+
+ void SetNTPCSubEvents(Int_t nsub);
+ void Set2TPCEPSubEvents(){SetNTPCSubEvents(2);}
+ void Set3TPCEPSubEvents(){SetNTPCSubEvents(3);}
void SetEventPlaneMethod(Int_t epmethod);
void SetTPCEPOnly(){SetEventPlaneMethod(kTPC);}
void SetVZEROEP(){SetEventPlaneMethod(kVZERO);}
void SetVZEROCEP(){SetEventPlaneMethod(kVZEROC);}
void SetTPCEP(){SetEventPlaneMethod(kTPCVZERO);}
void SetEventPlanesCompatibility(Float_t comp) {fEventPlanesComp=comp;}
- void SetSubEvents(Int_t subev){if(subev>=kFullTPC&&subev<=kSingleV0Side)fSubEvents=subev;}
+ //void SetSubEvents(Int_t subev){if(subev>=kFullTPC&&subev<=kSingleV0Side)fSubEvents=subev;}
Int_t GetEventPlaneMethod()const {return fEventPlaneMeth;}
+ Int_t GetNTPCSubEvents()const {return fSubEvents;}
Float_t GetEventPlanesCompatibility()const {return fEventPlanesComp;}
Float_t GetUpperMassLimit()const {return fUpmasslimit;}
Float_t GetLowerMassLimit()const {return fLowmasslimit;}
Int_t fMinCentr; //minimum centrality
Int_t fMaxCentr; //maximum centrality
Bool_t fEtaGap; // Eta gap feature for Eventplane from tracks; be careful that you do the correct settings in AddTaskEventPlane.C !!!!
- Int_t fSubEvents; //Sub-events definition for V0 EP
+ Int_t fSubEvents; //Sub-events definition for TPC EP
ClassDef(AliAnalysisTaskSEHFv2,2); // AliAnalysisTaskSE for the HF v2 analysis
};
#endif
+/* $Id$ */
+
//methods for the analysis of AliAnalysisTaskSEHFv2 output
//Authors: Chiara Bianchin cbianchi@pd.infn.it
// Giacomo Ortona ortona@to.infn.it
// Francesco Prino prino@to.infn.it
+//evPlane flag: -1=V0C,0=V0,1=V0A,2=TPC2subevs,3=TPC3subevs
//global variables to be set
Bool_t gnopng=kTRUE; //don't save in png format (only root and eps)
const Int_t nptbinsnew=3;
Int_t minPtBin[nptbinsnew]={-1,-1,-1};
Int_t maxPtBin[nptbinsnew]={-1,-1,-1};
Double_t mass;
+
//methods
-//Bool_t ReadFile(TList* &list,TH1F* &hstat,AliRDHFCuts* &cutobj,TString listname,TString partname,TString path="./",TString filename="AnalysisResults.root");
TList *LoadMassHistos(TList *inputlist,Int_t minCent,Int_t maxCent,Bool_t inoutanis);
Int_t FindPtBin(Int_t nbins, Float_t* array,Float_t value);
void FillSignalGraph(TList *histlist,TGraphAsymmErrors **gSignal,TGraphAsymmErrors **gSignalfs,Bool_t inoutanis);
void DmesonsFlowAnalysis(Bool_t inoutanis=kTRUE,Int_t minC=30,Int_t maxC=50,TString partname="Dplus",Int_t evPlane=2);
void DrawEventPlane(Int_t mincentr=0,Int_t maxcentr=0,TString filename="AnalysisResults.root",TString dirname="PWG3_D2H_HFv2",TString listname="coutputv2",Int_t evPlane=2);
-Double_t DrawEventPlane(TList *list,Int_t mincentr=0,Int_t maxcentr=0,Int_t evPlane=2);
+Double_t DrawEventPlane(TList *list,Int_t mincentr=0,Int_t maxcentr=0,Int_t evPlane=2,Double_t &error);
void DrawEventPlaneAllCentralities(TList *list,Int_t startCentr=20,Int_t endCentr=80,Int_t evPlane=2);
-//Aggiungere <pt> method
//methods implementation
//________________________________________________________________________________
TString dirname="PWG3_D2H_HFv2";
TString listname="coutputv2";
+ if(evPlane>3||evPlane<-1){
+ printf("Wrong EP flag %d \n",evPlane);
+ return;
+ }
+
AliRDHFCuts *cutsobj=0x0;
//Load input data from AliAnalysisTaskSEHFv2
TFile *f = TFile::Open(filename.Data());
Int_t nphi=nphibins;
if(inoutanis)nphi=2;
- Int_t nSubRes=1;
- if(evPlane<2)nSubRes=3;//3 sub-events method for VZERO EP
//EP resolution
TH1F* hevplresos=(TH1F*)list->FindObject(Form("hEvPlaneResocentr%d_%d",minC,minC+5));
TH1F* hevplresos2=0;
TH1F* hevplresos3=0;
- if(evPlane<2){
+ if(evPlane!=2){
hevplresos2=(TH1F*)list->FindObject(Form("hEvPlaneReso2centr%d_%d",minC,minC+5));
hevplresos3=(TH1F*)list->FindObject(Form("hEvPlaneReso3centr%d_%d",minC,minC+5));
}
Double_t resol=0;
for(Int_t icent=minC+5;icent<maxC;icent=icent+5){
hevplresos->Add((TH1F*)list->FindObject(Form("hEvPlaneResocentr%d_%d",icent,icent+5)));
- if(evPlane<2){
+ if(evPlane!=2){
hevplresos2->Add((TH1F*)list->FindObject(Form("hEvPlaneReso2centr%d_%d",icent,icent+5)));
hevplresos3->Add((TH1F*)list->FindObject(Form("hEvPlaneReso3centr%d_%d",icent,icent+5)));
}
}
- if(evPlane>=2){
+ if(evPlane==2){
resol=AliVertexingHFUtils::GetFullEvResol(hevplresos);
}else{
Double_t resolSub[3]={hevplresos->GetMean(),hevplresos2->GetMean(),hevplresos3->GetMean()};
- resol=TMath::Sqrt(resolSub[0]*resolSub[2]/resolSub[1]);
+ if(evPlane<=0)resol=TMath::Sqrt(resolSub[1]*resolSub[2]/resolSub[0]);
+ else if(evPlane==1) resol=TMath::Sqrt(resolSub[0]*resolSub[2]/resolSub[1]);
+ else if(evPlane==3) resol=TMath::Sqrt(resolSub[0]*resolSub[1]/resolSub[2]);
}
printf("average pt for pt bin \n");
if(!list){
printf("list %s not found in file, please check list name\n",listname.Data());return;
}
- DrawEventPlane(list,mincentr,maxcentr,evPlane);
+ DrawEventPlane(list,mincentr,maxcentr,evPlane,0);
}
//_______________________________________________________________________________________________________________________________
-Double_t DrawEventPlane(TList *list,Int_t mincentr,Int_t maxcentr,Int_t evPlane){
+Double_t DrawEventPlane(TList *list,Int_t mincentr,Int_t maxcentr,Int_t evPlane,Double_t &error){
//draw the histograms correlated to the event plane, returns event plane resolution
Double_t resolFull=0;
Int_t nSubRes=1;
- if(evPlane<2)nSubRes=3;//3 sub-events method for VZERO EP
+ if(evPlane!=2)nSubRes=3;//3 sub-events method for VZERO EP
TString namereso[3]={"Reso","Reso2","Reso3"};
TString suffixcentr=Form("centr%d_%d",mincentr,maxcentr);
TH2F* hevpls=(TH2F*)list->FindObject(Form("hEvPlanecentr%d_%d",mincentr,mincentr+5));
if(nSubRes>1){
hevplresos[1]->SetLineColor(2);
hevplresos[2]->SetLineColor(3);
- hevplresos[0]->SetTitle("cos2(#Psi_{V0A}-#Psi_{V0B})");
- hevplresos[1]->SetTitle("cos2(#Psi_{TPC}-#Psi_{V0A})");
- hevplresos[2]->SetTitle("cos2(#Psi_{TPC}-#Psi_{V0B})");
+ hevplresos[0]->SetTitle("cos2(#Psi_{TPC}-#Psi_{V0A})");
+ hevplresos[1]->SetTitle("cos2(#Psi_{TPC}-#Psi_{V0C})");
+ hevplresos[2]->SetTitle("cos2(#Psi_{V0A}-#Psi_{V0C})");
hevplresos[1]->Draw("SAME");
hevplresos[2]->Draw("SAME");
}
leg->SetLineColor(0);
leg->SetFillColor(0);
- if(nSubRes<=1) resolFull=AliVertexingHFUtils::GetFullEvResol(hevplresos[0]);
+ if(nSubRes<=1){
+ resolFull=AliVertexingHFUtils::GetFullEvResol(hevplresos[0]);
+ error = TMath::Abs(resolFull-AliVertexingHFUtils::GetSubEvResolLowLim(hevplresos[0]));
+ }
else{
Double_t resolSub[3];
- for(Int_t ires=0;ires<nSubRes;ires++)resolSub[ires]=hevplresos[ires]->GetMean();
- resolFull=TMath::Sqrt(resolSub[2]*resolSub[0]/resolSub[1]);
+ Double_t errors[3];
+ for(Int_t ires=0;ires<nSubRes;ires++){
+ resolSub[ires]=hevplresos[ires]->GetMean();
+ errors[ires]=hevplresos[ires]->GetMeanError();
+ }
+ Double_t lowlim[3];for(Int_t ie=0;ie<3;ie++)lowlim[ie]=TMath::Abs(resolSub[ie]-errors[ie]);
+ if(evPlane<=0){
+ resolFull=TMath::Sqrt(resolSub[1]*resolSub[2]/resolSub[0]);
+ error=resolFull-TMath::Sqrt(lowlim[2]*lowlim[1]/lowlim[0]);
+ }
+ else if(evPlane==1){
+ resolFull=TMath::Sqrt(resolSub[0]*resolSub[2]/resolSub[1]);
+ error=resolFull-TMath::Sqrt(lowlim[2]*lowlim[0]/lowlim[1]);
+ }
+ else if(evPlane==3){
+ resolFull=TMath::Sqrt(resolSub[0]*resolSub[1]/resolSub[2]);
+ error=resolFull-TMath::Sqrt(lowlim[0]*lowlim[1]/lowlim[2]);
+ }
}
TPaveText* pvreso=new TPaveText(0.1,0.1,0.6,0.2,"NDC");
//_______________________________________________________________________________________________________________________________
void DrawEventPlaneAllCentralities(TList *list,Int_t startCentr,Int_t endCentr,Int_t evPlane){
- TGraph* gresovscentr=new TGraph(0);
+ TGraphErrors* gresovscentr=new TGraphErrors(0);
gresovscentr->SetName("gresovscentr");
gresovscentr->SetTitle(Form("Resolution vs Centrality;centrality (%s);Resolution","%"));
Int_t iPoint=0;
Int_t step=5;//5% centrality step
// Int_t nCC=(endCentr-startCentr)/step;
-
+ Double_t errory=0;
for(Int_t i=startCentr;i<endCentr;i=i+step){
- Double_t resolFull=DrawEventPlane(list,i,i+step,evPlane);
+ Double_t resolFull=DrawEventPlane(list,i,i+step,evPlane,errory);
gresovscentr->SetPoint(iPoint,i+(Float_t)step/2.,resolFull);
+ gresovscentr->SetPointError(iPoint,5./2.,errory);
iPoint++;
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff