AliRoot Core  edcc906 (edcc906)
DrawHitsSDigits.C
Go to the documentation of this file.
1 //____________________________________________________________________
2 //
3 // $Id: DrawHitsSDigits.C 20907 2007-09-25 08:44:03Z cholm $
4 //
5 // Script that contains a class to draw eloss from hits, versus ADC
6 // counts from sdigits, using the AliFMDInputHits class in the util library.
7 //
8 // It draws the energy loss versus the p/(mq^2). It can be overlayed
9 // with the Bethe-Bloc curve to show how the simulation behaves
10 // relative to the expected.
11 //
12 // Use the script `Compile.C' to compile this class using ACLic.
13 //
14 #include <TH2D.h>
15 #include <AliFMDHit.h>
16 #include <AliFMDSDigit.h>
17 #include <AliFMDInput.h>
18 #include <AliFMDEdepMap.h>
19 #include <iostream>
20 #include <TStyle.h>
21 #include <TArrayF.h>
22 #include <AliLog.h>
23 #include <TMath.h>
24 
36 {
37 private:
38  TH2D* fElossVsAdc; // Histogram
40 public:
41  //__________________________________________________________________
42  DrawHitsSDigits(Int_t n=900, Double_t emin=1e-3, Double_t emax=10,
43  Int_t m=1100, Double_t amin=-0.5, Double_t amax=1099.5)
44  {
46  AddLoad(kHits);
47  TArrayF eloss(MakeLogScale(n, emin, emax));
48  TArrayF adcs(m+1);
49  adcs[0] = amin;
50  for (Int_t i = 1; i < m+1; i++) adcs[i] = adcs[i-1] + (amax-amin)/m;
51  fElossVsAdc = new TH2D("bad", "#Delta E vs. ADC",
52  eloss.fN-1, eloss.fArray, adcs.fN-1, adcs.fArray);
53  fElossVsAdc->SetXTitle("#Delta E/#Delta x [MeV/cm]");
54  fElossVsAdc->SetYTitle("ADC value");
55  }
56  //__________________________________________________________________
60  Bool_t Begin(Int_t ev)
61  {
62  fMap.Reset();
63  return AliFMDInput::Begin(ev);
64  }
65  //__________________________________________________________________
66  Bool_t ProcessHit(AliFMDHit* hit, TParticle*)
67  {
68  Info("ProcessHit", "Processing hit %p", hit);
69  // Cache the energy loss
70  if (!hit) return kFALSE;
71  UShort_t det = hit->Detector();
72  Char_t rng = hit->Ring();
73  UShort_t sec = hit->Sector();
74  UShort_t str = hit->Strip();
75  if (str > 511) {
76  AliWarning(Form("Bad strip number %d in hit", str));
77  return kTRUE;
78  }
79  fMap(det, rng, sec, str).fEdep += hit->Edep();
80  fMap(det, rng, sec, str).fN++;
81  hit->Print();
82  return kTRUE;
83  }
84  //__________________________________________________________________
85  Bool_t ProcessSDigit(AliFMDSDigit* sdigit)
86  {
87  Info("ProcessSDigit", "Processing sdigit %p", sdigit);
88  if (!sdigit) return kFALSE;
89  UShort_t det = sdigit->Detector();
90  Char_t rng = sdigit->Ring();
91  UShort_t sec = sdigit->Sector();
92  UShort_t str = sdigit->Strip();
93  if (str > 511) {
94  AliWarning(Form("Bad strip number %d in sdigit", str));
95  return kFALSE;
96  }
97  fElossVsAdc->Fill(fMap(det, rng, sec, str).fEdep, sdigit->Counts());
98  sdigit->Print();
99  return kTRUE;
100  }
101  //__________________________________________________________________
102  Bool_t Finish()
103  {
104  gStyle->SetPalette(1);
105  gStyle->SetOptTitle(0);
106  gStyle->SetCanvasColor(0);
107  gStyle->SetCanvasBorderSize(0);
108  gStyle->SetPadColor(0);
109  gStyle->SetPadBorderSize(0);
110  fElossVsAdc->SetStats(kFALSE);
111  fElossVsAdc->Draw("COLZ");
112  return kTRUE;
113  }
114 
115  ClassDef(DrawHitsSDigits,0);
116 };
117 
118 //____________________________________________________________________
119 //
120 // EOF
121 //
AliFMDEdepMap fMap
Bool_t Begin(Int_t ev)
Map of Energy deposited, hit information per strip. Contains a pair of energy deposited fEdep and num...
Definition: AliFMDEdepMap.h:33
TStyle * gStyle
Hit in the FMD.
FMD utility classes for reading FMD data.
UShort_t Sector() const
Char_t Ring() const
Per strip map of energy deposited and number of hits.
UShort_t Counts() const
Definition: AliFMDSDigit.h:134
Bool_t ProcessHit(AliFMDHit *hit, TParticle *)
virtual Bool_t Begin(Int_t event)
Draw hit energy loss versus sdigit ADC.
#define AliWarning(message)
Definition: AliLog.h:541
Float_t Edep() const
Definition: AliFMDHit.h:82
Char_t Ring() const
Definition: AliFMDHit.h:76
static TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max)
UShort_t Detector() const
Definition: AliFMDHit.h:74
DrawHitsSDigits(Int_t n=900, Double_t emin=1e-3, Double_t emax=10, Int_t m=1100, Double_t amin=-0.5, Double_t amax=1099.5)
UShort_t Strip() const
virtual void AddLoad(ETrees tree)
Definition: AliFMDInput.h:134
UShort_t Sector() const
Definition: AliFMDHit.h:78
void Print(Option_t *opt="") const
Definition: AliFMDHit.cxx:186
Digits for the FMD.
virtual void Reset()
Base class for reading in various FMD data. The class loops over all found events. For each event the specified data is read in. The class then loops over all elements of the read data, and process these with user defined code.
Definition: AliFMDInput.h:106
void Print(Option_t *opt="") const
UShort_t Strip() const
Definition: AliFMDHit.h:80
class for summable digits
Definition: AliFMDSDigit.h:27
AliFMDhit is the hit class for the FMD. Hits are the information that comes from a Monte Carlo at eac...
Definition: AliFMDHit.h:30
UShort_t Detector() const
Bool_t ProcessSDigit(AliFMDSDigit *sdigit)