#include <TROOT.h>
#include <TRandom.h>
+#include <TString.h>
#include <TSystem.h>
#include <TInterpreter.h>
#include <TChain.h>
#include <TList.h>
#include <TLorentzVector.h>
#include <TClonesArray.h>
-#include "TDatabasePDG.h"
+#include "TDatabasePDG.h"
#include "AliAnalysisTaskJetServices.h"
#include "AliAnalysisDataContainer.h"
#include "AliGenCocktailEventHeader.h"
#include "AliInputEventHandler.h"
#include "AliPhysicsSelection.h"
-
+#include "AliTriggerAnalysis.h"
#include "AliAnalysisHelperJetTasks.h"
AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(): AliAnalysisTaskSE(),
fUseAODInput(kFALSE),
+ fUsePhysicsSelection(kFALSE),
+ fMC(kFALSE),
+ fSelectionInfoESD(0),
+ fEventCutInfoESD(0),
fAvgTrials(1),
- fZVtxCut(8.),
- fRealData(kFALSE),
- fPhysicsSelection(0x0),
+ fVtxXMean(0),
+ fVtxYMean(0),
+ fVtxZMean(0),
+ fVtxRCut(1.),
+ fVtxZCut(8.),
+ fPtMinCosmic(5.),
+ fRIsolMinCosmic(3.),
+ fMaxCosmicAngle(0.01),
fh1Xsec(0x0),
fh1Trials(0x0),
fh1PtHard(0x0),
fh1PtHardTrials(0x0),
+ fh1SelectionInfoESD(0x0),
+ fh1EventCutInfoESD(0),
fh2TriggerCount(0x0),
fh2ESDTriggerCount(0x0),
fh2TriggerVtx(0x0),
fh2ESDTriggerVtx(0x0),
fh2ESDTriggerRun(0x0),
fh2VtxXY(0x0),
+ fh1NCosmicsPerEvent(0x0),
+ fTriggerAnalysis(0x0),
fHistList(0x0)
{
fRunRange[0] = fRunRange[1] = 0;
AliAnalysisTaskJetServices::AliAnalysisTaskJetServices(const char* name):
AliAnalysisTaskSE(name),
fUseAODInput(kFALSE),
+ fUsePhysicsSelection(kFALSE),
+ fMC(kFALSE),
+ fSelectionInfoESD(0),
+ fEventCutInfoESD(0),
fAvgTrials(1),
- fZVtxCut(8.),
- fRealData(kFALSE),
- fPhysicsSelection(0x0),
+ fVtxXMean(0),
+ fVtxYMean(0),
+ fVtxZMean(0),
+ fVtxRCut(1.),
+ fVtxZCut(8.),
+ fPtMinCosmic(5.),
+ fRIsolMinCosmic(3.),
+ fMaxCosmicAngle(0.01),
fh1Xsec(0x0),
fh1Trials(0x0),
fh1PtHard(0x0),
fh1PtHardTrials(0x0),
+ fh1SelectionInfoESD(0x0),
+ fh1EventCutInfoESD(0),
fh2TriggerCount(0x0),
fh2ESDTriggerCount(0x0),
fh2TriggerVtx(0x0),
fh2ESDTriggerVtx(0x0),
fh2ESDTriggerRun(0x0),
fh2VtxXY(0x0),
+ fh1NCosmicsPerEvent(0x0),
+ fTriggerAnalysis(0x0),
fHistList(0x0)
{
fRunRange[0] = fRunRange[1] = 0;
Bool_t oldStatus = TH1::AddDirectoryStatus();
TH1::AddDirectory(kFALSE);
-
- if(!fPhysicsSelection)
- fPhysicsSelection = new AliPhysicsSelection();
- fPhysicsSelection->SetName("AliPhysicsSelection_outputlist"); // to prevent conflict with object that is automatically streamed back
- //AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kDebug);
- fHistList->Add(fPhysicsSelection);
-
fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
fHistList->Add(fh1Xsec);
fHistList->Add(fh1PtHard);
fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",nBinPt,binLimitsPt);
fHistList->Add(fh1PtHardTrials);
+ fh1SelectionInfoESD = new TH1F("fh1SelectionInfoESD","Bit Masks that satisfy the selection info",
+ AliAnalysisHelperJetTasks::kTotalSelections,0.5,AliAnalysisHelperJetTasks::kTotalSelections+0.5);
+ fHistList->Add(fh1SelectionInfoESD);
+
+ fh1EventCutInfoESD = new TH1F("fh1EventCutInfoESD","Bit Masks that for the events after each step of cuts",
+ kTotalEventCuts,0.5,kTotalEventCuts+0.5);
+ fHistList->Add(fh1EventCutInfoESD);
// 3 decisions, 0 trigger X, X + SPD vertex, X + SPD vertex in range
// 3 triggers BB BE/EB EE
if(hn)hn->Sumw2();
}
+ fh1NCosmicsPerEvent = new TH1F("fh1NCosmicsPerEvent","Number of cosmic candidates per event",10,0.,10.);
+ fHistList->Add(fh1NCosmicsPerEvent),
+
TH1::AddDirectory(oldStatus);
}
//
// Initialization
//
-
- Printf(">>> AnalysisTaskJetServices::Init() debug level %d\n",fDebug);
if (fDebug > 1) printf("AnalysisTaskJetServices::Init() \n");
}
AliAODEvent *aod = 0;
AliESDEvent *esd = 0;
-
+
AliAnalysisHelperJetTasks::Selected(kTRUE,kFALSE); // set slection to false
+ fSelectionInfoESD = 0; // reset
+ fEventCutInfoESD = 0; // reset
+ AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD); // set slection to false
+
if(fUseAODInput){
aod = dynamic_cast<AliAODEvent*>(InputEvent());
esd = dynamic_cast<AliESDEvent*>(InputEvent());
}
+ fSelectionInfoESD |= kNoEventCut;
+ fEventCutInfoESD |= kNoEventCut;
+
+ Bool_t esdVtxValid = false;
+ Bool_t esdVtxIn = false;
+ Bool_t aodVtxValid = false;
+ Bool_t aodVtxIn = false;
+
if(esd){
+ // trigger analyisis
+ if(!fTriggerAnalysis){
+ fTriggerAnalysis = new AliTriggerAnalysis;
+ fTriggerAnalysis->SetAnalyzeMC(fMC);
+ fTriggerAnalysis->SetSPDGFOThreshhold(1);
+ }
+ // fTriggerAnalysis->FillTriggerClasses(esd);
+ Bool_t v0A = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0A);
+ Bool_t v0C = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0C);
+ Bool_t v0ABG = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0ABG);
+ Bool_t v0CBG = fTriggerAnalysis->IsOfflineTriggerFired(esd, AliTriggerAnalysis::kV0CBG);
+ Bool_t spdFO = fTriggerAnalysis->SPDFiredChips(esd, 0);;
+ if(v0A)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kV0A;
+ if(v0C)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kV0C;
+ if(!(v0ABG||v0CBG))fSelectionInfoESD |= AliAnalysisHelperJetTasks::kNoV0BG;
+ if(spdFO)fSelectionInfoESD |= AliAnalysisHelperJetTasks::kSPDFO;
+
Float_t run = (Float_t)esd->GetRunNumber();
const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
+ esdVtxValid =
+ esdVtxIn = IsVertexIn(vtxESD);
Float_t zvtx = -999;
Float_t xvtx = -999;
Float_t yvtx = -999;
- if(vtxESD->GetNContributors()>0){
+ if(esdVtxValid){
zvtx = vtxESD->GetZ();
yvtx = vtxESD->GetY();
xvtx = vtxESD->GetX();
}
+ // CKB this can be cleaned up a bit...
+
Int_t iTrig = -1;
if(esd->GetFiredTriggerClasses().Contains("CINT1B")
||esd->GetFiredTriggerClasses().Contains("CSMBB")
||esd->GetFiredTriggerClasses().Contains("MB1")
||esd->GetFiredTriggerClasses().Contains("CINT6B")){
iTrig = 0;
+ fSelectionInfoESD |= AliAnalysisHelperJetTasks::kBunchBunch;
}
else if(esd->GetFiredTriggerClasses().Contains("CINT1A")
||esd->GetFiredTriggerClasses().Contains("CSMBA")
||esd->GetFiredTriggerClasses().Contains("CINT6C")){
// empty bunch or bunch empty
iTrig = 1;
+ fSelectionInfoESD |= AliAnalysisHelperJetTasks::kBunchEmpty;
}
- if(esd->GetFiredTriggerClasses().Contains("CINT1-E")
+ else if(esd->GetFiredTriggerClasses().Contains("CINT1-E")
||esd->GetFiredTriggerClasses().Contains("CINT6-E")){
iTrig = 2;
+ fSelectionInfoESD |= AliAnalysisHelperJetTasks::kEmptyEmpty;
}
if(iTrig>=0){
iTrig *= 3;
fh2ESDTriggerRun->Fill(run,iTrig+1);
- if(vtxESD->GetNContributors()>0){
+ if(vtxESD->GetNContributors()>2){
fh2ESDTriggerRun->Fill(run,iTrig+2);
fh2VtxXY->Fill(xvtx,yvtx);
}
- if(TMath::Abs(zvtx)<fZVtxCut&&TMath::Abs(xvtx)<0.5&&TMath::Abs(yvtx)<0.5)fh2ESDTriggerRun->Fill(run,iTrig+3);
+ xvtx -= fVtxXMean;
+ yvtx -= fVtxYMean;
+ zvtx -= fVtxZMean;
+ Float_t r2 = xvtx *xvtx + yvtx *yvtx;
+ if(TMath::Abs(zvtx)<fVtxZCut&&r2<(fVtxRCut*fVtxRCut))fh2ESDTriggerRun->Fill(run,iTrig+3);
}
else{
fh2ESDTriggerRun->Fill(run,0);
}
+ // BKC
+ }
+
+
+ // Apply additional constraints
+ Bool_t esdEventSelected = IsEventSelected(esd);
+ Bool_t esdEventPileUp = IsEventPileUp(esd);
+ Bool_t esdEventCosmic = IsEventCosmic(esd);
+
+ Bool_t aodEventSelected = IsEventSelected(aod);
+ Bool_t physicsSelection = ((fInputHandler->IsEventSelected())&AliVEvent::kMB);
+ fEventCutInfoESD |= kPhysicsSelectionCut; // other alreay set via IsEventSelected
+ fh1EventCutInfoESD->Fill(fEventCutInfoESD);
+ if(esdEventSelected) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kVertexIn;
+ if(esdEventPileUp) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsPileUp;
+ if(esdEventCosmic) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kIsCosmic;
+ if(physicsSelection) fSelectionInfoESD |= AliAnalysisHelperJetTasks::kPhysicsSelection;
+
+
+ // here we have all selection information, fill histogram
+ for(unsigned int i = 1;i<(UInt_t)fh1SelectionInfoESD->GetNbinsX();i++){
+ if((i&fSelectionInfoESD)==i)fh1SelectionInfoESD->Fill(i);
}
-
- // loop over all possible trigger and
- for(int it = AliAnalysisHelperJetTasks::kAcceptAll;it < AliAnalysisHelperJetTasks::kTrigger;it++){
- Bool_t esdTrig = kFALSE;
- if(esd){
+ AliAnalysisHelperJetTasks::SelectInfo(kTRUE,fSelectionInfoESD);
+
+ if(esd&&aod&&fDebug){
+ if(esdEventSelected&&!aodEventSelected){
+ Printf("%s:%d Different Selection for ESD and AOD",(char*)__FILE__,__LINE__);
+ const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
+ const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
+ Printf("ESD Vtx %s %s %d",vtxESD->GetName(),vtxESD->GetTitle(),vtxESD->GetNContributors());
+ vtxESD->Print();
+ Printf("AOD Vtx %s %s %d",vtxAOD->GetName(),vtxAOD->GetTitle(),vtxAOD->GetNContributors());
+ vtxAOD->Print();
+ }
+ }
+
+ // loop over all possible triggers for esd
+
+ if(esd){
+ const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
+ // Printf(">> ESDvtx %s %s",vtxESD->GetName(),vtxESD->GetTitle());vtxESD->Print();
+ TString vtxName(vtxESD->GetName());
+ for(int it = AliAnalysisHelperJetTasks::kAcceptAll;it < AliAnalysisHelperJetTasks::kTrigger;it++){
+ Bool_t esdTrig = kFALSE;
esdTrig = AliAnalysisHelperJetTasks::IsTriggerFired(esd,(AliAnalysisHelperJetTasks::Trigger)it);
if(esdTrig)fh2ESDTriggerCount->Fill(it,kAllTriggered);
+ Bool_t cand = physicsSelection;
+ if(cand){
+ fh2ESDTriggerCount->Fill(it,kSelectedALICE);
+ }
+ if(!fUsePhysicsSelection)cand = AliAnalysisHelperJetTasks::IsTriggerFired(esd,AliAnalysisHelperJetTasks::kMB1);
+ if(vtxESD->GetNContributors()>2&&!vtxName.Contains("TPCVertex")){
+ if(esdTrig)fh2ESDTriggerCount->Fill(it,kTriggeredVertex);
+ Float_t zvtx = vtxESD->GetZ();
+ if(esdEventSelected&&esdTrig){
+ fh2ESDTriggerCount->Fill(it,kTriggeredVertexIn);
+ fh2ESDTriggerVtx->Fill(it,zvtx);
+ }
+ if(cand)fh2ESDTriggerCount->Fill(it,kSelectedALICEVertexValid);
+ }
+ if(cand&&esdEventSelected){
+ fh2ESDTriggerCount->Fill(it,kSelectedALICEVertexIn);
+ fh2ESDTriggerCount->Fill(it,kSelected);
+ AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
+ }
}
- Bool_t aodTrig = kFALSE;
- if(aod){
+ }
+
+ if(aod){
+ const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
+ aodVtxValid = IsVertexValid(vtxAOD);
+ aodVtxIn = IsVertexIn(vtxAOD);
+
+ for(int it = AliAnalysisHelperJetTasks::kAcceptAll;it < AliAnalysisHelperJetTasks::kTrigger;it++){
+ Bool_t aodTrig = kFALSE;
aodTrig = AliAnalysisHelperJetTasks::IsTriggerFired(aod,(AliAnalysisHelperJetTasks::Trigger)it);
+ Bool_t cand = aod->GetHeader()->GetOfflineTrigger()&AliVEvent::kMB;
if(aodTrig)fh2TriggerCount->Fill(it,kAllTriggered);
- }
-
- // Check wether we have also an SPD vertex
-
- if(aod){
- const AliAODVertex *vtxAOD = aod->GetPrimaryVertex();
- // Printf(">> AODvtx %s %s",vtxAOD->GetName(),vtxAOD->GetTitle());vtxAOD->Print();
-
- if(vtxAOD->GetNContributors()>0){
- if(aodTrig)fh2TriggerCount->Fill(it,kTriggeredSPDVertex);
+ if(aodVtxValid){
+ if(aodTrig)fh2TriggerCount->Fill(it,kTriggeredVertex);
Float_t zvtx = vtxAOD->GetZ();
- Float_t yvtx = vtxAOD->GetY();
- Float_t xvtx = vtxAOD->GetX();
- fh2TriggerVtx->Fill(it,zvtx);
- if(TMath::Abs(zvtx)<fZVtxCut&&aodTrig&&TMath::Abs(xvtx)<0.5&&TMath::Abs(yvtx)<0.5){
+ if(cand&&aodEventSelected){
+ fh2TriggerCount->Fill(it,kSelected);
+ }
+ if(aodTrig&&aodEventSelected){
+ fh2TriggerVtx->Fill(it,zvtx);
fh2TriggerCount->Fill(it,kTriggeredVertexIn);
}
- }
- }
- if(esd){
- const AliESDVertex *vtxESD = esd->GetPrimaryVertex();
- // Printf(">> ESDvtx %s %s",vtxESD->GetName(),vtxESD->GetTitle());vtxESD->Print();
- Bool_t cand = fPhysicsSelection->IsCollisionCandidate(esd);
- if(cand) fh2ESDTriggerCount->Fill(it,kSelectedALICE);
- if(vtxESD->GetNContributors()>0){
- if(esdTrig)fh2ESDTriggerCount->Fill(it,kTriggeredSPDVertex);
- Float_t zvtx = vtxESD->GetZ();
- Float_t yvtx = vtxESD->GetY();
- Float_t xvtx = vtxESD->GetX();
- fh2ESDTriggerVtx->Fill(it,zvtx);
- if(TMath::Abs(zvtx)<fZVtxCut&&esdTrig&&TMath::Abs(xvtx)<0.5&&TMath::Abs(yvtx)<0.5){
- fh2ESDTriggerCount->Fill(it,kTriggeredVertexIn);
- // here we select based on ESD info...
- if(fRealData){
- if(cand){
- fh2ESDTriggerCount->Fill(it,kSelected);
- AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
- }
- }
- else{
- fh2ESDTriggerCount->Fill(it,kSelected);
- AliAnalysisHelperJetTasks::Selected(kTRUE,kTRUE);// select this event
- }
+ if(cand){
+ fh2TriggerCount->Fill(it,kSelectedALICEVertexValid);
+ if(aodEventSelected)fh2TriggerCount->Fill(it,kSelectedALICEVertexIn);
}
-
-
}
-
}
-
}
-
-
if (fDebug > 1)printf(" AliAnalysisTaskJetServices: Analysing event # %5d\n", (Int_t) fEntry);
PostData(1, fHistList);
}
+Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliESDEvent* esd){
+ if(!esd)return kFALSE;
+ const AliESDVertex *vtx = esd->GetPrimaryVertex();
+ return IsVertexIn(vtx); // vertex in calls vertex valid
+}
-void AliAnalysisTaskJetServices::Terminate(Option_t */*option*/)
-{
- // Terminate analysis
- //
- TDirectory* owd = gDirectory;
+Bool_t AliAnalysisTaskJetServices::IsEventSelected(const AliAODEvent* aod) const {
+ if(!aod)return kFALSE;
+ const AliAODVertex *vtx = aod->GetPrimaryVertex();
+ return IsVertexIn(vtx); // VertexIn calls VertexValid
+}
- if (fDebug > 1) printf("AnalysisJetServices: Terminate() \n");
+Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliESDVertex* vtx) {
+
+ // We can check the type of the vertex by its name and title
+ // if vertexer failed title is empty (default c'tor) and ncontributors is 0
+ // vtx name title
+ // Tracks PrimaryVertex VertexerTracksNoConstraint
+ // Tracks PrimaryVertex VertexerTracksWithConstraint
+ // TPC TPCVertex VertexerTracksNoConstraint
+ // TPC TPCVertex VertexerTracksWithConstraint
+ // SPD SPDVertex vertexer: 3D
+ // SPD SPDVertex vertexer: Z
+
+ Int_t nCont = vtx->GetNContributors();
+
+ if(nCont>=1){
+ fEventCutInfoESD |= kContributorsCut1;
+ if(nCont>=2){
+ fEventCutInfoESD |= kContributorsCut2;
+ if(nCont>=3){
+ fEventCutInfoESD |= kContributorsCut3;
+ }
+ }
+ }
+
+ if(nCont<3)return kFALSE;
- fHistList = dynamic_cast<TList*> (GetOutputData(1));
- if (!fHistList)
- Printf("ERROR: fHistList not available");
+ // do not want tpc only primary vertex
+ TString vtxName(vtx->GetName());
+ if(vtxName.Contains("TPCVertex")){
+ fEventCutInfoESD |= kVertexTPC;
+ return kFALSE;
+ }
+ if(vtxName.Contains("SPDVertex"))fEventCutInfoESD |= kVertexSPD;
+ if(vtxName.Contains("PrimaryVertex"))fEventCutInfoESD |= kVertexGlobal;
+ TString vtxTitle(vtx->GetTitle());
+ if(vtxTitle.Contains("vertexer: Z")){
+ if(vtx->GetDispersion()>0.02)return kFALSE;
+ }
+ fEventCutInfoESD |= kSPDDispersionCut;
+ return kTRUE;
+}
- AliAnalysisDataContainer *cont = GetOutputSlot(1)->GetContainer();
- TString filename = cont->GetFileName();
- TFile *f = NULL;
- // Check first if the file is already opened
- f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
- if (f) {
- // Cd to file
- f->cd();
- // Check for a folder request
- TString dir = cont->GetFolderName();
- if (!dir.IsNull()) {
- if (!f->GetDirectory(dir)) f->mkdir(dir);
- f->cd(dir);
- }
+
+Bool_t AliAnalysisTaskJetServices::IsVertexValid ( const AliAODVertex* vtx) const {
+
+ // We can check the type of the vertex by its name and title
+ // if vertexer failed title is empty (default c'tor) and ncontributors is 0
+ // vtx name title
+ // Tracks PrimaryVertex VertexerTracksNoConstraint
+ // TPC TPCVertex VertexerTracksNoConstraint
+ // SPD SPDVertex vertexer: 3D
+ // SPD SPDVertex vertexer: Z
+
+ if(vtx->GetNContributors()<3)return kFALSE;
+ // do not want tpc only primary vertex
+ TString vtxName(vtx->GetName());
+ if(vtxName.Contains("TPCVertex"))return kFALSE;
+
+ // no dispersion yet...
+ /*
+ TString vtxTitle(vtx->GetTitle());
+ if(vtxTitle.Contains("vertexer: Z")){
+ if(vtx->GetDispersion()>0.02)return kFALSE;
}
+ */
+ return kTRUE;
+}
+
+
+Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliESDVertex* vtx) {
+
+ if(!IsVertexValid(vtx))return kFALSE;
+
+ Float_t zvtx = vtx->GetZ();
+ Float_t yvtx = vtx->GetY();
+ Float_t xvtx = vtx->GetX();
+
+ xvtx -= fVtxXMean;
+ yvtx -= fVtxYMean;
+ zvtx -= fVtxZMean;
- if (fHistList)
- {
- fPhysicsSelection = dynamic_cast<AliPhysicsSelection*> (fHistList->FindObject("AliPhysicsSelection_outputlist"));
+
+
+ if(TMath::Abs(zvtx)>fVtxZCut){
+ return kFALSE;
}
+ fEventCutInfoESD |= kVertexZCut;
+ Float_t r2 = yvtx*yvtx+xvtx*xvtx;
+ if(r2>(fVtxRCut*fVtxRCut)){
+ return kFALSE;
+ }
+ fEventCutInfoESD |= kVertexRCut;
+ return kTRUE;
+}
+
+
+Bool_t AliAnalysisTaskJetServices::IsVertexIn (const AliAODVertex* vtx) const {
- if (fPhysicsSelection)
- {
- fPhysicsSelection->SaveHistograms("physics_selection");
- fPhysicsSelection->Print();
+ if(!IsVertexValid(vtx))return kFALSE;
+
+ Float_t zvtx = vtx->GetZ();
+ Float_t yvtx = vtx->GetY();
+ Float_t xvtx = vtx->GetX();
+
+ xvtx -= fVtxXMean;
+ yvtx -= fVtxYMean;
+ zvtx -= fVtxZMean;
+
+ Float_t r2 = yvtx*yvtx+xvtx*xvtx;
+
+ Bool_t vertexIn = TMath::Abs(zvtx)<fVtxZCut&&r2<(fVtxRCut*fVtxRCut);
+ return vertexIn;
+}
+
+Bool_t AliAnalysisTaskJetServices::IsEventPileUp(const AliESDEvent* esd) const{
+ if(!esd)return kFALSE;
+ return esd->IsPileupFromSPD();
+}
+
+Bool_t AliAnalysisTaskJetServices::IsEventCosmic(const AliESDEvent* esd) const {
+ if(!esd)return kFALSE;
+ // add track cuts for which we look for cosmics...
+
+ Bool_t isCosmic = kFALSE;
+ Int_t nTracks = esd->GetNumberOfTracks();
+ Int_t nCosmicCandidates = 0;
+
+ for (Int_t iTrack1 = 0; iTrack1 < nTracks; iTrack1++) {
+ AliESDtrack* track1 = (AliESDtrack*)esd->GetTrack(iTrack1);
+ if (!track1) continue;
+ UInt_t status1 = track1->GetStatus();
+ //If track is ITS stand alone track, skip the track
+ if (((status1 & AliESDtrack::kITSin) == 0 || (status1 & AliESDtrack::kTPCin))) continue;
+ if(track1->Pt()<fPtMinCosmic) continue;
+ //Start 2nd track loop to look for correlations
+ for (Int_t iTrack2 = iTrack1+1; iTrack2 < nTracks; iTrack2++) {
+ AliESDtrack* track2 = (AliESDtrack*)esd->GetTrack(iTrack2);
+ if(!track2) continue;
+ UInt_t status2 = track2->GetStatus();
+ //If track is ITS stand alone track, skip the track
+ if (((status2 & AliESDtrack::kITSin) == 0 || (status2 & AliESDtrack::kTPCin))) continue;
+ if(track2->Pt()<fPtMinCosmic) continue;
+ //Check if back-to-back
+ Double_t mom1[3],mom2[3];
+ track1->GetPxPyPz(mom1);
+ track2->GetPxPyPz(mom2);
+ TVector3 momv1(mom1[0],mom1[1],mom1[2]);
+ TVector3 momv2(mom2[0],mom2[1],mom2[2]);
+ Float_t theta = (float)(momv1.Phi()-momv2.Phi());
+ if(theta<-0.5*TMath::Pi()) theta+=2.*TMath::Pi();
+
+ Float_t deltaPhi = track1->Phi()-track2->Phi();
+ if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
+
+ Float_t rIsol = (float)(TMath::Sqrt( deltaPhi*deltaPhi+(track1->Eta()-track2->Eta())*(track1->Eta()-track2->Eta()) ));
+ if(rIsol<fRIsolMinCosmic) continue;
+
+ if(TMath::Abs(TMath::Pi()-theta)<fMaxCosmicAngle) {
+ nCosmicCandidates+=1;
+ isCosmic = kTRUE;
+ }
+
}
+ }
+
+ fh1NCosmicsPerEvent->Fill((float)nCosmicCandidates);
+
+ return isCosmic;
+}
- owd->cd();
+
+
+void AliAnalysisTaskJetServices::Terminate(Option_t */*option*/)
+{
+ // Terminate analysis
+ //
}