1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 ////////////////////////////////////////////////////////////////////////////////
18 // This class contains all code which is used to compute any of the values
19 // which can be of interest within a resonance analysis. Besides the obvious
20 // invariant mass, it allows to compute other utility values on all possible
21 // targets, in order to allow a wide spectrum of binning and checks.
22 // When needed, this object can also define a binning in the variable which
23 // it is required to compute, which is used for initializing axes of output
24 // histograms (see AliRsnFunction).
25 // The value computation requires this object to be passed the object whose
26 // informations will be used. This object can be of any allowed input type
27 // (track, pair, event), then this class must inherit from AliRsnTarget.
28 // Then, when value computation is attempted, a check on target type is done
29 // and computation is successful only if expected target matches that of the
31 // In some cases, the value computation can require a support external object,
32 // which must then be passed to this class. It can be of any type inheriting
35 // authors: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
36 // M. Vala (martin.vala@cern.ch)
38 ////////////////////////////////////////////////////////////////////////////////
40 #include "Riostream.h"
48 #include "AliPIDResponse.h"
50 #include "AliRsnDaughter.h"
51 #include "AliRsnEvent.h"
53 #include "AliRsnMiniMonitor.h"
55 ClassImp(AliRsnMiniMonitor)
57 //__________________________________________________________________________________________________
58 AliRsnMiniMonitor::AliRsnMiniMonitor() :
71 //__________________________________________________________________________________________________
72 AliRsnMiniMonitor::AliRsnMiniMonitor(const char *name, EType type, Int_t cutID) :
81 // Default constructor
85 //__________________________________________________________________________________________________
86 AliRsnMiniMonitor::AliRsnMiniMonitor(const AliRsnMiniMonitor ©) :
90 fCharge(copy.fCharge),
91 fListID(copy.fListID),
99 //__________________________________________________________________________________________________
100 AliRsnMiniMonitor &AliRsnMiniMonitor::operator=(const AliRsnMiniMonitor ©)
103 // Assignment operator
106 TNamed::operator=(copy);
110 fCutID = copy.fCutID;
111 fCharge = copy.fCharge;
112 fListID = copy.fListID;
118 //__________________________________________________________________________________________________
119 Bool_t AliRsnMiniMonitor::Init(const char *name, TList *list)
122 // Initialize this output histogram and put into the passed list
133 AliError("No list!");
138 TH1 *histogram = 0x0;
142 sname += "_TrackPt_cut";
144 histogram = new TH1F(sname.Data(), "", 100, 0.0, 10.0);
147 sname += "_TPCsignal";
148 histogram = new TH2F(sname.Data(), "", 500, 0.0, 5.0, 1000, 0.0, 1000.0);
150 case ktimeTOFvsPPion:
151 sname += "_TOFsignalPi";
152 histogram = new TH2F(sname.Data(), "", 500, 0.0, 5.0, 1000, 0.0, 200000.0);
153 case ktimeTOFvsPKaon:
154 sname += "_TOFsignalK";
155 histogram = new TH2F(sname.Data(), "", 500, 0.0, 5.0, 1000, 0.0, 200000.0);
157 case ktimeTOFvsPProton:
158 sname += "_TOFsignalP";
159 histogram = new TH2F(sname.Data(), "", 500, 0.0, 5.0, 1000, 0.0, 200000.0);
162 AliError("Wrong enum type");
167 if (histogram && fList) {
169 fList->Add(histogram);
170 fListID = fList->IndexOf(histogram);
171 AliInfo(Form("Histogram '%s' added to list in slot #%d", histogram->GetName(), fListID));
178 //_____________________________________________________________________________
179 Bool_t AliRsnMiniMonitor::Fill(AliRsnDaughter *track, AliRsnEvent *event)
182 // Fill the histogram
185 // retrieve object from list
187 AliError("List pointer is NULL");
190 TObject *obj = fList->At(fListID);
192 AliError("List object is NULL");
196 Double_t valueX, valueY;
197 AliVTrack *vtrack = track->Ref2Vtrack();
199 AliPIDResponse *pid = event->GetPIDResponse();
204 AliWarning("Required vtrack for this value");
207 if (fCharge == '+' && vtrack->Charge() <= 0) return kFALSE;
208 if (fCharge == '-' && vtrack->Charge() >= 0) return kFALSE;
209 if (fCharge == '0' && vtrack->Charge() != 0) return kFALSE;
210 valueX = vtrack->Pt();
211 ((TH1F *)obj)->Fill(valueX);
215 AliWarning("Required vtrack for this value");
218 valueX = vtrack->GetTPCmomentum();
219 valueY = vtrack->GetTPCsignal();
220 ((TH2F *)obj)->Fill(valueX, valueY);
222 case ktimeTOFvsPPion:
224 AliWarning("Required vtrack for this value");
227 valueX = vtrack->P();
228 //valueY = vtrack->GetTOFsignal();
230 if (pid) valueY = pid->NumberOfSigmasTOF(vtrack, AliPID::kPion);
231 ((TH2F *)obj)->Fill(valueX, valueY);
233 case ktimeTOFvsPKaon:
235 AliWarning("Required vtrack for this value");
238 valueX = vtrack->P();
239 //valueY = vtrack->GetTOFsignal();
241 if (pid) valueY = pid->NumberOfSigmasTOF(vtrack, AliPID::kKaon);
242 ((TH2F *)obj)->Fill(valueX, valueY);
244 case ktimeTOFvsPProton:
246 AliWarning("Required vtrack for this value");
249 valueX = vtrack->P();
250 //valueY = vtrack->GetTOFsignal();
252 if (pid) valueY = pid->NumberOfSigmasTOF(vtrack, AliPID::kProton);
253 ((TH2F *)obj)->Fill(valueX, valueY);
256 AliError("Invalid value type");