]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ZDC/AliZDCQADataMakerSim.cxx
Deleted erroneously added simlinks
[u/mrichter/AliRoot.git] / ZDC / AliZDCQADataMakerSim.cxx
1 /**************************************************************************\r
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3  *                                                                        *\r
4  * Author: The ALICE Off-line Project.                                    *\r
5  * Contributors are mentioned in the code where appropriate.              *\r
6  *                                                                        *\r
7  * Permission to use, copy, modify and distribute this software and its   *\r
8  * documentation strictly for non-commercial purposes is hereby granted   *\r
9  * without fee, provided that the above copyright notice appears in all   *\r
10  * copies and that both the copyright notice and this permission notice   *\r
11  * appear in the supporting documentation. The authors make no claims     *\r
12  * about the suitability of this software for any purpose. It is          *\r
13  * provided "as is" without express or implied warranty.                  *\r
14  **************************************************************************/\r
15 \r
16 // --- Standard library ---\r
17 #include <Riostream.h>\r
18 // --- ROOT system ---\r
19 #include <TClonesArray.h>\r
20 #include <TFile.h>     \r
21 #include <TH1F.h> \r
22 #include <TH2F.h>\r
23 #include <TBranch.h>\r
24 #include <TTree.h>\r
25 // --- AliRoot header files ---\r
26 #include "AliLog.h"\r
27 #include "AliQAChecker.h"\r
28 #include "AliZDCQADataMakerSim.h"\r
29 #include "AliZDCHit.h"\r
30 #include "AliZDCDigit.h"\r
31 #include "AliZDCRawStream.h"\r
32 \r
33 ClassImp(AliZDCQADataMakerSim)\r
34            \r
35 //____________________________________________________________________________ \r
36   AliZDCQADataMakerSim::AliZDCQADataMakerSim() : \r
37       AliQADataMakerSim(AliQAv1::GetDetName(AliQAv1::kZDC), "ZDC Quality Assurance Data Maker"),\r
38       fDigit(0)\r
39 {\r
40   // ctor\r
41 }\r
42 \r
43 //____________________________________________________________________________ \r
44 AliZDCQADataMakerSim::AliZDCQADataMakerSim(const AliZDCQADataMakerSim& qadm) :\r
45     AliQADataMakerSim(), \r
46     fDigit(0) \r
47 {\r
48   //copy ctor \r
49   SetName((const char*)qadm.GetName()); \r
50   SetTitle((const char*)qadm.GetTitle()); \r
51 }\r
52 \r
53 //__________________________________________________________________\r
54 AliZDCQADataMakerSim& AliZDCQADataMakerSim::operator = (const AliZDCQADataMakerSim& qadm )\r
55 {\r
56   // Equal operator.\r
57   this->~AliZDCQADataMakerSim();\r
58   new(this) AliZDCQADataMakerSim(qadm);\r
59   return *this;\r
60 }\r
61  \r
62 //____________________________________________________________________________ \r
63 void AliZDCQADataMakerSim::InitHits()\r
64 {\r
65   // create Hits histograms in Hits subdir\r
66   //\r
67   const Bool_t expert   = kTRUE ; \r
68   const Bool_t image    = kTRUE ; \r
69 \r
70   TH1F * hHitsZNCTot = new TH1F("hHitsZNCTot", "Signal in ZNC; N_{phe}", 100, 0., 6000.);\r
71   TH1F * hHitsZNATot = new TH1F("hHitsZNATot", "Signal in ZNA; N_{phe}", 100, 0., 6000.);\r
72   TH1F * hHitsZPCTot = new TH1F("hHitsZPCTot", "Signal in ZPC; N_{phe}", 100, 0., 6000.);\r
73   TH1F * hHitsZPATot = new TH1F("hHitsZPATot", "Signal in ZPA; N_{phe}", 100, 0., 6000.);\r
74   Add2HitsList(hHitsZNCTot, 0, !expert, image);\r
75   Add2HitsList(hHitsZNATot, 1, !expert, image);\r
76   Add2HitsList(hHitsZPCTot, 2, !expert, image);\r
77   Add2HitsList(hHitsZPATot, 3, !expert, image);\r
78   //\r
79   TH1F * hHitsSumQZNC = new TH1F("hHitsSumQZNC", "Signal in 4 ZNC PMQ; N_{phe}",100, 0., 4000.);\r
80   TH1F * hHitsSumQZNA = new TH1F("hHitsSumQZNA", "Signal in 4 ZNA PMQ; N_{phe}",100, 0., 4000.);\r
81   TH1F * hHitsSumQZPC = new TH1F("hHitsSumQZPC", "Signal in 4 ZPC PMQ; N_{phe}",100, 0., 4000.);\r
82   TH1F * hHitsSumQZPA = new TH1F("hHitsSumQZPA", "Signal in 4 ZPA PMQ; N_{phe}",100, 0., 4000.);\r
83   Add2HitsList(hHitsSumQZNC, 4, expert, !image);\r
84   Add2HitsList(hHitsSumQZNA, 5, expert, !image);\r
85   Add2HitsList(hHitsSumQZPC, 6, expert, !image);\r
86   Add2HitsList(hHitsSumQZPA, 7, expert, !image);\r
87   //\r
88   TH1F * hHitsPMCZNC = new TH1F("hHitsPMCZNC", "Signal in ZNC PMC; N_{phe}",100, 0., 4000.);\r
89   TH1F * hHitsPMCZNA = new TH1F("hHitsPMCZNA", "Signal in ZNA PMC; N_{phe}",100, 0., 4000.);\r
90   TH1F * hHitsPMCZPC = new TH1F("hHitsPMCZPC", "Signal in ZPC PMC; N_{phe}",100, 0., 4000.);\r
91   TH1F * hHitsPMCZPA = new TH1F("hHitsPMCZPA", "Signal in ZPA PMC; N_{phe}",100, 0., 4000.);\r
92   Add2HitsList(hHitsPMCZNC, 8, expert, !image);\r
93   Add2HitsList(hHitsPMCZNA, 9, expert, !image);\r
94   Add2HitsList(hHitsPMCZPC, 10, expert, !image);\r
95   Add2HitsList(hHitsPMCZPA, 11, expert, !image);\r
96   \r
97   TH2F * hHitsZNCh  = new TH2F("hHitsZNCh", "Hits centroid in ZNC;Centroid position [cm];Counts", 100, -5.,5.,100,-5.,5.);\r
98   TH2F * hHitsZNAh  = new TH2F("hHitsZNAh", "Hits centroid in ZNA;Centroid position [cm];Counts", 100, -5.,5.,100,-5.,5.);\r
99   Add2HitsList(hHitsZNCh, 12, !expert, image);\r
100   Add2HitsList(hHitsZNAh, 13, !expert, image);\r
101   // NB -> For the moment no check is performesd on ZP centroids\r
102   TH2F * hHitsZPCh  = new TH2F("hHitsZPCh", "Hits centroid in ZPC;Centroid position [cm];Counts", 100,-12.,12.,100,-12.,12.); \r
103   TH2F * hHitsZPAh  = new TH2F("hHitsZPAh", "Hits centroid in ZPA;Centroid position [cm];Counts", 100,-12.,12.,100,-12.,12.); \r
104   Add2HitsList(hHitsZPCh, 14, !expert, image);\r
105   Add2HitsList(hHitsZPAh, 15, !expert, image);\r
106 }\r
107 \r
108 \r
109 //____________________________________________________________________________ \r
110 void AliZDCQADataMakerSim::InitDigits()\r
111 {\r
112   // create Digits histograms in Digits subdir\r
113   //\r
114   const Bool_t expert   = kTRUE ; \r
115   const Bool_t image    = kTRUE ; \r
116   \r
117   // ------------------- HIGH GAIN CHAIN ---------------------------\r
118   TH1F * hDigZNCTot = new TH1F("hDigZNCTot", "Signal in ZNC;Amplitude [ADC counts];Counts", 100, 0., 6000.);\r
119   TH1F * hDigZNATot = new TH1F("hDigZNATot", "Signal in ZNA;Amplitude [ADC counts];Counts", 100, 0., 6000.);\r
120   TH1F * hDigZPCTot = new TH1F("hDigZPCTot", "Signal in ZPC;Amplitude [ADC counts];Counts", 100, 0., 6000.);\r
121   TH1F * hDigZPATot = new TH1F("hDigZPATot", "Signal in ZPA;Amplitude [ADC counts];Counts", 100, 0., 6000.);\r
122   Add2DigitsList(hDigZNCTot, 0, !expert, image);\r
123   Add2DigitsList(hDigZNATot, 1, !expert, image);\r
124   Add2DigitsList(hDigZPCTot, 2, !expert, image);\r
125   Add2DigitsList(hDigZPATot, 3, !expert, image);\r
126   //\r
127   TH1F * hDigSumQZNC = new TH1F("hDigSumQZNC", "Signal in 4 ZNC PMQ;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
128   TH1F * hDigSumQZNA = new TH1F("hDigSumQZNA", "Signal in 4 ZNA PMQ;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
129   TH1F * hDigSumQZPC = new TH1F("hDigSumQZPC", "Signal in 4 ZPC PMQ;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
130   TH1F * hDigSumQZPA = new TH1F("hDigSumQZPA", "Signal in 4 ZPA PMQ;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
131   Add2DigitsList(hDigSumQZNC, 4, expert, !image);\r
132   Add2DigitsList(hDigSumQZNA, 5, expert, !image);\r
133   Add2DigitsList(hDigSumQZPC, 6, expert, !image);\r
134   Add2DigitsList(hDigSumQZPA, 7, expert, !image);\r
135   //\r
136   TH1F * hDigPMCZNC = new TH1F("hDigPMCZNC", "Signal in ZNC PMC;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
137   TH1F * hDigPMCZNA = new TH1F("hDigPMCZNA", "Signal in ZNA PMC;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
138   TH1F * hDigPMCZPC = new TH1F("hDigPMCZPC", "Signal in ZPC PMC;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
139   TH1F * hDigPMCZPA = new TH1F("hDigPMCZPA", "Signal in ZPA PMC;Amplitude [ADC counts];Counts",100, 0., 4000.);\r
140   Add2DigitsList(hDigPMCZNC, 8, expert, !image);\r
141   Add2DigitsList(hDigPMCZNA, 9, expert, !image);\r
142   Add2DigitsList(hDigPMCZPC, 10, expert, !image);\r
143   Add2DigitsList(hDigPMCZPA, 11, expert, !image);\r
144   // \r
145   // ------------------- LOW GAIN CHAIN ---------------------------\r
146   TH1F * hDigZNCTotlg = new TH1F("hDigZNCTotlg", "Digit lg signal in ZNC", 100, 0., 6000.);\r
147   TH1F * hDigZNATotlg = new TH1F("hDigZNATotlg", "Digit lg signal in ZNA", 100, 0., 6000.);\r
148   TH1F * hDigZPCTotlg = new TH1F("hDigZPCTotlg", "Digit lg signal in ZPC", 100, 0., 6000.);\r
149   TH1F * hDigZPATotlg = new TH1F("hDigZPATotlg", "Digit lg signal in ZPA", 100, 0., 6000.);\r
150   Add2DigitsList(hDigZNCTotlg, 12, !expert, image);\r
151   Add2DigitsList(hDigZNATotlg, 13, !expert, image);\r
152   Add2DigitsList(hDigZPCTotlg, 14, !expert, image);\r
153   Add2DigitsList(hDigZPATotlg, 15, !expert, image);\r
154   //\r
155   TH1F * hDigSumQZNClg = new TH1F("hDigSumQZNClg", "Signal in 4 ZNC PMQlg",100, 0., 4000.);\r
156   TH1F * hDigSumQZNAlg = new TH1F("hDigSumQZNAlg", "Signal in 4 ZNA PMQlg",100, 0., 4000.);\r
157   TH1F * hDigSumQZPClg = new TH1F("hDigSumQZPClg", "Signal in 4 ZPC PMQlg",100, 0., 4000.);\r
158   TH1F * hDigSumQZPAlg = new TH1F("hDigSumQZPAlg", "Signal in 4 ZPA PMQlg",100, 0., 4000.);\r
159   Add2DigitsList(hDigSumQZNClg, 16, expert, !image);\r
160   Add2DigitsList(hDigSumQZNAlg, 17, expert, !image);\r
161   Add2DigitsList(hDigSumQZPClg, 18, expert, !image);\r
162   Add2DigitsList(hDigSumQZPAlg, 19, expert, !image);\r
163   //\r
164   TH1F * hDigPMCZNClg = new TH1F("hDigPMCZNClg", "Signal in ZNC PMClg",100, 0., 4000.);\r
165   TH1F * hDigPMCZNAlg = new TH1F("hDigPMCZNAlg", "Signal in ZNA PMClg",100, 0., 4000.);\r
166   TH1F * hDigPMCZPClg = new TH1F("hDigPMCZPClg", "Signal in ZPC PMClg",100, 0., 4000.);\r
167   TH1F * hDigPMCZPAlg = new TH1F("hDigPMCZPAlg", "Signal in ZPA PMClg",100, 0., 4000.);\r
168   Add2DigitsList(hDigPMCZNClg, 20, expert, !image);\r
169   Add2DigitsList(hDigPMCZNAlg, 21, expert, !image);\r
170   Add2DigitsList(hDigPMCZPClg, 22, expert, !image);\r
171   Add2DigitsList(hDigPMCZPAlg, 23, expert, !image);\r
172 \r
173 }\r
174 \r
175 //____________________________________________________________________________\r
176 void AliZDCQADataMakerSim::MakeHits()\r
177 {\r
178   //filling QA histos for Hits\r
179   //\r
180 \r
181   // Check id histograms already created for this Event Specie\r
182   if ( ! GetHitsData(0) )\r
183     InitHits() ;\r
184   \r
185   TIter next(fHitsArray); \r
186   AliZDCHit * hit; \r
187   Float_t adcSum_ZNC=0., adcSum_ZNA=0., adcSum_ZPC=0., adcSum_ZPA=0.;\r
188   Float_t adcSumQ_ZNC=0., adcSumQ_ZNA=0., adcSumQ_ZPC=0., adcSumQ_ZPA=0.;\r
189   while((hit = dynamic_cast<AliZDCHit *>(next()))){\r
190     if(hit->GetVolume(0)==1){\r
191        adcSumQ_ZNC += hit->GetLightPMQ();\r
192        adcSum_ZNC  += hit->GetLightPMC() + hit->GetLightPMQ();\r
193        //\r
194        GetHitsData(8)->Fill(hit->GetLightPMQ());\r
195        //\r
196        GetHitsData(12)->Fill(hit->GetXImpact(),hit->GetYImpact());        \r
197     }\r
198     else if(hit->GetVolume(0)==4){\r
199        adcSumQ_ZNA += hit->GetLightPMQ();\r
200        adcSum_ZNA  += hit->GetLightPMC() + hit->GetLightPMQ();\r
201        //\r
202        GetHitsData(9)->Fill(hit->GetLightPMQ());\r
203        //\r
204        GetHitsData(13)->Fill(hit->GetXImpact(), hit->GetYImpact());\r
205     }\r
206     else if(hit->GetVolume(0)==2){\r
207        adcSumQ_ZNC += hit->GetLightPMQ();\r
208        adcSum_ZNC  += hit->GetLightPMC() + hit->GetLightPMQ();\r
209        //\r
210        GetHitsData(10)->Fill(hit->GetLightPMQ());\r
211        //\r
212        GetHitsData(14)->Fill(hit->GetXImpact(), hit->GetYImpact());\r
213     }\r
214     else if(hit->GetVolume(0)==5){\r
215        adcSumQ_ZNC += hit->GetLightPMQ();\r
216        adcSum_ZNC  += hit->GetLightPMC() + hit->GetLightPMQ();\r
217        //\r
218        GetHitsData(11)->Fill(hit->GetLightPMQ());\r
219        //\r
220        GetHitsData(15)->Fill(hit->GetXImpact(), hit->GetYImpact());\r
221     }\r
222     //\r
223     GetHitsData(0)->Fill(adcSum_ZNC);\r
224     GetHitsData(1)->Fill(adcSum_ZNA);\r
225     GetHitsData(2)->Fill(adcSum_ZPC);\r
226     GetHitsData(3)->Fill(adcSum_ZPA);\r
227     //\r
228     GetHitsData(4)->Fill(adcSumQ_ZNC);\r
229     GetHitsData(5)->Fill(adcSumQ_ZNA);\r
230     GetHitsData(6)->Fill(adcSumQ_ZPC);\r
231     GetHitsData(7)->Fill(adcSumQ_ZPA);\r
232   }\r
233 }\r
234 \r
235 //___________________________________________________________________________\r
236 void AliZDCQADataMakerSim::MakeHits(TTree * hitTree)\r
237 {\r
238   // make QA data from Hit Tree\r
239   //\r
240   if(!hitTree){\r
241     AliError("Hit Tree not found!"); \r
242     return;\r
243   }\r
244   //\r
245 \r
246   TBranch * branch = hitTree->GetBranch("ZDC") ;\r
247 \r
248   if(!branch){\r
249     AliError("ZDC branch in Hit Tree not found!"); \r
250     return;\r
251   } \r
252   else{\r
253     if (fHitsArray) \r
254       fHitsArray->Clear() ;                    \r
255     char** add = (char**) (branch->GetAddress());\r
256     if(add){\r
257         fHitsArray = (TClonesArray*)(*add);\r
258     } \r
259     else{\r
260         if(!fHitsArray) fHitsArray = new TClonesArray("AliZDCHit", 1000);\r
261         branch->SetAddress(&fHitsArray);\r
262     }\r
263     Int_t ntracks = (Int_t) hitTree->GetEntries();\r
264     if (ntracks<=0) return;\r
265     //\r
266     for(Int_t itrack=0; itrack<ntracks; itrack++){\r
267         \r
268         branch->GetEntry(itrack);\r
269         //\r
270         MakeHits(); \r
271         fHitsArray->Clear();\r
272     }   \r
273   }\r
274 }\r
275 \r
276 //___________________________________________________________________________\r
277 void AliZDCQADataMakerSim::MakeDigits(TTree *digitTree )\r
278 {\r
279   // makes data from Digit Tree\r
280   TBranch * branch = digitTree->GetBranch("ZDC");\r
281   if(!branch){\r
282     AliError("ZDC branch in Digit Tree not found"); \r
283     return;\r
284   } \r
285   \r
286   // Check id histograms already created for this Event Specie\r
287   if ( ! GetDigitsData(0) )\r
288     InitDigits() ;\r
289   \r
290   branch->SetAddress(&fDigit);\r
291   \r
292   Int_t ndig = digitTree->GetEntries();\r
293    \r
294   Float_t adcSum_ZNC=0., adcSum_ZNA=0., adcSum_ZPC=0., adcSum_ZPA=0.;\r
295   Float_t adcSumQ_ZNC=0., adcSumQ_ZNA=0., adcSumQ_ZPC=0., adcSumQ_ZPA=0.;\r
296   Float_t adcSum_ZNC_lg=0., adcSum_ZNA_lg=0., adcSum_ZPC_lg=0., adcSum_ZPA_lg=0.;\r
297   Float_t adcSumQ_ZNC_lg=0., adcSumQ_ZNA_lg=0., adcSumQ_ZPC_lg=0., adcSumQ_ZPA_lg=0.;\r
298   //\r
299   for(Int_t i = 0; i < ndig; i++){\r
300       digitTree->GetEntry(i);\r
301       if(fDigit->GetSector(0)==1){\r
302           adcSum_ZNC += fDigit->GetADCValue(0);\r
303           adcSum_ZNC_lg += fDigit->GetADCValue(1);\r
304           //\r
305           if(fDigit->GetSector(1)!=0){\r
306               adcSumQ_ZNC += fDigit->GetADCValue(0);\r
307               adcSumQ_ZNC_lg+= fDigit->GetADCValue(1);\r
308           }\r
309           else{\r
310               GetDigitsData(8)->Fill(fDigit->GetADCValue(0));\r
311               GetDigitsData(20)->Fill(fDigit->GetADCValue(1));\r
312           }\r
313       }\r
314       else if(fDigit->GetSector(0)==2){\r
315           adcSum_ZPC += fDigit->GetADCValue(0);\r
316           adcSum_ZPC_lg += fDigit->GetADCValue(1);\r
317           //\r
318           if(fDigit->GetSector(1)!=0){\r
319               adcSumQ_ZPC += fDigit->GetADCValue(0);\r
320               adcSumQ_ZPC_lg+= fDigit->GetADCValue(1);\r
321           }\r
322           else{\r
323               GetDigitsData(10)->Fill(fDigit->GetADCValue(0));\r
324               GetDigitsData(22)->Fill(fDigit->GetADCValue(1));\r
325           }\r
326       }\r
327       else if(fDigit->GetSector(0)==4){\r
328           adcSum_ZNA += fDigit->GetADCValue(0);\r
329           adcSum_ZNA_lg += fDigit->GetADCValue(1);\r
330           //\r
331           if(fDigit->GetSector(1)!=0){\r
332               adcSumQ_ZNA += fDigit->GetADCValue(0);\r
333               adcSumQ_ZNA_lg+= fDigit->GetADCValue(1);\r
334           }\r
335           else{\r
336               GetDigitsData(9)->Fill(fDigit->GetADCValue(0));\r
337               GetDigitsData(21)->Fill(fDigit->GetADCValue(1));\r
338           }\r
339       }\r
340       else if(fDigit->GetSector(0)==5){\r
341           adcSum_ZPA += fDigit->GetADCValue(0);\r
342           adcSum_ZPA_lg += fDigit->GetADCValue(1);\r
343           //\r
344           if(fDigit->GetSector(1)!=0){\r
345               adcSumQ_ZPA += fDigit->GetADCValue(0);\r
346               adcSumQ_ZPA_lg+= fDigit->GetADCValue(1);\r
347           }\r
348           else{\r
349               GetDigitsData(11)->Fill(fDigit->GetADCValue(0));\r
350               GetDigitsData(23)->Fill(fDigit->GetADCValue(1));\r
351           }\r
352       }\r
353   }\r
354   //\r
355   GetDigitsData(0)->Fill(adcSum_ZNC);\r
356   GetDigitsData(1)->Fill(adcSum_ZNA);\r
357   GetDigitsData(2)->Fill(adcSum_ZPC);\r
358   GetDigitsData(3)->Fill(adcSum_ZPA);\r
359   //\r
360   GetDigitsData(4)->Fill(adcSumQ_ZNC);\r
361   GetDigitsData(5)->Fill(adcSumQ_ZNA);\r
362   GetDigitsData(6)->Fill(adcSumQ_ZPC);\r
363   GetDigitsData(7)->Fill(adcSumQ_ZPA);\r
364   //\r
365   GetDigitsData(12)->Fill(adcSum_ZNC_lg);\r
366   GetDigitsData(13)->Fill(adcSum_ZNA_lg);\r
367   GetDigitsData(14)->Fill(adcSum_ZPC_lg);\r
368   GetDigitsData(15)->Fill(adcSum_ZPA_lg);\r
369   //\r
370   GetDigitsData(16)->Fill(adcSumQ_ZNC_lg);\r
371   GetDigitsData(17)->Fill(adcSumQ_ZNA_lg);\r
372   GetDigitsData(18)->Fill(adcSumQ_ZPC_lg);\r
373   GetDigitsData(19)->Fill(adcSumQ_ZPA_lg);\r
374 }\r
375 \r
376 //____________________________________________________________________________\r
377 void AliZDCQADataMakerSim::StartOfDetectorCycle()\r
378 {\r
379   //Detector specific actions at start of cycle\r
380   \r
381 }\r
382 \r
383 //____________________________________________________________________________ \r
384 void AliZDCQADataMakerSim::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)\r
385 {\r
386   // Detector specific actions at end of cycle\r
387   // do the QA checking\r
388   AliQAChecker::Instance()->Run(AliQAv1::kZDC, task, list);  \r
389 }\r