Added cut set for kstar analysis (fbellini)
[u/mrichter/AliRoot.git] / PWGLF / RESONANCES / AliRsnCutPIDNSigma.cxx
1 //
2 // Class AliRsnCutPIDNSigma
3 //
4 // General implementation of a single cut strategy, which can be:
5 // - a value contained in a given interval  [--> IsBetween()   ]
6 // - a value equal to a given reference     [--> MatchesValue()]
7 //
8 // In all cases, the reference value(s) is (are) given as data members
9 // and each kind of cut requires a given value type (Int, UInt, Double),
10 // but the cut check procedure is then automatized and chosen thanks to
11 // an enumeration of the implemented cut types.
12 // At the end, the user (or any other point which uses this object) has
13 // to use the method IsSelected() to check if this cut has been passed.
14 //
15 // authors: Martin Vala (martin.vala@cern.ch)
16 //          Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
17 //
18
19 #include "AliPIDResponse.h"
20 #include "AliESDpid.h"
21 #include "AliAODpidUtil.h"
22
23 #include "AliRsnCutPIDNSigma.h"
24
25 ClassImp(AliRsnCutPIDNSigma)
26
27 //_________________________________________________________________________________________________
28 AliRsnCutPIDNSigma::AliRsnCutPIDNSigma() :
29    AliRsnCut("cut", AliRsnTarget::kDaughter),
30    fSpecies(AliPID::kUnknown),
31    fDetector(kDetectors),
32    fRejectUnmatched(kTRUE),
33    fTrackNSigma(0.0),
34    fTrackMom(0.0),
35    fMyPID(0x0),
36    fRanges("AliRsnPIDRange", 0)
37 {
38 //
39 // Main constructor.
40 //
41 }
42
43 //_________________________________________________________________________________________________
44 AliRsnCutPIDNSigma::AliRsnCutPIDNSigma
45 (const char *name, AliPID::EParticleType species, EDetector det) :
46    AliRsnCut(name, AliRsnTarget::kDaughter),
47    fSpecies(species),
48    fDetector(det),
49    fRejectUnmatched(kTRUE),
50    fTrackNSigma(0.0),
51    fTrackMom(0.0),
52    fMyPID(0x0),
53    fRanges("AliRsnPIDRange", 0)
54 {
55 //
56 // Main constructor.
57 //
58 }
59
60 //_________________________________________________________________________________________________
61 AliRsnCutPIDNSigma::AliRsnCutPIDNSigma
62 (const AliRsnCutPIDNSigma &copy) :
63    AliRsnCut(copy),
64    fSpecies(copy.fSpecies),
65    fDetector(copy.fDetector),
66    fRejectUnmatched(copy.fRejectUnmatched),
67    fTrackNSigma(0.0),
68    fTrackMom(0.0),
69    fMyPID(copy.fMyPID),
70    fRanges(copy.fRanges)
71 {
72 //
73 // Copy constructor.
74 //
75 }
76
77 //_________________________________________________________________________________________________
78 AliRsnCutPIDNSigma &AliRsnCutPIDNSigma::operator=(const AliRsnCutPIDNSigma &copy)
79 {
80 //
81 // Assignment operator
82 //
83
84    AliRsnCut::operator=(copy);
85    if (this == &copy)
86       return *this;
87    fSpecies = copy.fSpecies;
88    fDetector = copy.fDetector;
89    fRejectUnmatched = copy.fRejectUnmatched;
90    fMyPID = copy.fMyPID;
91    fRanges = copy.fRanges;
92
93    return (*this);
94 }
95
96 //__________________________________________________________________________________________________
97 void AliRsnCutPIDNSigma::InitMyPID(Bool_t isMC, Bool_t isESD)
98 {
99 //
100 // Initialize manual PID object
101 //
102
103    if (isESD)
104       fMyPID = new AliESDpid(isMC);
105    else
106       fMyPID = new AliAODpidUtil(isMC);
107 }
108
109 //_________________________________________________________________________________________________
110 Bool_t AliRsnCutPIDNSigma::IsSelected(TObject *object)
111 {
112 //
113 // Cut checker.
114 // As usual, there are 'kFALSE' exit points whenever one of the conditions is not passed,
115 // and at the end, it returns kTRUE since it bypassed all possible exit points.
116 //
117
118    // coherence check
119    if (!TargetOK(object)) return kFALSE;
120
121    // check initialization of PID object
122    // if manual PID is used, use that, otherwise get from source event
123    AliPIDResponse *pid = 0x0;
124    if (fMyPID)
125       pid = fMyPID;
126    else
127       pid = fEvent->GetPIDResponse();
128    if (!pid) {
129       AliFatal("NULL PID response");
130       return kFALSE;
131    }
132
133    // convert input object into AliVTrack
134    // if this fails, the cut cannot be checked
135    AliVTrack *vtrack = fDaughter->Ref2Vtrack();
136    if (!vtrack) {
137       AliDebugClass(2, "Referenced daughter is not a track");
138       return kFALSE;
139    }
140
141    // check matching, if required
142    // a) if detector is not matched and matching is required, reject the track
143    // b) if detector is not matched and matching is not required, accept blindly the track
144    //    since missing the matching causes one not to be able to rely that detector
145    if (!MatchDetector(vtrack)) {
146       AliDebugClass(2, Form("Detector not matched. fRejectUnmatched = %s --> track is %s", (fRejectUnmatched ? "true" : "false"), (fRejectUnmatched ? "rejected" : "accepted")));
147       return (!fRejectUnmatched);
148    }
149
150    // get reference momentum
151    fTrackMom = (fDetector == kTPC) ? vtrack->GetTPCmomentum() : vtrack->P();
152
153    // get number of sigmas
154    switch (fDetector) {
155       case kITS:
156          fTrackNSigma = TMath::Abs(pid->NumberOfSigmasITS(vtrack, fSpecies));
157          break;
158       case kTPC:
159          fTrackNSigma = TMath::Abs(pid->NumberOfSigmasTPC(vtrack, fSpecies));
160          break;
161       case kTOF:
162          fTrackNSigma = TMath::Abs(pid->NumberOfSigmasTOF(vtrack, fSpecies));
163          break;
164       default:
165          AliError("Bad detector chosen. Rejecting track");
166          return kFALSE;
167    }
168
169    // loop on all ranges, and use the one which contains this momentum
170    // if none is found, the cut is not passed
171    Bool_t accept = kFALSE;
172    Int_t  i, goodRange = -1, nRanges = fRanges.GetEntriesFast();
173    for (i = 0; i < nRanges; i++) {
174       AliRsnPIDRange *range = (AliRsnPIDRange *)fRanges[i];
175       if (!range) continue;
176       if (!range->IsInRange(fTrackMom)) continue;
177       else {
178          goodRange = i;
179          accept = range->CutPass(fTrackNSigma);
180          AliDebugClass(2, Form("[%s] NSigma = %.3f, max = %.3f, track %s", GetName(), fTrackNSigma, range->NSigmaCut(), (accept ? "accepted" : "rejected")));
181          break;
182       }
183    }
184    if (goodRange < 0) {
185       AliDebugClass(2, Form("[%s] No good range found. Rejecting track", GetName()));
186       return kFALSE;
187    } else {
188       AliDebugClass(2, Form("[%s] Mom = %.3f, good range found (#%d), track was %s", GetName(), fTrackMom, goodRange, (accept ? "accepted" : "rejected")));
189       return accept;
190    }
191 }
192
193 //_________________________________________________________________________________________________
194 void AliRsnCutPIDNSigma::Print(const Option_t *) const
195 {
196 //
197 // Print information on this cut
198 //
199
200    Char_t mom[200], det[100], match[200];
201
202    if (fRejectUnmatched)
203       snprintf(match, 200, "Unmatched tracks are rejected");
204    else
205       snprintf(match, 200, "No check on track matching");
206
207    switch (fDetector) {
208       case kITS: snprintf(det, 3, "ITS"); break;
209       case kTPC: snprintf(det, 3, "TPC"); break;
210       case kTOF: snprintf(det, 3, "TOF"); break;
211       default  : snprintf(det, 3, "undefined");
212    }
213
214    AliInfo(Form("Cut name          : %s", GetName()));
215    AliInfo(Form("--> PID detector  : %s", det));
216    AliInfo(Form("--> match criteria: %s", match));
217    AliInfo(Form("--> momentum range: %s", mom));
218 }
219
220 //_________________________________________________________________________________________________
221 Bool_t AliRsnCutPIDNSigma::MatchITS(const AliVTrack *vtrack) const
222 {
223 //
224 // Checks if the track has the status flags required for an ITS standalone track
225 //
226
227    if ((vtrack->GetStatus() & AliESDtrack::kITSin)  == 0) return kFALSE;
228    if ((vtrack->GetStatus() & AliESDtrack::kITSpid) == 0) return kFALSE;
229
230    return kTRUE;
231 }
232 //_________________________________________________________________________________________________
233 Bool_t AliRsnCutPIDNSigma::MatchTPC(const AliVTrack *vtrack) const
234 {
235 //
236 // Checks if the track has the status flags required for a TPC track
237 //
238
239    if ((vtrack->GetStatus() & AliESDtrack::kTPCin) == 0) return kFALSE;
240
241    return kTRUE;
242 }
243
244 //_________________________________________________________________________________________________
245 Bool_t AliRsnCutPIDNSigma::MatchTOF(const AliVTrack *vtrack) const
246 {
247 //
248 // Checks if the track has the status flags required for an ITS standalone track
249 //
250
251    if ((vtrack->GetStatus() & AliESDtrack::kTOFout) == 0) return kFALSE;
252    if ((vtrack->GetStatus() & AliESDtrack::kTIME)   == 0) return kFALSE;
253
254    return kTRUE;
255 }
256
257 //_________________________________________________________________________________________________
258 Bool_t AliRsnCutPIDNSigma::MatchDetector(const AliVTrack *vtrack) const
259 {
260 //
261 // Checks if the track has matched the required detector.
262 // If no valid detector is specified, kFALSE is always returned.
263 //
264
265    switch (fDetector) {
266       case kITS: return MatchITS(vtrack);
267       case kTPC: return MatchTPC(vtrack);
268       case kTOF: return MatchTOF(vtrack);
269       default  : return kFALSE;
270    }
271 }
272
273 //_________________________________________________________________________________________________
274 void AliRsnCutPIDNSigma::AddPIDRange(Double_t nsigma, Double_t pmin, Double_t pmax)
275 {
276 //
277 // Add a new slot for checking PID
278 //
279
280    Int_t n = fRanges.GetEntries();
281
282    new (fRanges[n]) AliRsnPIDRange(nsigma, pmin, pmax);
283 }
284
285 //_________________________________________________________________________________________________
286 void AliRsnCutPIDNSigma::SinglePIDRange(Double_t nsigma)
287 {
288 //
289 // Clear all slots and sets a unique one
290 //
291
292    fRanges.Delete();
293
294    new (fRanges[0]) AliRsnPIDRange(nsigma, 0.0, 1E20);
295 }