More docs
[u/mrichter/AliRoot.git] / FMD / scripts / DrawHitsDigits.C
1 //____________________________________________________________________
2 //
3 // $Id$
4 //
5 // Script that contains a class to draw eloss from hits, versus ADC
6 // counts from digits, 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 <AliFMDDigit.h>
17 #include <AliFMDInput.h>
18 #include <AliFMDEdepMap.h>
19 #include <iostream>
20 #include <TStyle.h>
21 #include <TArrayF.h>
22
23 /** @class DrawHitsDigits
24     @brief Draw hit energy loss versus digit ADC
25     @code 
26     Root> .L Compile.C
27     Root> Compile("DrawHitsDigits.C")
28     Root> DrawHitsDigits c
29     Root> c.Run();
30     @endcode
31     @ingroup FMD_script
32  */
33 class DrawHitsDigits : public AliFMDInputHits
34 {
35 private:
36   TH2D* fElossVsAdc; // Histogram 
37   AliFMDEdepMap fMap;
38 public:
39   //__________________________________________________________________
40   TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max) 
41   {
42     TArrayF bins(n+1);
43     Float_t dp   = n / TMath::Log10(max / min);
44     Float_t pmin = TMath::Log10(min);
45     bins[0]      = min;
46     for (Int_t i = 1; i < n+1; i++) {
47       Float_t p = pmin + i / dp;
48       bins[i]   = TMath::Power(10, p);
49     }
50     return bins;
51   }
52   //__________________________________________________________________
53   DrawHitsDigits(Int_t n=900, Double_t emin=1e-3, Double_t emax=10, 
54                  Int_t m=1100, Double_t amin=-0.5, Double_t amax=1099.5) 
55   { 
56     AddLoad(kDigits);
57     TArrayF eloss(MakeLogScale(n, emin, emax));
58     TArrayF adcs(m+1);
59     adcs[0] = amin;
60     for (Int_t i = 1; i < m+1; i++) adcs[i] = adcs[i-1] + (amax-amin)/m;
61     fElossVsAdc = new TH2D("bad", "#Delta E vs. ADC", 
62                            eloss.fN-1, eloss.fArray, adcs.fN-1, adcs.fArray);
63     fElossVsAdc->SetXTitle("#Delta E/#Delta x [MeV/cm]");
64     fElossVsAdc->SetYTitle("ADC value");
65   }
66   //__________________________________________________________________
67   /** Begining of event
68       @param ev Event number
69       @return @c false on error */
70   Bool_t Begin(Int_t ev) 
71   {
72     fMap.Reset();
73     return AliFMDInputHits::Begin(ev);
74   }
75   //__________________________________________________________________
76   Bool_t ProcessHit(AliFMDHit* hit, TParticle*) 
77   {
78     // Cache the energy loss 
79     if (!hit) return kFALSE;
80     UShort_t det = hit->Detector();
81     Char_t   rng = hit->Ring();
82     UShort_t sec = hit->Sector();
83     UShort_t str = hit->Strip();
84     if (str > 511) {
85       AliWarning(Form("Bad strip number %d in hit", str));
86       return kTRUE;
87     }
88     fMap(det, rng, sec, str).fEdep += hit->Edep();
89     fMap(det, rng, sec, str).fN++;
90     return kTRUE;
91   }
92   //__________________________________________________________________
93   Bool_t ProcessDigit(AliFMDDigit* digit)
94   {
95     if (!digit) continue;
96     UShort_t det = digit->Detector();
97     Char_t   rng = digit->Ring();
98     UShort_t sec = digit->Sector();
99     UShort_t str = digit->Strip();
100     if (str > 511) {
101       AliWarning(Form("Bad strip number %d in digit", str));
102       continue;
103     }
104     fElossVsAdc->Fill(fMap(det, rng, sec, str).fEdep, digit->Counts());
105     return kTRUE;
106   }
107   //__________________________________________________________________
108   Bool_t Finish()
109   {
110     gStyle->SetPalette(1);
111     gStyle->SetOptTitle(0);
112     gStyle->SetCanvasColor(0);
113     gStyle->SetCanvasBorderSize(0);
114     gStyle->SetPadColor(0);
115     gStyle->SetPadBorderSize(0);
116     fElossVsAdc->SetStats(kFALSE);
117     fElossVsAdc->Draw("COLZ");
118     return kTRUE;
119   }
120
121   ClassDef(DrawHitsDigits,0);
122 };
123
124 //____________________________________________________________________
125 //
126 // EOF
127 //