New classes plus some renamed for a better user friendly interface
[u/mrichter/AliRoot.git] / PWG2 / RESONANCES / AliRsnSimpleFcnInvMass.cxx
1 //
2 // Class AliRsnSimpleFcnInvMass
3 //
4 // This is the most fundamental AliRsnSimpleFunction,
5 // which computes the invariant mass spectrum of a resonance,
6 // by correlating pairs of tracks from an event (or mixing, for BG)
7 //
8
9 #include <TH1.h>
10 #include <TArrayI.h>
11
12 #include "AliLog.h"
13
14 #include "AliRsnDaughter.h"
15 #include "AliRsnEvent.h"
16 #include "AliRsnSimpleFcnInvMass.h"
17
18 ClassImp(AliRsnSimpleFcnInvMass)
19
20 //________________________________________________________________________________________
21 AliRsnSimpleFcnInvMass::AliRsnSimpleFcnInvMass() :
22   AliRsnSimpleFunction(),
23   fUseMCValues(kFALSE),
24   fPIDMethod(AliRsnDaughter::kNoPID)
25 {
26 //
27 // Constructor.
28 // Only default initializations.
29 //
30 }
31
32 //________________________________________________________________________________________
33 AliRsnSimpleFcnInvMass::AliRsnSimpleFcnInvMass
34 (const char *name, AliRsnDaughter::EPIDMethod method, AliRsnPairDef *pd,
35  AliRsnHistoDef *hd, AliRsnCutMgr *cuts, Option_t *option) :
36   AliRsnSimpleFunction(name, pd, hd, cuts, option),
37   fUseMCValues(kFALSE),
38   fPIDMethod(method)
39 {
40 //
41 // Constructor.
42 // Only default initializations.
43 //
44 }
45
46 //________________________________________________________________________________________
47 Bool_t AliRsnSimpleFcnInvMass::ProcessOne(AliRsnEvent* event)
48 {
49 //
50 // Process a single event to build the invariant mass histogram
51 // from the correlation of track pairs, according to the AliRsnPairDef settings.
52 // The collection of used tracks depends on the PID method chosen.
53 //
54
55     if (!event) {
56         AliError("Argument cannot be NULL.");
57         return kFALSE;
58     }
59     
60     // check PID method
61     if (fPIDMethod < AliRsnDaughter::kNoPID || fPIDMethod >= AliRsnDaughter::kMethods) {
62         AliError("PID method not properly initialized");
63         return kFALSE;
64     }
65     
66     // check event cut
67     if (!CutPass(event)) return kFALSE;
68
69     // assign pointers to the list of indexes to be used
70     TArrayI *listCharged1 = 0x0, *listCharged2 = 0x0;
71     if (fPIDMethod == AliRsnDaughter::kNoPID || fPIDMethod == AliRsnDaughter::kWeighted) {
72         listCharged1 = event->GetCharged(fPairDef->GetCharge(0));
73         listCharged2 = event->GetCharged(fPairDef->GetCharge(1));
74     }
75     else {
76         listCharged1 = event->GetTracksArray(fPIDMethod, fPairDef->GetCharge(0), fPairDef->GetType(0));
77         listCharged2 = event->GetTracksArray(fPIDMethod, fPairDef->GetCharge(1), fPairDef->GetType(1));
78     }
79     if (!listCharged1 || !listCharged2) return 0.0;
80     TArrayI &list1 = *listCharged1;
81     TArrayI &list2 = *listCharged2;
82
83     Stat_t count = 0;
84     Int_t i1, i2, start2;
85     AliRsnDaughter *track1 = 0x0, *track2 = 0x0;
86
87     // loop on particle of type 1
88     for (i1 = 0; i1 < list1.GetSize(); i1++) {
89         track1 = event->GetTrack(list1[i1]);
90         if (!CutPass(track1)) continue;
91         // loop on particle of type 2
92         // in case we are building a like-sign histogram with particles
93         // of the same type, we must avoid that each pair is used twice
94         start2 = 0;
95         if (listCharged1 == listCharged2) start2 = i1 + 1;
96         for (i2 = start2; i2 < list2.GetSize(); i2++) {
97             track2 = event->GetTrack(list2[i2]);
98             if (!CutPass(track2)) continue;
99             if (Add(track1, track2)) count++;
100         }
101     }
102
103     return (count > (Stat_t)0);
104 }
105
106 //________________________________________________________________________________________
107 Bool_t AliRsnSimpleFcnInvMass::ProcessTwo
108 (AliRsnEvent* event1, AliRsnEvent* event2)
109 {
110 //
111 // Process a single event to build the invariant mass histogram
112 // from the correlation of track pairs, according to the AliRsnPairDef settings.
113 // The collection of used tracks depends on the PID method chosen.
114 //
115
116     if (!event1 || !event2) {
117         AliError("Arguments cannot be NULL.");
118         return kFALSE;
119     }
120     
121     // check event cut
122     if (!CutPass(event1)) return kFALSE;
123     if (!CutPass(event2)) return kFALSE;
124
125     // assign pointers to the list of indexes to be used
126     TArrayI *listCharged1 = 0x0, *listCharged2 = 0x0;
127     if (fPIDMethod == AliRsnDaughter::kNoPID || fPIDMethod == AliRsnDaughter::kWeighted) {
128         listCharged1 = event1->GetCharged(fPairDef->GetCharge(0));
129         listCharged2 = event2->GetCharged(fPairDef->GetCharge(1));
130     }
131     else {
132         listCharged1 = event1->GetTracksArray(fPIDMethod, fPairDef->GetCharge(0), fPairDef->GetType(0));
133         listCharged2 = event2->GetTracksArray(fPIDMethod, fPairDef->GetCharge(1), fPairDef->GetType(1));
134     }
135     if (!listCharged1 || !listCharged2) return 0.0;
136     TArrayI &list1 = *listCharged1;
137     TArrayI &list2 = *listCharged2;
138
139     Stat_t count = 0;
140     Int_t i1, i2;
141     AliRsnDaughter *track1 = 0x0, *track2 = 0x0;
142
143     // loop on particle of type 1
144     for (i1 = 0; i1 < list1.GetSize(); i1++) {
145         track1 = event1->GetTrack(list1[i1]);
146         if (!CutPass(track1)) continue;
147         // loop on particle of type 2
148         for (i2 = 0; i2 < list2.GetSize(); i2++) {
149             track2 = event2->GetTrack(list2[i2]);
150             if (!CutPass(track2)) continue;
151             if (Add(track1, track2)) count++;
152         }
153     }
154
155     return (count > (Stat_t)0);
156 }
157
158 //________________________________________________________________________________________
159 Bool_t AliRsnSimpleFcnInvMass::Add
160 (AliRsnDaughter *track1, AliRsnDaughter *track2)
161 {
162 //
163 // Add a histogram entry from two tracks, if they pass the cuts.
164 // The order matters, because track1 is processed using first element
165 // in the AliRsnPairDef definition, and track2 is processed using second ones.
166 // In case the study is done only for true pairs, this is checked automatically.
167 //
168
169     // set the static variable pointing to PID method
170     // to the current one in use by this object
171     AliRsnDaughter::SetPIDMethod(fPIDMethod);
172
173     // setup pair and check pair cuts
174     fPair.SetPair(track1, track2);
175     if (!CutPass(&fPair)) return kFALSE;
176
177     // computation variables
178     Double_t mass1[2] = {fPairDef->GetMass(0), fPairDef->GetMass(1)};
179     Double_t mass2[2] = {fPairDef->GetMass(1), fPairDef->GetMass(0)};
180     
181     // define the value of the entry weight according to PID method selected:
182     // - in case of "weighted" computations, which is done with all tracks
183     //   using the PID computed probabilities, the weight will be the product
184     //   of the ones related to the used PID for the tracks
185     Bool_t   useWeight = kFALSE;
186     Double_t entryWeight[2];
187     entryWeight[0] = fPairDef->ComputeWeight(track1, track2);
188     entryWeight[1] = fPairDef->ComputeWeight(track2, track1);
189     if (fPIDMethod != AliRsnDaughter::kWeighted) {
190         useWeight = kTRUE;
191         entryWeight[0] = 1.0;
192         entryWeight[0] = 1.0;
193     }
194     
195     // when filling the histogram, if the defined pair is a like-sign one,
196     // and if the two PID involved are different, it is necessary
197     // to fill it twice, exchanging the tracks with respect to pair definition
198     // in all other cases (unlike sign, like sign with same PID) this is not done
199     Int_t i, lastAdd = 0;
200     if (fPairDef->IsLikeSign() && !fPairDef->HasEqualTypes()) lastAdd = 1;
201     
202     // make computation
203     Double_t invmass;
204     for (i = 0; i <= lastAdd; i++) {
205         if (fUseMCValues) invmass = fPair.GetInvMassMC(mass1[i], mass2[i]);
206         else invmass = fPair.GetInvMass(mass1[i], mass2[i]);
207         
208         if (useWeight) fHisto1D->Fill(invmass, entryWeight[i]);
209         else fHisto1D->Fill(invmass);
210     }
211
212     return kTRUE;
213 }