macro dir
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliEmcalPhysicsSelection.cxx
CommitLineData
b0a53615 1// $Id$
2//
980821ba 3// Emcal physics selection class.
b0a53615 4//
cd231d42 5// Author: C.Loizides
be94a07d 6
7#include "AliEmcalPhysicsSelection.h"
bb33f819 8#include "AliAODEvent.h"
be94a07d 9#include "AliESDEvent.h"
0abef324 10#include "AliLog.h"
be94a07d 11
12ClassImp(AliEmcalPhysicsSelection)
13
fb81cdb7 14AliEmcalPhysicsSelection::AliEmcalPhysicsSelection() :
15 AliPhysicsSelection(),
16 fMarkFastOnly(0),
17 fMarkLedEvent(0),
18 fSkipFastOnly(0),
19 fSkipLedEvent(0),
2c091b4a 20 fCellMinE(-1),
21 fClusMinE(-1),
22 fTrackMinPt(-1),
23 fTriggers(0),
1dc3ded9 24 fZvertex(-1),
25 fZvertexDiff(0),
26 fCentMin(-1),
27 fCentMax(-1),
6a7fb48f 28 fMinCellTrackScale(-1),
29 fMaxCellTrackScale(-1),
fb81cdb7 30 fIsFastOnly(0),
31 fIsLedEvent(0),
32 fIsGoodEvent(0),
33 fCellMaxE(0),
2c091b4a 34 fClusMaxE(0),
35 fTrackMaxPt(0)
fb81cdb7 36{
37 // Default constructor.
38}
39
be94a07d 40//__________________________________________________________________________________________________
41UInt_t AliEmcalPhysicsSelection::GetSelectionMask(const TObject* obj)
42{
bb33f819 43 // Calculate selection mask.
44
45 const AliVEvent *ev = dynamic_cast<const AliVEvent*>(obj);
09422fad 46 if (!ev)
bb33f819 47 return 0;
09422fad 48
49 AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
bb33f819 50
51 UInt_t res = 0;
52 const AliESDEvent *eev = dynamic_cast<const AliESDEvent*>(obj);
53 const AliAODEvent *aev = 0;
54 if (eev) {
09422fad 55 am->LoadBranch("AliESDHeader.");
56 am->LoadBranch("AliESDRun.");
57 TString title(eev->GetHeader()->GetTitle());
58 if (1&&(title.Length()>0)) {
59 res = eev->GetHeader()->GetUniqueID();
60 res &= 0x4FFFFFFF;
61 } else {
62 res = IsCollisionCandidate(eev);
63 }
bb33f819 64 } else {
65 aev = dynamic_cast<const AliAODEvent*>(obj);
66 res = aev->GetHeader()->GetOfflineTrigger();
67 }
be94a07d 68
1dc3ded9 69 // return 0, if 0 found
f4fa7802 70 if (res==0)
71 return 0;
72
1dc3ded9 73 // return 0, if ptrs are not set
74 if ((eev==0) && (aev==0))
75 return 0;
76
2c091b4a 77 if (fTriggers) { // only process given triggers
1dc3ded9 78 if ((res & fTriggers) == 0)
2c091b4a 79 return res;
80 }
81
fb81cdb7 82 fIsFastOnly = kFALSE;
83 fIsGoodEvent = kFALSE;
84 fIsLedEvent = kFALSE;
f4fa7802 85 fCellMaxE = -1;
86 fClusMaxE = -1;
87 fTrackMaxPt = -1;
fb81cdb7 88
1dc3ded9 89 if ((res & AliVEvent::kAnyINT) ||
90 (res & AliVEvent::kSemiCentral) ||
91 (res & AliVEvent::kCentral) ||
92 (res & AliVEvent::kEMC1) ||
93 (res & AliVEvent::kEMC7) ||
94 (res & AliVEvent::kEMCEJE) ||
bb33f819 95 (res & AliVEvent::kEMCEGA))
fb81cdb7 96 fIsGoodEvent = kTRUE;
09422fad 97 else {
98 return 0;
99 }
fb81cdb7 100
1dc3ded9 101 if (fZvertexDiff || (fZvertex>0)) {
102 Double_t vzPRI = +999;
103 Double_t vzSPD = -999;
104 const AliVVertex *pv = 0;
09422fad 105 if (eev) {
106 am->LoadBranch("PrimaryVertex.");
1dc3ded9 107 pv = eev->GetPrimaryVertex();
09422fad 108 } else {
1dc3ded9 109 pv = aev->GetPrimaryVertex();
09422fad 110 }
1dc3ded9 111 if (pv && pv->GetNContributors()>0) {
112 vzPRI = pv->GetZ();
113 }
114 const AliVVertex *sv = 0;
09422fad 115 if (eev) {
116 am->LoadBranch("SPDVertex.");
1dc3ded9 117 sv = eev->GetPrimaryVertexSPD();
09422fad 118 } else {
1dc3ded9 119 sv = aev->GetPrimaryVertexSPD();
09422fad 120 }
1dc3ded9 121 if (sv && sv->GetNContributors()>0) {
122 vzSPD = sv->GetZ();
123 }
124 Double_t dvertex = TMath::Abs(vzPRI-vzSPD);
1b91882c 125 // skip events with dvertex>1mm if requested
1dc3ded9 126 // https://indico.cern.ch/getFile.py/access?contribId=4&resId=0&materialId=slides&confId=189624
127 // also check on vertex z if requested
1b91882c 128 if (fZvertexDiff && (dvertex>0.1))
1dc3ded9 129 fIsGoodEvent = kFALSE;
ee139088 130 if ((fZvertex>0) && (TMath::Abs(vzPRI)>fZvertex))
1dc3ded9 131 fIsGoodEvent = kFALSE;
132 }
133
134 if ((fCentMin>-1) && (fCentMax>-1)) {
135 Double_t v0mcent = -1;
136 AliCentrality *centin = 0;
137 if (eev) {
6cf0dbe4 138 TTree *tree = am->GetTree();
139 if (tree->GetBranch("Centrality."))
140 am->LoadBranch("Centrality.");
1dc3ded9 141 centin = dynamic_cast<AliCentrality*>(eev->FindListObject("Centrality"));
142 } else {
143 centin = const_cast<AliAODEvent*>(aev)->GetCentrality();
144 }
145 if (centin)
146 v0mcent = centin->GetCentralityPercentileUnchecked("V0M");
77a933e7 147 if ((v0mcent<fCentMin) || (v0mcent>fCentMax))
1dc3ded9 148 fIsGoodEvent = kFALSE;
149 }
150
2c091b4a 151 AliVCaloCells *cells = ev->GetEMCALCells();
fb81cdb7 152 const Short_t nCells = cells->GetNumberOfCells();
153
154 // mark LHC11a fast only partition if requested
155 if (res & AliVEvent::kFastOnly) {
156 fIsFastOnly = kTRUE;
157 if (fMarkFastOnly||fSkipFastOnly) {
158 if (nCells>0) {
159 AliFatal(Form("Number of cells %d, even though EMCAL should not be in fast only partition.",nCells));
be94a07d 160 }
fb81cdb7 161 fIsGoodEvent = kFALSE;
be94a07d 162 }
163 }
164
2c091b4a 165 if (fCellMinE>0) {
166 if (eev)
167 am->LoadBranch("EMCALCells.");
168 for(Int_t iCell=0; iCell<nCells; ++iCell) {
169 Short_t cellId = cells->GetCellNumber(iCell);
170 Double_t cellE = cells->GetCellAmplitude(cellId);
171 if (cellE>fCellMaxE)
172 fCellMaxE = cellE;
173 }
fb81cdb7 174 }
2c091b4a 175
176 if (fClusMinE>0) {
177 if (eev)
178 am->LoadBranch("CaloClusters");
179
180 const Int_t nCaloClusters = ev->GetNumberOfCaloClusters();
181 for(Int_t iClus = 0; iClus<nCaloClusters; ++iClus) {
182 AliVCluster *cl = ev->GetCaloCluster(iClus);
183 if (!cl->IsEMCAL())
184 continue;
185 Double_t e = cl->E();
186 if (e>fClusMaxE)
187 fClusMaxE = e;
188 }
189 }
190
6a7fb48f 191 TClonesArray *trks = 0;
192 Int_t Ntracks = 0;
193
194 if (fTrackMinPt>0 || fMinCellTrackScale > 0 || fMaxCellTrackScale > 0) {
2c091b4a 195 if (eev) {
196 am->LoadBranch("PicoTracks");
197 trks = dynamic_cast<TClonesArray*>(eev->FindListObject("PicoTracks"));
198 if (!trks) {
199 am->LoadBranch("Tracks");
200 trks = dynamic_cast<TClonesArray*>(eev->FindListObject("Tracks"));
201 }
202 } else {
203 trks = dynamic_cast<TClonesArray*>(aev->FindListObject("tracks"));
204 }
6a7fb48f 205 if (trks)
206 Ntracks = trks->GetEntriesFast();
207 }
208
209 Int_t nAccTracks = 0;
210
211 if (fTrackMinPt>0 || fMinCellTrackScale > 0 || fMaxCellTrackScale > 0) {
212 for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
213 AliVTrack *track = static_cast<AliVTrack*>(trks->At(iTracks));
214 if (!track)
215 continue;
216 if (aev) { // a bit ugly since cuts are hard coded for now
217 AliAODTrack *aodtrack = static_cast<AliAODTrack*>(track);
218 if (!aodtrack->TestFilterBit(256) && !aodtrack->TestFilterBit(512))
219 continue;
f4fa7802 220 }
6a7fb48f 221 nAccTracks++;
222 Double_t pt = track->Pt();
223 if (pt>fTrackMaxPt)
224 fTrackMaxPt = pt;
2c091b4a 225 }
fb81cdb7 226 }
227
6a7fb48f 228 if (fMinCellTrackScale > 0 || fMaxCellTrackScale > 0) {
229 if (nCells < fMinCellTrackScale * nAccTracks || nCells > fMaxCellTrackScale * nAccTracks)
230 fIsGoodEvent = kFALSE;
231 }
232
fb81cdb7 233 // bad cell criterion for LHC11a from
234 // https://indico.cern.ch/materialDisplay.py?contribId=4&materialId=slides&confId=147067
235 const Int_t runN = ev->GetRunNumber();
236 if ((runN>=144871) && (runN<=146860)) {
2c091b4a 237
238 if (eev)
239 am->LoadBranch("EMCALCells.");
240
241 // count cells above threshold
242 Int_t nCellCount[12] = {0,0,0,0,0,0,0,0,0,0,0,0};
243 for(Int_t iCell=0; iCell<nCells; ++iCell) {
244 Short_t cellId = cells->GetCellNumber(iCell);
245 Double_t cellE = cells->GetCellAmplitude(cellId);
246 Int_t sm = cellId / (24*48);
247 if (cellE>0.1)
248 ++nCellCount[sm];
249 }
250
fb81cdb7 251 if (nCellCount[4] > 100)
252 fIsLedEvent = kTRUE;
253 else {
254 if ((runN>=146858) && (runN<=146860)) {
255 if ((res&AliVEvent::kMB) && (nCellCount[3]>=21))
256 fIsLedEvent = kTRUE;
257 else if ((res&AliVEvent::kEMC1) && (nCellCount[3]>=35))
258 fIsLedEvent = kTRUE;
259 }
260 }
261 if (fIsLedEvent) {
262 fIsGoodEvent = kFALSE;
263 }
264 }
265
77a933e7 266 if (fIsGoodEvent) {
267 if (fCellMaxE>fCellMinE)
268 res |= kEmcalHC;
269 if ((fClusMaxE>fClusMinE) || (fTrackMaxPt>fTrackMinPt))
270 res |= kEmcalHT;
fb81cdb7 271 res |= kEmcalOk;
77a933e7 272 }
fb81cdb7 273
274 if ((fSkipLedEvent && fIsLedEvent) ||
275 (fSkipFastOnly && fIsFastOnly))
276 res = 0;
277
be94a07d 278 return res;
279}