]>
Commit | Line | Data |
---|---|---|
faf80567 | 1 | //____________________________________________________________________ |
2 | // | |
3 | // $Id: DrawBothDigits.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 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 <TH2F.h> | |
15 | #include <AliFMDDigit.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 | ||
24 | /** @class DrawBothDigits | |
25 | @brief Draw hit energy loss versus digit ADC | |
26 | @code | |
27 | Root> .L Compile.C | |
28 | Root> Compile("DrawBothDigits.C") | |
29 | Root> DrawBothDigits c | |
30 | Root> c.Run(); | |
31 | @endcode | |
32 | @ingroup FMD_script | |
33 | */ | |
34 | class DrawBothDigits : public AliFMDInput | |
35 | { | |
36 | private: | |
37 | TH2F* fTrackNos; // Histogram | |
38 | AliFMDEdepMap fCache; | |
39 | public: | |
40 | //__________________________________________________________________ | |
41 | DrawBothDigits(Int_t max=300) | |
42 | : AliFMDInput("galice.root") | |
43 | { | |
44 | AddLoad(kDigits); | |
45 | AddLoad(kSDigits); | |
46 | fTrackNos = new TH2F("trackNos", "Track numbers", | |
47 | max+1, -1.5, max-.5, max+1, -1.5, max-.5); | |
48 | fTrackNos->SetXTitle("Digit track"); | |
49 | fTrackNos->SetYTitle("SDigit track"); | |
50 | } | |
51 | //__________________________________________________________________ | |
52 | Bool_t Begin(Int_t evno) | |
53 | { | |
54 | fCache.Reset(); | |
55 | return AliFMDInput::Begin(evno); | |
56 | } | |
57 | //__________________________________________________________________ | |
58 | Bool_t ProcessSDigit(AliFMDSDigit* sdigit) | |
59 | { | |
60 | if (!sdigit) return kTRUE; | |
61 | AliFMDEdepHitPair& entry = fCache(sdigit->Detector(), | |
62 | sdigit->Ring(), | |
63 | sdigit->Sector(), | |
64 | sdigit->Strip()); | |
65 | entry.fLabels.Set(sdigit->GetNTrack()); | |
66 | Info("ProcessSDigit", "Got %d SDigit tracks", sdigit->GetNTrack()); | |
67 | for (size_t i = 0; i < sdigit->GetNTrack(); i++) | |
68 | entry.fLabels.fArray[i] = sdigit->GetTrack(i); | |
69 | return kTRUE; | |
70 | } | |
71 | //__________________________________________________________________ | |
72 | Bool_t ProcessDigit(AliFMDDigit* digit) | |
73 | { | |
74 | if (!digit) return kTRUE; | |
75 | AliFMDEdepHitPair& entry = fCache(digit->Detector(), | |
76 | digit->Ring(), | |
77 | digit->Sector(), | |
78 | digit->Strip()); | |
79 | TArrayI& stracks = entry.fLabels; | |
80 | Info("ProcessDigit", "Got %d SDigit tracks, and %d Digit tracks", | |
81 | stracks.fN, digit->GetNTrack()); | |
82 | for (Int_t i = 0; (i < stracks.fN || i < digit->GetNTrack()); i++) { | |
83 | Int_t strack = (i < stracks.fN ? stracks.fArray[i] : -1); | |
84 | Int_t dtrack = (i < digit->GetNTrack() ? digit->GetTrack(i) : -1); | |
85 | fTrackNos->Fill(strack, dtrack); | |
86 | } | |
87 | return kTRUE; | |
88 | } | |
89 | //__________________________________________________________________ | |
90 | Bool_t Finish() | |
91 | { | |
92 | gStyle->SetPalette(1); | |
93 | gStyle->SetOptTitle(0); | |
94 | gStyle->SetCanvasColor(0); | |
95 | gStyle->SetCanvasBorderSize(0); | |
96 | gStyle->SetPadColor(0); | |
97 | gStyle->SetPadBorderSize(0); | |
98 | fTrackNos->SetStats(kFALSE); | |
99 | // fTrackNos->Scale(1. / fAdc->GetEntries()); | |
100 | fTrackNos->Draw("colz"); | |
101 | return kTRUE; | |
102 | } | |
103 | ||
104 | ClassDef(DrawBothDigits,0); | |
105 | }; | |
106 | ||
107 | //____________________________________________________________________ | |
108 | // | |
109 | // EOF | |
110 | // |