]>
Commit | Line | Data |
---|---|---|
d61f73d9 | 1 | Int_t AliTOFanalyzeSDigitsV2(Int_t ndump=0, Int_t iEvNum=-1, Int_t nEvent=0) |
ca4508a6 | 2 | { |
d61f73d9 | 3 | |
4 | ///////////////////////////////////////////////////////////////////////// | |
ca4508a6 | 5 | // |
6 | // Analyzes the TOF sdigits and fills QA-histograms | |
d61f73d9 | 7 | // iEvNum=-1 and nEvent=0 means all events in the file |
8 | // | |
9 | // root[0] .L AliTOFanalyzeSDigitsV2.C | |
10 | // root[1] AliTOFanalyzeSDigitsV2() | |
11 | // | |
12 | // If you want analyze only the sdigits in the 2th event | |
13 | // (existing int the header file), see in the following: | |
14 | // | |
15 | // root[0] .L AliTOFanalyzeSDigitsV2.C | |
16 | // root[1] AliTOFanalyzeSDigitsV2(0,2,1) | |
32bd3eb5 | 17 | // |
18 | // Updated to the new I/O by: A. De Caro, C. Zampolli | |
d61f73d9 | 19 | // |
20 | // Report problems to: decaro@sa.infn.it | |
21 | // | |
22 | ///////////////////////////////////////////////////////////////////////// | |
ca4508a6 | 23 | |
32bd3eb5 | 24 | Int_t rc = 0; |
ca4508a6 | 25 | |
26 | // adc and tdc | |
da3d3acd | 27 | TH1F *htdc = new TH1F("htdc","TDC [bin]",5000,0.,150000.); |
28 | TH1F *hadc = new TH1F("hadc","ADC [bin]",100,0., 3000.); | |
ca4508a6 | 29 | // TOF sdigit volumes |
da3d3acd | 30 | TH1F *hsector = new TH1F("hsector","Sector",18,0.,18.); |
31 | TH1F *hplate = new TH1F("hplate","Plate ", 5,0., 5.); | |
32 | TH1F *hstrip = new TH1F("hstrip","Strip ",20,0.,20.); | |
33 | TH1F *hpadz = new TH1F("hpadz","Pad along z ",2,0.,2.); | |
34 | TH1F *hpadx = new TH1F("hpadx","Pad along x",48,0.,48.); | |
ca4508a6 | 35 | // ADC-TDC correlation |
32bd3eb5 | 36 | TH2F *h2tdcVSadc = new TH2F("h2tdcVSadc","TDC [bin] VS ADC [bin]",500,0.,150000.,100,0.,3000.); |
37 | ||
38 | // Dynamically link some shared libs | |
39 | if (gClassTable->GetID("AliRun") < 0) { | |
40 | gROOT->LoadMacro("loadlibs.C"); | |
41 | loadlibs(); | |
42 | } | |
43 | ||
44 | if (gAlice) | |
45 | { | |
46 | delete gAlice->GetRunLoader(); | |
47 | delete gAlice; | |
48 | gAlice = 0x0; | |
49 | } | |
50 | ||
f5a857b2 | 51 | AliRunLoader *rl = AliRunLoader::Open("galice.root",AliConfig::GetDefaultEventFolderName(),"read"); |
32bd3eb5 | 52 | if (!rl) |
53 | { | |
54 | cerr<<"Can't load RunLoader from file"<<"!\n"; | |
55 | rc = 1; | |
56 | return rc; | |
57 | } | |
58 | ||
59 | rl->LoadgAlice(); | |
60 | gAlice = rl->GetAliRun(); | |
61 | ||
62 | if (!gAlice) | |
63 | { | |
64 | cerr << "<AliTOFanalyzeSDigits> AliRun object not found on file\n "; | |
65 | rc = 2; | |
66 | return rc; | |
67 | } | |
ca4508a6 | 68 | |
d61f73d9 | 69 | rl->LoadHeader(); |
70 | ||
32bd3eb5 | 71 | AliLoader *tofl = rl->GetLoader("TOFLoader"); |
72 | AliTOF *tof = (AliTOF *) rl->GetAliRun()->GetDetector("TOF"); | |
73 | ||
74 | if (tof==0x0 || tofl==0x0) | |
75 | { | |
76 | cerr << "<AliTOFanalyzeSDigits> no TOF detector found" << endl; | |
77 | rc = 3; | |
78 | return rc; | |
79 | } | |
80 | ||
d61f73d9 | 81 | cout << "First " << ndump << " SDigits found in TOF TreeS branch have: \n"; |
ca4508a6 | 82 | |
d61f73d9 | 83 | Int_t upperLimit; |
84 | Int_t bottomLimit; | |
85 | ||
86 | if (iEvNum<0) bottomLimit=0; | |
87 | else bottomLimit = iEvNum; | |
88 | ||
89 | if (nEvent == 0) upperLimit = (Int_t)(rl->GetNumberOfEvents()); | |
90 | else upperLimit = nEvent+bottomLimit; | |
ca4508a6 | 91 | |
32bd3eb5 | 92 | AliTOFSDigit *tofsdigit; |
93 | ||
d61f73d9 | 94 | for (Int_t ievent = bottomLimit; ievent < upperLimit; ievent++) { |
95 | ||
32bd3eb5 | 96 | rl->GetEvent(ievent); |
d61f73d9 | 97 | printf ("Processing event %d \n", ievent); |
ca4508a6 | 98 | |
32bd3eb5 | 99 | // Get the pointer SDigit tree |
d61f73d9 | 100 | tofl->LoadSDigits("read"); |
32bd3eb5 | 101 | TTree *TS=tofl->TreeS(); |
102 | tof->SetTreeAddress(); | |
ca4508a6 | 103 | |
32bd3eb5 | 104 | if(!TS) |
105 | { | |
d61f73d9 | 106 | cout << "<AliTOFanalyzeSDigits> No TreeS found \n"; |
32bd3eb5 | 107 | rc = 4; |
108 | return rc; | |
109 | } | |
ca4508a6 | 110 | |
32bd3eb5 | 111 | TClonesArray * TOFsdigits = new TClonesArray("AliTOFSDigit",1000); |
d61f73d9 | 112 | //TOFsdigits = tof->SDigits(); |
113 | tof->SDigits(); | |
32bd3eb5 | 114 | TOFsdigits = TS->GetBranch("TOF")->SetAddress(&TOFsdigits); |
115 | ||
116 | Int_t nEntries = TS->GetEntries(); | |
117 | ||
118 | for (Int_t iEntry = 0; iEntry < nEntries; iEntry ++) | |
119 | { | |
120 | tof->ResetDigits(); | |
121 | TS->GetEvent(iEntry); | |
122 | Int_t ndig = TOFsdigits->GetEntriesFast(); | |
123 | cout << "<AliTOFanalyzeSDigits> found " << ndig | |
124 | << " TOF sdigits for event " << ievent << endl; | |
ca4508a6 | 125 | |
32bd3eb5 | 126 | for (Int_t k=0; k<ndig; k++) { |
127 | tofsdigit = (AliTOFSDigit*) TOFsdigits->UncheckedAt(k); | |
128 | Float_t firstTDC = tofsdigit->GetTdc(0); | |
129 | Float_t firstADC = tofsdigit->GetAdc(0); | |
130 | htdc->Fill(firstTDC); | |
131 | hadc->Fill(firstADC); | |
132 | // TOF sdigit volumes | |
da3d3acd | 133 | Int_t sector = tofsdigit->GetSector(); // range [0-17] |
134 | Int_t plate = tofsdigit->GetPlate(); // range [0- 4] | |
135 | Int_t strip = tofsdigit->GetStrip(); // range [0-19] | |
136 | Int_t padz = tofsdigit->GetPadz(); // range [0- 1] | |
137 | Int_t padx = tofsdigit->GetPadx(); // range [0-47] | |
32bd3eb5 | 138 | // it is QA, then I perform QA! |
da3d3acd | 139 | Bool_t isSDigitBad = (sector<0 || sector>17 || plate<0 || plate >4 || padz<0 || padz>1 || padx<0 || padx>47); |
32bd3eb5 | 140 | |
da3d3acd | 141 | |
32bd3eb5 | 142 | if (isSDigitBad) { |
d61f73d9 | 143 | cout << "<AliTOFanalyzeSDigits> strange sdigit found \n"; |
32bd3eb5 | 144 | rc = 4; |
145 | return rc; | |
146 | } | |
147 | ||
148 | if(k<ndump){ | |
d61f73d9 | 149 | cout << k << "-th | Sector " << sector << " | Plate " << plate << " | Strip " << strip << " | PadZ " << padz << " | PadX " << padx << endl; |
150 | cout << k << "-th | ADC " << firstADC << " [bin] | TDC " << firstTDC << " [bin] \n"; | |
151 | cout << "---------------------------------------------------- \n"; | |
32bd3eb5 | 152 | } |
153 | ||
154 | // filling sdigit volume histos | |
155 | hsector->Fill(sector); | |
156 | hplate->Fill(plate); | |
157 | hstrip->Fill(strip); | |
158 | hpadx->Fill(padx); | |
159 | hpadz->Fill(padz); | |
160 | h2tdcVSadc->Fill(firstTDC,firstADC); | |
161 | ||
ca4508a6 | 162 | } |
32bd3eb5 | 163 | } |
ca4508a6 | 164 | |
32bd3eb5 | 165 | tofl->UnloadSDigits(); |
ca4508a6 | 166 | |
32bd3eb5 | 167 | } // end loop on events |
168 | ||
169 | rl->UnloadHeader(); | |
170 | rl->UnloadgAlice(); | |
171 | ||
ca4508a6 | 172 | TFile *fout = new TFile("TOF_sdigitsQA.root","RECREATE"); |
173 | htdc->Write(); | |
174 | hadc->Write(); | |
ca4508a6 | 175 | h2tdcVSadc->Write(); |
176 | hsector->Write(); | |
177 | hplate->Write(); | |
178 | hstrip->Write(); | |
179 | hpadz->Write(); | |
180 | hpadx->Write(); | |
32bd3eb5 | 181 | fout->Close(); |
182 | ||
183 | delete htdc; | |
184 | delete hadc; | |
185 | delete h2tdcVSadc; | |
186 | delete hsector; | |
187 | delete hplate; | |
188 | delete hstrip; | |
189 | delete hpadz; | |
190 | delete hpadx; | |
ca4508a6 | 191 | |
32bd3eb5 | 192 | if (gAlice) |
193 | { | |
194 | delete gAlice->GetRunLoader(); | |
195 | delete gAlice; | |
196 | gAlice = 0x0; | |
197 | } | |
ca4508a6 | 198 | |
199 | return rc; | |
200 | ||
ca4508a6 | 201 | } |
32bd3eb5 | 202 |