AliPhysics  84d6c4a (84d6c4a)
AliFMDESDFixer.cxx
Go to the documentation of this file.
1 #include "AliFMDESDFixer.h"
2 #include "AliESDFMD.h"
3 #include "AliFMDStripIndex.h"
4 #include "AliForwardUtil.h"
6 #include "AliFMDCorrNoiseGain.h"
7 #include "AliLog.h"
8 #include <TH1.h>
9 #include <TList.h>
10 #include <TObjArray.h>
11 #include <TROOT.h>
12 #include <TMath.h>
13 #include <TSystem.h>
14 #include <TInterpreter.h>
15 #include <TVector3.h>
16 #include <iostream>
17 #include <iomanip>
18 
19 
20 //____________________________________________________________________
22  : TObject(),
23  fRecoFactor(1),
24  fMaxNoiseCorr(0.05),
25  fRecalculateEta(true),
26  fXtraDead(0),
27  fHasXtraDead(false),
28  fInvalidIsEmpty(false),
29  fNoiseChange(0),
30  fEtaChange(0),
31  fDeadChange(0)
32 {}
33 
34 //____________________________________________________________________
36  : TObject(),
37  fRecoFactor(1),
38  fMaxNoiseCorr(0.05),
39  fRecalculateEta(false),
40  fXtraDead(AliFMDStripIndex::Pack(3,'O',19,511)+1),
41  fHasXtraDead(false),
42  fInvalidIsEmpty(false),
43  fNoiseChange(0),
44  fEtaChange(0),
45  fDeadChange(0)
46 {}
47 
48 
49 //____________________________________________________________________
51  : TObject(),
58  fNoiseChange(0),
59  fEtaChange(0),
60  fDeadChange(0)
61 {}
62 
63 //____________________________________________________________________
66 {
67  if (&o == this) return *this;
68 
72  fXtraDead = o.fXtraDead;
75  fNoiseChange = 0;
76  fEtaChange = 0;
77  fDeadChange = 0;
78 
79  return *this;
80 };
81 
82 //____________________________________________________________________
83 void
85 {
86  TList* d = new TList;
87  d->SetOwner();
88  d->SetName(GetName());
89  l->Add(d);
90 
91  d->Add(AliForwardUtil::MakeParameter("recoFactor", fRecoFactor));
92  d->Add(AliForwardUtil::MakeParameter("recalcEta", fRecalculateEta));
93  d->Add(AliForwardUtil::MakeParameter("invalidIsEmpty", fInvalidIsEmpty));
94 
95  fNoiseChange = new TH1D("noiseChange", "#delta#Delta due to noise",30,0,.3);
96  fNoiseChange->SetDirectory(0);
97  fNoiseChange->SetFillColor(kYellow+1);
98  fNoiseChange->SetFillStyle(3001);
99  d->Add(fNoiseChange);
100 
101  fEtaChange = new TH1D("etaChange", "#delta#eta",40,-1,1);
102  fEtaChange->SetDirectory(0);
103  fEtaChange->SetFillColor(kCyan+1);
104  fEtaChange->SetFillStyle(3001);
105  d->Add(fEtaChange);
106 
107  fDeadChange = new TH1D("deadChange", "#deltaN_{dead}",100,0,51200);
108  fDeadChange->SetDirectory(0);
109  fDeadChange->SetFillColor(kMagenta+1);
110  fDeadChange->SetFillStyle(3001);
111  d->Add(fDeadChange);
112 
113  TObjArray* extraDead = new TObjArray;
114  extraDead->SetOwner();
115  extraDead->SetName("extraDead");
116  fXtraDead.Compact();
117  UInt_t firstBit = fXtraDead.FirstSetBit();
118  UInt_t nBits = fXtraDead.GetNbits();
119  for (UInt_t i = firstBit; i < nBits; i++) {
120  if (!fXtraDead.TestBitNumber(i)) continue;
121  UShort_t dd, s, t;
122  Char_t r;
123  AliFMDStripIndex::Unpack(i, dd, r, s, t);
124  extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
125  dd, r, s, t),
126  UShort_t(i)));
127  }
128  d->Add(extraDead);
129  fHasXtraDead = nBits > 0;
130 
131 }
132 
133 //____________________________________________________________________
134 void
136 {
137  if (d < 1 || d > 3) {
138  Warning("AddDead", "Invalid detector FMD%d", d);
139  return;
140  }
141  Bool_t inner = (r == 'I' || r == 'i');
142  if (d == 1 && !inner) {
143  Warning("AddDead", "Invalid ring FMD%d%c", d, r);
144  return;
145  }
146  if ((inner && s >= 20) || (!inner && s >= 40)) {
147  Warning("AddDead", "Invalid sector FMD%d%c[%02d]", d, r, s);
148  return;
149  }
150  if ((inner && t >= 512) || (!inner && t >= 256)) {
151  Warning("AddDead", "Invalid strip FMD%d%c[%02d,%03d]", d, r, s, t);
152  return;
153  }
154 
155  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
156  // Int_t i = 0;
157  fXtraDead.SetBitNumber(id, true);
158 }
159 //____________________________________________________________________
160 void
162  UShort_t s1, UShort_t s2,
163  UShort_t t1, UShort_t t2)
164 {
165  // Add a dead region spanning from FMD<d><r>[<s1>,<t1>] to
166  // FMD<d><r>[<s2>,<t2>] (both inclusive)
167  for (Int_t s = s1; s <= s2; s++)
168  for (Int_t t = t1; t <= t2; t++)
169  AddDead(d, r, s, t);
170 }
171 //____________________________________________________________________
172 void
174 {
175  if (!script || script[0] == '\0') return;
176 
177  const char* scr = gSystem->Which(gROOT->GetMacroPath(), script);
178  if (!scr) {
179  AliWarningF("%s not found in %s", script, gROOT->GetMacroPath());
180  return;
181  }
182  AliInfoF("Reading additional dead strips from %s", scr);
183 
184  gROOT->Macro(Form("%s((AliFMDESDFixer*)%p);", scr, this));
185 
186  gInterpreter->UnloadFile(scr);
187  delete scr;
188 }
189 
190 //____________________________________________________________________
191 Bool_t
193 {
194  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
195  return fXtraDead.TestBitNumber(id);
196 }
197 
198 //____________________________________________________________________
199 Int_t
201 {
202  if (!IsUseNoiseCorrection())
203  // If the reconstruction factor was high (4 or more), do nothing
204  return 0;
205 
206  Int_t target = 0;
207  if (AliESDFMD::Class_Version() < 4) {
208  // IF we running with older STEER - we fix it here
209  target = 4;
210  } else {
211 #if 1
212  if (!esd.TestBit(1 << 14)) {
213  // If the bit isn't set, do nothing
214  return 0;
215  }
216 #else
217  // Uncommented until Peter commits patch to STEER/ESD
218  if (!esd.NeedNoiseFix()) {
219  // If the bit isn't set, do nothing
220  return 0;
221  }
222 #endif
223  target = Int_t(esd.GetNoiseFactor());
224  }
225  // Get the target factor - even thought the method below returns a
226  // floating point value, we know that the noise factor is always
227  // integer, so we coerce it to be the same here.
228  target -= fRecoFactor;
229 
230  // If the target factor is the same or smaller than the assumed
231  // factor, we have nothing to do here, and we return immediately
232  if (target <= 0) return 0;
233 
234  // Get the scaled noise from the correction mananger
236  return 0;
237 
238  return target;
239 }
240 
241 //____________________________________________________________________
242 #define ETA2COS(ETA) \
243  TMath::Cos(2*TMath::ATan(TMath::Exp(-TMath::Abs(ETA))))
244 
245 //____________________________________________________________________
246 void
247 AliFMDESDFixer::Fix(AliESDFMD& esd, const TVector3& ip)
248 {
249 
250  const AliFMDCorrNoiseGain* ng = 0;
251  Int_t tgtFactor = FindTargetNoiseFactor(esd, false);
252  if (tgtFactor > 0)
254 
255  if (!ng && !fHasXtraDead && !fRecalculateEta && !fInvalidIsEmpty)
256  // We have nothing to do!
257  return;
258 
259  UShort_t nDead = 0;
260  for (UShort_t d = 1; d <= 3; d++) {
261  UShort_t nQ = d == 1 ? 1 : 2;
262 
263  for (UShort_t q = 0; q < nQ; q++) {
264  Char_t r = (q == 0 ? 'I' : 'O');
265  UShort_t nS = (q == 0 ? 20 : 40);
266  UShort_t nT = (q == 0 ? 512 : 256);
267 
268  for (UShort_t s = 0; s < nS; s++) {
269  for (UShort_t t = 0; t < nT; t++) {
270  Double_t mult = esd.Multiplicity(d,r,s,t);
271  Double_t eta = esd.Eta(d,r,s,t);
272  Double_t cosTheta = 0;
273 
274  if (CheckDead(d,r,s,t,mult)) nDead++;
275 
276  // Possibly re-calculate eta
277  if (fRecalculateEta) RecalculateEta(d,r,s,t,ip,eta,mult,cosTheta);
278 
279  // Possibly correct for poor treatment of ZS in reconstruction.
280  if (ng && mult != AliESDFMD::kInvalidMult) {
281  if (cosTheta <= 0) cosTheta = ETA2COS(eta);
282  if (!NoiseCorrect(tgtFactor,ng->Get(d,r,s,t), cosTheta, mult))
283  nDead++;
284  }
285 
286  // Write out final values to object
287  if (mult >= AliESDFMD::kInvalidMult) mult = AliESDFMD::kInvalidMult;
288  esd.SetMultiplicity(d,r,s,t,mult);
289  esd.SetEta(d,r,s,t,eta);
290  } // for t
291  } // for s
292  } // for q
293  } // for d
294  fDeadChange->Fill(nDead);
295 }
296 
297 //____________________________________________________________________
298 Bool_t
300  Double_t& mult)
301 {
302  // Correct for zero's being flagged as invalid
303  if (mult == AliESDFMD::kInvalidMult && fInvalidIsEmpty) mult = 0;
304 
305  // Take into account what we're defined as dead
306  if (IsDead(d,r,s,t)) {
307  mult = AliESDFMD::kInvalidMult;
308  return true;
309  }
310  return false;
311 }
312 
313 //____________________________________________________________________
314 void
316  const TVector3& ip, Double_t& eta,
317  Double_t& mult,
318  Double_t& cosTheta)
319 {
320  Double_t oldEta = eta, newEta = eta, newPhi=0;
321  // Double_t newEta = AliForwardUtil::GetEtaFromStrip(d,r,s,t, zvtx);
322  // eta = newEta;
323  cosTheta = ETA2COS(eta);
324  if (!AliForwardUtil::GetEtaPhi(d, r, s, t, ip, newEta, newPhi)) return;
325  if (TMath::Abs(eta) < 1) {
326  ::Warning("RecalculateEta",
327  "FMD%d%c[%2d,%3d] (%f,%f,%f) eta=%f phi=%f (was %f)",
328  d, r, s, t, ip.X(), ip.Y(), ip.Z(), newEta, newPhi, oldEta);
329  return;
330  }
331 
332  eta = newEta;
333  fEtaChange->Fill(newEta-oldEta);
334 
335  if (mult == AliESDFMD::kInvalidMult) return;
336 
337  Double_t newCos = ETA2COS(newEta);
338  Double_t oldCos = ETA2COS(oldEta);
339  Double_t corr = newCos / oldCos;
340  cosTheta = newCos;
341  mult *= corr;
342 }
343 //____________________________________________________________________
344 Bool_t
346  Double_t& mult)
347 {
348  if (corr > fMaxNoiseCorr || corr <= 0) {
349  mult = AliESDFMD::kInvalidMult;
350  return false;
351  }
352  Double_t add = corr * target * cosTheta;
353  fNoiseChange->Fill(add);
354 
355  mult += add;
356  return true;
357 }
358 
359 #define PF(N,V,...) \
360  AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
361 #define PFB(N,FLAG) \
362  do { \
363  AliForwardUtil::PrintName(N); \
364  std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
365  } while(false)
366 #define PFV(N,VALUE) \
367  do { \
368  AliForwardUtil::PrintName(N); \
369  std::cout << (VALUE) << std::endl; } while(false)
370 
371 //____________________________________________________________________
372 void
374 {
376  gROOT->IncreaseDirLevel();
377  PFB("Consider invalid null", fInvalidIsEmpty);
378  PFB("Has extra dead", fHasXtraDead || fXtraDead.GetNbits() > 0);
379  PFV("Reco noise factor", fRecoFactor);
380  PFV("Max noise corr", fMaxNoiseCorr);
381  PFB("Recalc. eta", fRecalculateEta);
382  gROOT->DecreaseDirLevel();
383 }
384 
385 //____________________________________________________________________
386 //
387 // EOF
388 //
Bool_t NoiseCorrect(Int_t f, Double_t c, Double_t cosTheta, Double_t &mult)
void RecalculateEta(UShort_t d, Char_t r, UShort_t s, UShort_t t, const TVector3 &ip, Double_t &mult, Double_t &eta, Double_t &cosTheta)
AliFMDESDFixer & operator=(const AliFMDESDFixer &)
double Double_t
Definition: External.C:58
void Print(Option_t *option="") const
void AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t)
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)
TSystem * gSystem
char Char_t
Definition: External.C:18
Float_t Get(UShort_t d, Char_t r, UShort_t s, UShort_t t) const
const char * GetName() const
Int_t FindTargetNoiseFactor(const AliESDFMD &esd, Bool_t check=true) const
Bool_t CheckDead(UShort_t d, Char_t r, UShort_t s, UShort_t t, Double_t &m)
Bool_t fInvalidIsEmpty
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
Double_t fMaxNoiseCorr
Various utilities used in PWGLF/FORWARD.
Definition: External.C:212
static UInt_t Pack(UShort_t det, Char_t rng, UShort_t sec, UShort_t str)
virtual Bool_t IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const
#define PFB(N, FLAG)
void AddDeadRegion(UShort_t d, Char_t r, UShort_t s1, UShort_t s2, UShort_t t1, UShort_t t2)
const AliFMDCorrNoiseGain * GetNoiseGain() const
void CreateOutputObjects(TList *l)
static void PrintTask(const TObject &o)
#define PFV(N, VALUE)
static TObject * MakeParameter(const char *name, UShort_t value)
#define ETA2COS(ETA)
Bool_t fRecalculateEta
unsigned short UShort_t
Definition: External.C:28
const char Option_t
Definition: External.C:48
bool Bool_t
Definition: External.C:53
void Fix(AliESDFMD &esd, const TVector3 &ip)
static AliForwardCorrectionManager & Instance()
static void Unpack(UInt_t id, UShort_t &det, Char_t &rng, UShort_t &sec, UShort_t &str)
Bool_t IsUseNoiseCorrection() const