]>
Commit | Line | Data |
---|---|---|
088f8e79 | 1 | // |
2 | // $Id$ | |
3 | // | |
4 | // Script that contains a class to draw eloss from hits, versus ADC | |
5 | // counts from digits, using the AliFMDInputHits class in the util library. | |
6 | // | |
7 | // It draws the energy loss versus the p/(mq^2). It can be overlayed | |
8 | // with the Bethe-Bloc curve to show how the simulation behaves | |
9 | // relative to the expected. | |
10 | // | |
11 | // Use the script `Compile.C' to compile this class using ACLic. | |
12 | // | |
13 | #include <TH2D.h> | |
14 | #include <AliFMDHit.h> | |
15 | #include <AliFMDDigit.h> | |
16 | #include <AliFMDInput.h> | |
17 | #include <AliFMDEdepMap.h> | |
18 | #include <iostream> | |
19 | #include <TStyle.h> | |
20 | #include <TArrayF.h> | |
21 | ||
22 | class DrawHitsDigits : public AliFMDInputHits | |
23 | { | |
24 | private: | |
25 | TH2D* fElossVsAdc; // Histogram | |
26 | AliFMDEdepMap fMap; | |
27 | public: | |
28 | TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max) | |
29 | { | |
30 | TArrayF bins(n+1); | |
31 | Float_t dp = n / TMath::Log10(max / min); | |
32 | Float_t pmin = TMath::Log10(min); | |
33 | bins[0] = min; | |
34 | for (Int_t i = 1; i < n+1; i++) { | |
35 | Float_t p = pmin + i / dp; | |
36 | bins[i] = TMath::Power(10, p); | |
37 | } | |
38 | return bins; | |
39 | } | |
40 | DrawHitsDigits(Int_t n=900, Double_t emin=1e-3, Double_t emax=10, | |
41 | Int_t m=1100, Double_t amin=-0.5, Double_t amax=1099.5) | |
42 | { | |
43 | AddLoad(kDigits); | |
44 | TArrayF eloss(MakeLogScale(n, emin, emax)); | |
45 | TArrayF adcs(m+1); | |
46 | adcs[0] = amin; | |
47 | for (Int_t i = 1; i < m+1; i++) adcs[i] = adcs[i-1] + (amax-amin)/m; | |
48 | fElossVsAdc = new TH2D("bad", "#Delta E vs. ADC", | |
49 | eloss.fN-1, eloss.fArray, adcs.fN-1, adcs.fArray); | |
50 | fElossVsAdc->SetXTitle("#Delta E/#Delta x [MeV/cm]"); | |
51 | fElossVsAdc->SetYTitle("ADC value"); | |
52 | } | |
53 | Bool_t Begin(Int_t ev) | |
54 | { | |
55 | fMap.Reset(); | |
56 | return AliFMDInputHits::Begin(ev); | |
57 | } | |
58 | Bool_t ProcessHit(AliFMDHit* hit, TParticle*) | |
59 | { | |
60 | // Cache the energy loss | |
61 | if (!hit) return kFALSE; | |
62 | UShort_t det = hit->Detector(); | |
63 | Char_t rng = hit->Ring(); | |
64 | UShort_t sec = hit->Sector(); | |
65 | UShort_t str = hit->Strip(); | |
8f6ee336 | 66 | if (str > 511) { |
67 | AliWarning(Form("Bad strip number %d in hit", str)); | |
68 | return kTRUE; | |
69 | } | |
088f8e79 | 70 | fMap(det, rng, sec, str).fEdep += hit->Edep(); |
71 | fMap(det, rng, sec, str).fN++; | |
72 | return kTRUE; | |
73 | } | |
74 | Bool_t Event() | |
75 | { | |
76 | if (!AliFMDInputHits::Event()) return kFALSE; | |
77 | Int_t nEv = fTreeD->GetEntries(); | |
78 | for (Int_t i = 0; i < nEv; i++) { | |
79 | Int_t digitRead = fTreeD->GetEntry(i); | |
80 | if (digitRead <= 0) continue; | |
81 | Int_t nDigit = fArrayD->GetEntries(); | |
82 | if (nDigit <= 0) continue; | |
83 | for (Int_t j = 0; j < nDigit; j++) { | |
84 | AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j)); | |
85 | if (!digit) continue; | |
86 | UShort_t det = digit->Detector(); | |
87 | Char_t rng = digit->Ring(); | |
88 | UShort_t sec = digit->Sector(); | |
89 | UShort_t str = digit->Strip(); | |
8f6ee336 | 90 | if (str > 511) { |
91 | AliWarning(Form("Bad strip number %d in digit", str)); | |
92 | continue; | |
93 | } | |
088f8e79 | 94 | fElossVsAdc->Fill(fMap(det, rng, sec, str).fEdep, digit->Counts()); |
95 | } | |
96 | } | |
97 | return kTRUE; | |
98 | } | |
99 | Bool_t Finish() | |
100 | { | |
101 | gStyle->SetPalette(1); | |
102 | gStyle->SetOptTitle(0); | |
103 | gStyle->SetCanvasColor(0); | |
104 | gStyle->SetCanvasBorderSize(0); | |
105 | gStyle->SetPadColor(0); | |
106 | gStyle->SetPadBorderSize(0); | |
107 | fElossVsAdc->SetStats(kFALSE); | |
108 | fElossVsAdc->Draw("COLZ"); | |
109 | return kTRUE; | |
110 | } | |
8f6ee336 | 111 | |
112 | ClassDef(DrawHitsDigits,0); | |
088f8e79 | 113 | }; |
114 | ||
115 | //____________________________________________________________________ | |
116 | // | |
117 | // EOF | |
118 | // |