AliPhysics  1168478 (1168478)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliFMDMCTrackDensity.cxx
Go to the documentation of this file.
1 #include "AliFMDMCTrackDensity.h"
2 #include "AliESDFMD.h"
3 #include "AliTrackReference.h"
4 #include <TMath.h>
5 #include "AliFMDStripIndex.h"
6 #include "AliFMDEncodedEdx.h"
7 #include "AliLog.h"
8 #include <TH2D.h>
9 #include <TH1D.h>
10 #include <TList.h>
11 #include <TROOT.h>
12 #include <iostream>
13 
14 //____________________________________________________________________
15 void
17 {
18  angle = 0;
19  oldDetector = 0;
20  oldRing = '\0';
21  oldSector = 1024;
22  oldStrip = 1024;
23  startStrip = 1024;
24  nRefs = 0;
25  nStrips = 0;
26  longest = 0x0;
27  if (alsoCount) count = 0;
28 }
29 
30 //____________________________________________________________________
33 {
34  if (&o == this) return *this;
35  angle = o.angle;
36  oldDetector = o.oldDetector;
37  oldRing = o.oldRing;
38  oldSector = o.oldSector;
39  oldStrip = o.oldStrip;
40  startStrip = o.startStrip;
41  nRefs = o.nRefs;
42  nStrips = o.nStrips;
43  count = o.count;
44  longest = o.longest;
45  return *this;
46 }
47 
48 //____________________________________________________________________
51  fState(),
53  fLowCutvalue(0),
54  fNr(0),
55  fNt(0),
56  fNc(0),
57  fNcr(0),
58  fOutput(0)
59 {
60  // Default constructor
61 }
62 
63 //____________________________________________________________________
65  : AliBaseMCTrackDensity("fmdMCTrackDensity"),
66  fState(),
67  fMaxConsequtiveStrips(3),
68  fLowCutvalue(0),
69  fNr(0),
70  fNt(0),
71  fNc(0),
72  fNcr(0),
73  fOutput(0)
74 {
75  // Normal constructor constructor
76 }
77 
78 //____________________________________________________________________
81  fState(o.fState),
82  fMaxConsequtiveStrips(o.fMaxConsequtiveStrips),
83  fLowCutvalue(o.fLowCutvalue),
84  fNr(o.fNr),
85  fNt(o.fNt),
86  fNc(o.fNc),
87  fNcr(o.fNcr),
88  fOutput(o.fOutput)
89 {
90  // Normal constructor constructor
91 }
92 
93 //____________________________________________________________________
96 {
97  // Assignment operator
98  if (&o == this) return *this;
102  fNr = o.fNr;
103  fNt = o.fNt;
104  fNc = o.fNc;
105  fNcr = o.fNcr;
106  fState = o.fState;
107  fOutput = o.fOutput;
108 
109  return *this;
110 }
111 
112 //____________________________________________________________________
113 void
115 {
117  TList* ll = static_cast<TList*>(l->FindObject(GetTitle()));
118  if (!ll) ll = l;
119 
120  fNr = new TH1D("clusterRefs", "# track references per cluster",
121  21, -.5, 20.5);
122  fNr->SetXTitle("# of track references/cluster");
123  fNr->SetFillColor(kRed+1);
124  fNr->SetFillStyle(3001);
125  fNr->SetDirectory(0);
126  ll->Add(fNr);
127 
128  fNt = new TH1D("clusterSize", "cluster length in strips", 21, -.5, 20.5);
129  fNt->SetXTitle("Cluster size [strips]");
130  fNt->SetFillColor(kBlue+1);
131  fNt->SetFillStyle(3001);
132  fNt->SetDirectory(0);
133  ll->Add(fNt);
134 
135  fNc = new TH1D("nClusters", "# clusters per track", 21, -.5, 20.5);
136  fNc->SetXTitle("# clusters");
137  fNc->SetFillColor(kGreen+1);
138  fNc->SetFillStyle(3001);
139  fNc->SetDirectory(0);
140  ll->Add(fNc);
141 
142  fNcr = new TH2D("clusterVsRefs", "# clusters vs. # refs",
143  21, -.5, 20.5, 21, -.5, 20.5);
144  fNcr->SetXTitle("# References");
145  fNcr->SetYTitle("# Clusters");
146  fNcr->SetOption("COLZ");
147  fNcr->SetDirectory(0);
148  ll->Add(fNcr);
149 
150 
151 }
152 //____________________________________________________________________
153 Int_t
155 {
156  return AliTrackReference::kFMD;
157 }
158 
159 //____________________________________________________________________
160 void
162 {
163  fState.Clear(true);
164 }
165 
166 //____________________________________________________________________
167 void
169 {
170  fNc->Fill(fState.count);
171  fNcr->Fill(nRefs, fState.count);
172  fState.Clear(true);
173 }
174 
175 //____________________________________________________________________
176 AliTrackReference*
177 AliFMDMCTrackDensity::ProcessRef(AliMCParticle* particle,
178  const AliMCParticle* mother,
179  AliTrackReference* ref)
180 {
181  // Process track references of a track
182  //
183  // Note: If particle refers to a primary, then particle and mother
184  // refers to the same particle (the address are the same)
185  //
186 
187  // Get the detector coordinates
188  UShort_t d, s, t;
189  Char_t r;
190  AliFMDStripIndex::Unpack(ref->UserId(), d, r, s, t);
191  Double_t edep, length, dEdep, dLength;
192  AliFMDEncodedEdx::Decode((ref->UserId() >> 19), edep, length, dEdep, dLength);
193 
194  Double_t normaldEdx=0.0;
195  if(length>0.0)
196  normaldEdx=(edep/length)/4.406; // 4.406 mip in Si per 1 cm
197 
198  // Calculate distance of previous reference to base of cluster
199  UShort_t nT = TMath::Abs(t - fState.startStrip) + 1;
200 
201  // Now check if we should flush to output
202  Bool_t used = false;
203 
204  // If this is a new detector/ring, then reset the other one
205  // Check if we have a valid old detectorm ring, and sector
206  if (fState.oldDetector > 0 &&
207  fState.oldRing != '\0' &&
208  fState.oldSector != 1024) {
209  // New detector, new ring, or new sector
210  if (d != fState.oldDetector ||
211  r != fState.oldRing ||
212  s != fState.oldSector) {
213  if (fDebug) Info("Process", "New because new sector");
214  used = true;
215  }
216  else if (nT > fMaxConsequtiveStrips) {
217  if (fDebug) Info("Process", "New because too long: %d (%d,%d,%d)",
219  used = true;
220  }
221  }
222  if (used) {
223  if (fDebug)
224  Info("Process", "I=%p L=%p D=%d (was %d), R=%c (was %c), "
225  "S=%2d (was %2d) t=%3d (was %3d) nT=%3d/%4d",
226  ref, fState.longest,
227  d, fState.oldDetector,
228  r, fState.oldRing,
229  s, fState.oldSector,
230  t, fState.oldStrip,
232  // Int_t nnT = TMath::Abs(fState.oldStrip - fState.startStrip) + 1;
233  StoreParticle(particle, mother, fState.longest);
234  fState.Clear(false);
235  }
236 
237  if(normaldEdx<fLowCutvalue)
238  return 0x0;
239  // If base of cluster not set, set it here.
240  if (fState.startStrip == 1024) fState.startStrip = t;
241 
242  // Calculate distance of previous reference to base of cluster
243  fState.nStrips = TMath::Abs(t - fState.startStrip) + 1;
244 
245  // Count number of track refs in this sector
246  fState.nRefs++;
247 
248  fState.oldDetector = d;
249  fState.oldRing = r;
250  fState.oldSector = s;
251  fState.oldStrip = t;
252 
253  // Debug output
254  if (fDebug) {
255  if (t == fState.startStrip)
256  Info("Process", "New cluster starting at FMD%d%c[%2d,%3d]",
257  d, r, s, t);
258  else
259  Info("Process", "Adding to cluster starting at FMD%d%c[%2d,%3d], "
260  "length=%3d (now in %3d, previous %3d)",
262  }
263 
264  // The longest passage is determined through the angle
265  Double_t ang = GetTrackRefTheta(ref);
266  if (ang > fState.angle) {
267  fState.longest = ref;
268  fState.angle = ang;
269  }
270 
271  return fState.longest;
272 }
273 
274 //____________________________________________________________________
275 Double_t
276 AliFMDMCTrackDensity::StoreParticle(AliMCParticle* particle,
277  const AliMCParticle* mother,
278  AliTrackReference* ref) const
279 {
280  Double_t w =
281  AliBaseMCTrackDensity::StoreParticle(particle, mother, ref);
282  if (w <= 0) return w;
283 
284  // Get the detector coordinates
285  UShort_t d, s, t;
286  Char_t r;
287  AliFMDStripIndex::Unpack(ref->UserId(), d, r, s, t);
288 
289  // Check if we have value already
290  Double_t old = fOutput->Multiplicity(d,r,s,t);
291 
292  // If invalid, force it valid
293  if (old == AliESDFMD::kInvalidMult) old = 0;
294 
295  // Increment count
296 
297  fOutput->SetMultiplicity(d,r,s,t,old+w);
298 
299  // Fill histograms
300  fNr->Fill(fState.nRefs);
301  fNt->Fill(fState.nStrips);
302 
303  fState.count++;
304 
305  return w;
306 }
307 
308 //____________________________________________________________________
309 Bool_t
311  const AliMCEvent& event,
312  const TVector3& ip,
313  AliESDFMD& output,
314  TH2D* primary)
315 {
316  //
317  // Filter the input kinematics and track references, using
318  // some of the ESD information
319  //
320  // Parameters:
321  // input Input ESD event
322  // event Input MC event
323  // vz Vertex position
324  // output Output ESD-like object
325  // primary Per-event histogram of primaries
326  //
327  // Return:
328  // True on succes, false otherwise
329  //
330  output.Clear();
331  fOutput = &output;
332 
333  // Copy eta values to output
334  for (UShort_t ed = 1; ed <= 3; ed++) {
335  UShort_t nq = (ed == 1 ? 1 : 2);
336  for (UShort_t eq = 0; eq < nq; eq++) {
337  Char_t er = (eq == 0 ? 'I' : 'O');
338  UShort_t ns = (eq == 0 ? 20 : 40);
339  UShort_t nt = (eq == 0 ? 512 : 256);
340  for (UShort_t es = 0; es < ns; es++) {
341  for (UShort_t et = 0; et < nt; et++) {
342  Double_t eta, phi;
343  AliForwardUtil::GetEtaPhi(ed, er, es, et, ip, eta, phi);
344  output.SetEta(ed, er, es, et, eta); // input.Eta(ed, er, es, et);
345  // output.SetPhi(ed, er, es, et, phi*TMath::RadToDeg());
346  }
347  }
348  }
349  }
350 
351  return ProcessTracks(event, ip, primary);
352 }
353 
354 #define PFV(N,VALUE) \
355  do { \
356  AliForwardUtil::PrintName(N); \
357  std::cout << (VALUE) << std::endl; } while(false)
358 //____________________________________________________________________
359 void
361 {
363  gROOT->IncreaseDirLevel();
364  PFV("Max cluster size", fMaxConsequtiveStrips);
365  gROOT->DecreaseDirLevel();
366 }
367 
368 //____________________________________________________________________
369 //
370 // EOF
371 //
#define PFV(N, VALUE)
State & operator=(const State &o)
double Double_t
Definition: External.C:58
void Clear(Bool_t alsoCount=false)
Longest track through.
static Bool_t GetEtaPhi(UShort_t det, Char_t ring, UShort_t sec, UShort_t str, const TVector3 &ip, Double_t &eta, Double_t &phi)
virtual void CreateOutputObjects(TList *list)
char Char_t
Definition: External.C:18
Bool_t Calculate(const AliESDFMD &esd, const AliMCEvent &event, const TVector3 &ip, AliESDFMD &output, TH2D *primary)
AliTrackReference * ProcessRef(AliMCParticle *particle, const AliMCParticle *mother, AliTrackReference *ref)
void EndTrackRefs(Int_t nRefs)
void CreateOutputObjects(TList *list)
virtual Double_t StoreParticle(AliMCParticle *particle, const AliMCParticle *mother, AliTrackReference *ref) const
virtual void Print(Option_t *option="") const
Double_t GetTrackRefTheta(const AliTrackReference *ref) const
int Int_t
Definition: External.C:63
Definition: External.C:228
Definition: External.C:212
static void Decode(UInt_t bits, Double_t &edep, Double_t &length)
Class to encode the energy loss and path length of a particle into track reference bits...
AliBaseMCTrackDensity & operator=(const AliBaseMCTrackDensity &o)
void Print(Option_t *option="") const
Double_t StoreParticle(AliMCParticle *particle, const AliMCParticle *mother, AliTrackReference *ref) const
UShort_t fMaxConsequtiveStrips
State.
unsigned short UShort_t
Definition: External.C:28
const char Option_t
Definition: External.C:48
Bool_t ProcessTracks(const AliMCEvent &event, const TVector3 &ip, TH2D *primary)
bool Bool_t
Definition: External.C:53
AliFMDMCTrackDensity & operator=(const AliFMDMCTrackDensity &o)
struct AliFMDMCTrackDensity::State fState
static void Unpack(UInt_t id, UShort_t &det, Char_t &rng, UShort_t &sec, UShort_t &str)