]> git.uio.no Git - u/mrichter/AliRoot.git/blob - FMD/AliFMDPedestalDA.cxx
Starting a collection of QA/Comparison macros
[u/mrichter/AliRoot.git] / FMD / AliFMDPedestalDA.cxx
1 /**************************************************************************
2  * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /** @file    AliFMDPedestalDA.cxx
17     @author  Hans Hjersing Dalsgaard <canute@nbi.dk>
18     @date    Mon Mar 10 09:46:05 2008
19     @brief   Derived class for the pedestal detector algorithm.
20 */
21 //
22 // This class implements the virtual functions of the AliFMDBaseDA
23 // class.  The most important of these functions, FillChannels(..) and
24 // Analyse(...) collect and analyse the data of each channel. The
25 // resulting pedestal and noise values are written to a comma
26 // separated values (csv) file on the go. The csv files produced in
27 // this way are the basic input to the AliFMDPreprocessor.
28 //
29
30 #include "AliFMDPedestalDA.h"
31 #include "AliFMDAltroMapping.h"
32 #include "AliFMDParameters.h"
33 #include "AliFMDCalibPedestal.h"
34 #include "AliFMDDigit.h"
35 #include "AliLog.h"
36 #include <iostream>
37 #include <fstream>
38 #include <iomanip>
39 #include <TFile.h>
40 #include <TF1.h>
41 #include <TObject.h>
42 #include <TMath.h>
43 #include <TSystem.h>
44 #include <TDatime.h>
45 #include <TH2.h>
46
47 //_____________________________________________________________________
48 ClassImp(AliFMDPedestalDA)
49
50 //_____________________________________________________________________
51 AliFMDPedestalDA::AliFMDPedestalDA() 
52   : AliFMDBaseDA(),
53     fCurrentChannel(1),
54     fPedSummary("PedestalSummary","pedestals",51200,0,51200),
55     fNoiseSummary("NoiseSummary","noise",51200,0,51200),
56     fZSfileFMD1(),
57     fZSfileFMD2(),
58     fZSfileFMD3(), 
59     fMinTimebin(3 * 4 * 3 * 16), // 3 ddls, 4 FECs, 3 Altros, 16 channels
60     fMaxTimebin(3 * 4 * 3 * 16), // 3 ddls, 4 FECs, 3 Altros, 16 channels
61     fSummaryFMD1i(0),
62     fSummaryFMD2i(0),
63     fSummaryFMD2o(0),
64     fSummaryFMD3i(0),
65     fSummaryFMD3o(0)
66 {
67   // Default constructor 
68   Rotate("peds.csv", 3);
69   fOutputFile.open("peds.csv");
70   Rotate("ddl3072.csv", 10);
71   fZSfileFMD1.open("ddl3072.csv");
72   Rotate("ddl3073.csv", 10);
73   fZSfileFMD2.open("ddl3073.csv");
74   Rotate("ddl3074.csv", 10);
75   fZSfileFMD3.open("ddl3074.csv");  
76   fDiagnosticsFilename = "diagnosticsPedestal.root";
77 }
78
79 //_____________________________________________________________________
80 AliFMDPedestalDA::AliFMDPedestalDA(const AliFMDPedestalDA & pedDA) : 
81   AliFMDBaseDA(pedDA),
82   fCurrentChannel(1),
83   fPedSummary("PedestalSummary","pedestals",51200,0,51200),
84   fNoiseSummary("NoiseSummary","noise",51200,0,51200),
85   fZSfileFMD1(),
86   fZSfileFMD2(),
87   fZSfileFMD3(),
88   fMinTimebin(pedDA.fMinTimebin),
89   fMaxTimebin(pedDA.fMaxTimebin),
90   fSummaryFMD1i(pedDA.fSummaryFMD1i),
91   fSummaryFMD2i(pedDA.fSummaryFMD2i),
92   fSummaryFMD2o(pedDA.fSummaryFMD2o),
93   fSummaryFMD3i(pedDA.fSummaryFMD3i),
94   fSummaryFMD3o(pedDA.fSummaryFMD3o)
95 {
96   // Copy constructor 
97 }
98
99 //_____________________________________________________________________
100 AliFMDPedestalDA::~AliFMDPedestalDA() 
101
102   // Destructor.
103 }
104
105 //_____________________________________________________________________
106 void AliFMDPedestalDA::Init() 
107
108   // Initialise 
109   SetRequiredEvents(1000);
110   fMinTimebin.Reset(1024);
111   fMaxTimebin.Reset(-1);
112
113
114 }
115
116 //_____________________________________________________________________
117 void AliFMDPedestalDA::AddChannelContainer(TObjArray* sampleArray, 
118                                            UShort_t det, 
119                                            Char_t   ring, 
120                                            UShort_t sec, 
121                                            UShort_t strip) 
122 {
123   // Add a channel to the containers.
124   //
125   // Parameters:
126   //     sectorArray  Array of sectors
127   //     det          Detector 
128   //     ring         Ring
129   //     sec          Sector 
130   //     strip        Strip
131   AliFMDParameters* pars        = AliFMDParameters::Instance();
132   UInt_t            samples     = pars->GetSampleRate(det, ring, sec, strip);
133   for (UInt_t sample = 0; sample < samples; sample++) {
134     TString name(Form("FMD%d%c[%02d,%03d]_%d", det,ring,sec,strip,sample));
135     TH1S* hSample = new TH1S(name.Data(),name.Data(), 1024,-.5,1023.5);
136     hSample->SetXTitle("ADC");
137     hSample->SetYTitle("Events");
138     hSample->SetDirectory(0);
139     sampleArray->AddAtAndExpand(hSample, sample);
140   }
141 }
142
143 //_____________________________________________________________________
144 void AliFMDPedestalDA::AddSectorSummary(TObjArray* sectorArray, 
145                                         UShort_t det, 
146                                         Char_t   ring, 
147                                         UShort_t sec, 
148                                         UShort_t nStr) 
149 {
150   TH1F* sumPed = new TH1F("Pedestals", 
151                           Form("Summary of pedestals in FMD%d%c[%02d]", 
152                                det, ring, sec), 
153                           nStr, -.5, nStr-.5);
154   sumPed->SetXTitle("Strip");
155   sumPed->SetYTitle("Pedestal [ADC]");
156   sumPed->SetDirectory(0);
157   
158   TH1F* sumNoise = static_cast<TH1F*>(sumPed->Clone("Noise"));
159   sumNoise->SetYTitle("Noise [ADC]");
160   sumNoise->SetDirectory(0);
161   
162   Int_t n = sectorArray->GetEntriesFast();
163   sectorArray->AddAtAndExpand(sumPed,   n + kPedestalOffset - 1);
164   sectorArray->AddAtAndExpand(sumNoise, n + kNoiseOffset - 1);
165 }
166
167 //_____________________________________________________________________
168 void AliFMDPedestalDA::FillChannels(AliFMDDigit* digit) 
169 {
170   // Fill ADC values from a digit into the corresponding histogram.
171   //
172   // Parameters:
173   //    digit Digit to fill ADC values for.
174   UShort_t det   = digit->Detector();
175   Char_t   ring  = digit->Ring();
176   UShort_t sec   = digit->Sector();
177   UShort_t strip = digit->Strip();
178   
179   AliFMDParameters* pars     = AliFMDParameters::Instance();
180   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
181   for (UInt_t sample = 0; sample < samples; sample++) {
182     TH1S* hSample = GetChannel(det, ring, sec, strip, sample);
183     if (!hSample) continue;
184     
185     hSample->Fill(digit->Count(sample));
186   }
187   
188 }
189
190 //_____________________________________________________________________
191 void AliFMDPedestalDA::MakeSummary(UShort_t det, Char_t ring)
192 {
193   //Create summary hists for FMD pedestals
194   // std::cout << "Making summary for FMD" << det << ring << " ..." 
195   //           << std::endl;
196   switch (det) { 
197   case 1: 
198     fSummaryFMD1i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
199     break;
200   case 2:
201     switch (ring) { 
202     case 'I': case 'i':
203       fSummaryFMD2i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
204       break;
205     case 'O': case 'o':
206       fSummaryFMD2o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
207       break;
208     }
209     break;
210   case 3:
211     switch (ring) { 
212     case 'I': case 'i':
213       fSummaryFMD3i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
214       break;
215     case 'O': case 'o':
216       fSummaryFMD3o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
217       break;
218     }
219     break;
220   }
221 }
222
223 //_____________________________________________________________________
224 void AliFMDPedestalDA::Analyse(UShort_t det, 
225                                Char_t   ring, 
226                                UShort_t sec, 
227                                UShort_t strip) 
228 {
229   // Analyse a strip.  That is, compute the mean and spread of the ADC
230   // spectra for all strips.  Also output on files the values. 
231   //
232   // Parameters:
233   //     det   Detector
234   //     ring  Ring
235   //     sec   Sector 
236   //     strip Strip.
237   AliFMDParameters* pars     = AliFMDParameters::Instance();
238   TH2* summary = 0;
239   switch (det) { 
240   case 1: summary = fSummaryFMD1i; break;
241   case 2: 
242     switch (ring) { 
243     case 'I':  summary = fSummaryFMD2i; break;
244     case 'O':  summary = fSummaryFMD2o; break;
245     }
246     break;
247   case 3:
248     switch (ring) { 
249     case 'I':  summary = fSummaryFMD3i; break;
250     case 'O':  summary = fSummaryFMD3o; break;
251     }
252     break;
253   }
254 #if 0
255   static bool first = true;
256   if (summary && first) { 
257     std::cout << "Filling summary " << summary->GetName() << std::endl;
258     first = false;
259   }
260 #endif
261   // Float_t           factor   = pars->GetPedestalFactor();
262   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
263   for (UShort_t sample = 0; sample < samples; sample++) {
264   
265     TH1S* hChannel = GetChannel(det, ring, sec, strip,sample);
266     if(!hChannel || hChannel->GetEntries() == 0) {
267       //AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
268       //            det,ring,sec,strip));
269       return;
270     }
271     
272     AliDebug(50, Form("Fitting FMD%d%c[%02d,%03d] with %d entries",
273                       det,ring,sec,strip, int(hChannel->GetEntries())));
274     TF1 fitFunc("fitFunc","gausn",0,300);
275     fitFunc.SetParameters(100,100,1);
276     hChannel->Fit("fitFunc","Q0+","",10,200);
277     
278     Float_t mean = hChannel->GetMean();
279     Float_t rms  = hChannel->GetRMS();
280     
281     
282     
283     hChannel->GetXaxis()->SetRangeUser(mean-5*rms,mean+5*rms);
284     
285     mean = hChannel->GetMean();
286     rms  = hChannel->GetRMS();
287   
288     
289     UShort_t ddl, board, altro, channel;
290     UShort_t timebin;
291     
292     pars->Detector2Hardware(det,ring,sec,strip,sample,
293                             ddl,board,altro,channel,timebin);
294     Int_t idx = HWIndex(ddl, board, altro, channel);
295     if (idx >= 0) { 
296       fMinTimebin[idx] = TMath::Min(Short_t(timebin),   fMinTimebin[idx]);
297       fMaxTimebin[idx] = TMath::Max(Short_t(timebin+1), fMaxTimebin[idx]);
298     }
299     
300     std::ostream* zsFile = 0;
301     switch(det) {
302     case 1:  zsFile = &fZSfileFMD1; break;
303     case 2:  zsFile = &fZSfileFMD2; break;
304     case 3:  zsFile = &fZSfileFMD3; break;
305     default: AliWarning("Unknown sample!"); break;
306       
307     }
308     *zsFile << board   << ',' 
309             << altro   << ',' 
310             << channel << ',' 
311             << timebin << ','
312             << mean    << ',' 
313             << rms     << "\n";
314     
315     Float_t chi2ndf = 0;
316     
317     
318     if(fitFunc.GetNDF())
319       chi2ndf = fitFunc.GetChisquare() / fitFunc.GetNDF();
320     
321     
322     Int_t sampleToWrite = 2;
323     if      (samples == 2) sampleToWrite = 1;
324     else if (samples <  2) sampleToWrite = 0;
325     
326     hChannel->GetXaxis()->SetRange(1,1024);
327
328     if(sample != sampleToWrite) continue;
329     
330     
331     fOutputFile << det                         << ','
332                 << ring                        << ','
333                 << sec                         << ','
334                 << strip                       << ','
335                 << mean                        << ','
336                 << rms                         << ','
337                 << fitFunc.GetParameter(1)     << ','
338                 << fitFunc.GetParameter(2)     << ','
339                 << chi2ndf                     <<"\n";
340
341     if (summary) { 
342       Int_t bin = summary->FindBin(sec, strip);
343       summary->SetBinContent(bin, mean);
344       summary->SetBinError(bin, rms);
345     }
346
347     if(fSaveHistograms  ) {
348       TH1F* sumPed   = GetSectorSummary(det, ring, sec, true);
349       TH1F* sumNoise = GetSectorSummary(det, ring, sec, false);
350 #if 0
351       gDirectory->cd(GetSectorPath(det, ring, sec, kTRUE));
352       TH1F* sumPed   = dynamic_cast<TH1F*>(gDirectory->Get("Pedestals"));
353       TH1F* sumNoise = dynamic_cast<TH1F*>(gDirectory->Get("Noise"));
354       Int_t nStr = (ring == 'I' ? 512 : 256);
355       if (!sumPed) {
356         sumPed = new TH1F("Pedestals", 
357                           Form("Summary of pedestals in FMD%d%c[%02d]", 
358                                det, ring, sec), 
359                           nStr, -.5, nStr-.5);
360         sumPed->SetXTitle("Strip");
361         sumPed->SetYTitle("Pedestal [ADC]");
362         sumPed->SetDirectory(gDirectory);
363       }
364       if (!sumNoise) { 
365         sumNoise = new TH1F("Noise", 
366                             Form("Summary of noise in FMD%d%c[%02d]", 
367                                  det, ring, sec), 
368                             nStr, -.5, nStr-.5);
369         sumNoise->SetXTitle("Strip");
370         sumNoise->SetYTitle("Noise [ADC]");
371         
372         sumNoise->SetDirectory(gDirectory);
373       }
374 #endif
375       sumPed->SetBinContent(strip+1, mean);
376       sumPed->SetBinError(strip+1, rms);
377       sumNoise->SetBinContent(strip+1, rms);
378
379 #if 0
380       if(sumNoise->GetEntries() == nStr)
381         sumNoise->Write(sumNoise->GetName(),TObject::kOverwrite);
382       if(sumPed->GetEntries() == nStr)
383         sumPed->Write(sumPed->GetName(),TObject::kOverwrite);
384 #endif
385  
386       fPedSummary.SetBinContent(fCurrentChannel,mean);
387       
388       fNoiseSummary.SetBinContent(fCurrentChannel,rms);
389       fCurrentChannel++;
390
391 #if 0      
392       gDirectory->cd(GetStripPath(det, ring, sec, strip, kTRUE));
393       hChannel->GetXaxis()->SetRange(1,1024);
394       
395       hChannel->Write();
396 #endif
397     }
398   }
399 }
400
401 //_____________________________________________________________________
402 void AliFMDPedestalDA::Terminate(TFile* diagFile) 
403 {
404   // Called at the end of a job.  Fills in missing time-bins and
405   // closes output files
406   if(fSaveHistograms && diagFile) {
407     diagFile->cd();
408     
409     fPedSummary.Write();
410     fNoiseSummary.Write();
411   }
412   AliFMDAltroMapping* map = AliFMDParameters::Instance()->GetAltroMap();
413   for (Int_t i = 0; i < 3; i++) { 
414     std::ofstream& out = (i == 0 ? fZSfileFMD1 : 
415                           i == 1 ? fZSfileFMD2 : 
416                           fZSfileFMD3);
417     if (out.is_open() && fSeenDetectors[i]) { 
418       FillinTimebins(out, map->Detector2DDL(i+1));
419     }
420     if (!fSeenDetectors[i]) {
421       TString n(Form("ddl%d.csv",3072+map->Detector2DDL(i+1)));
422       gSystem->Unlink(n.Data());
423     }
424   }
425   
426 }
427
428 //_____________________________________________________________________
429 void AliFMDPedestalDA::FillinTimebins(std::ofstream& out, UShort_t /*ddl*/)
430 {
431   // 
432   // Fill missing timebins
433   //
434 #if 0
435   unsigned short  boards[] = { 0x0, 0x1, 0x10, 0x11, 0xFFFF };
436   unsigned short* board    = boards;
437   while ((*boards) != 0xFFFF) { 
438     for (UShort_t altro = 0; altro < 3; altro++) { 
439       for (UShort_t channel = 0; channel < 16; channel++) { 
440         Int_t idx = HWIndex(ddl, *board, altro, channel);
441         if (idx < 0) { 
442           AliWarning(Form("Invalid index for %4d/0x%02x/0x%x/0x%x: %d", 
443                           ddl, *board, altro, channel, idx));
444           continue;
445         }
446         Short_t min = fMinTimebin[idx];
447         Short_t max = fMaxTimebin[idx];
448         
449         // Channel not seen at all. 
450         if (min > 1023 || max < 0) continue;
451
452         out << "# Extra timebins for 0x" << std::hex 
453             << board << ',' << altro << ',' << channel 
454             << " got time-bins " << min << " to " << max-1 
455             << std::dec << std::endl;
456         
457         for (UShort_t t = 15; t < min; t++) 
458           // Write a phony line 
459           out << board << "," << altro << "," << channel << "," 
460               << t << "," << 1023 << "," << 0 << std::endl;
461
462         for (UShort_t t = max; t < 1024; t++) 
463           // Write a phony line 
464           out << board << "," << altro << "," << channel << "," 
465               << t << "," << 1023 << "," << 0 << std::endl;
466       } // channel loop
467     } // altro loop
468   } // board loop
469   // Write trailer, and close 
470 #endif
471   out.write("# EOF\n", 6);
472   out.close();
473 }
474
475 //_____________________________________________________________________
476 void AliFMDPedestalDA::WriteHeaderToFile() 
477 {
478   //
479   // Write headers to output files 
480   //
481   AliFMDParameters* pars       = AliFMDParameters::Instance();
482   fOutputFile.write(Form("# %s \n",pars->GetPedestalShuttleID()),13);
483   TDatime now;
484   fOutputFile << "# This file created from run # " << fRunno 
485               << " @ " << now.AsString() << std::endl;
486   fOutputFile.write("# Detector, "
487                     "Ring, "
488                     "Sector, "
489                     "Strip, "
490                     "Pedestal, "
491                     "Noise, "
492                     "Mu, "
493                     "Sigma, "
494                     "Chi2/NDF \n", 71);
495
496   std::ostream* zss[] = { &fZSfileFMD1, &fZSfileFMD2, &fZSfileFMD3, 0 };
497   for (size_t i = 0; i < 3; i++)  {
498     *(zss[i]) << "# FMD " << (i+1) << " pedestals \n"
499               << "# board, "
500               << "altro, "
501               << "channel, "
502               << "timebin, "
503               << "pedestal, "
504               << "noise\n";
505     *(zss[i]) << "# This file created from run # " << fRunno 
506               << " @ " << now.AsString() << std::endl;
507   }
508 }
509
510 //_____________________________________________________________________
511 TH1S* AliFMDPedestalDA::GetChannel(UShort_t det, 
512                                    Char_t   ring, 
513                                    UShort_t sec, 
514                                    UShort_t strip,
515                                    UInt_t   sample) 
516 {
517   // Get the histogram corresponding to a strip sample.
518   // 
519   // Parameters:
520   //     det    Detector
521   //     ring   Ring
522   //     sec    Sector
523   //     strip  Strip
524   //     sample Sample
525   //
526   // Return:
527   //     ADC spectra of a strip.
528   TObjArray* sampleArray = GetStripArray(det, ring, sec, strip);
529   if (!sampleArray) return 0;
530   TH1S*      hSample     = static_cast<TH1S*>(sampleArray->At(sample));
531   if (!hSample) {
532     AliErrorF("No channel histogram for FMD%d%c[%02d,%03d]_%d", 
533               det, ring, sec, strip, sample);
534     sampleArray->ls();
535     AliErrorF("Path is %s <- %s <- %s <- %s", 
536               sampleArray->GetName(),
537               GetSectorArray(det, ring, sec)->GetName(), 
538               GetRingArray(det, ring)->GetName(),
539               GetDetectorArray(det)->GetName());
540     
541   }
542   return hSample;
543   
544 }
545 //_____________________________________________________________________
546 TH1F* AliFMDPedestalDA::GetSectorSummary(UShort_t det, 
547                                          Char_t   ring, 
548                                          UShort_t sec, 
549                                          Bool_t   pedNotNoise) 
550 {
551   TObjArray* secArray    = GetSectorArray(det, ring, sec);
552   Int_t      n           = secArray->GetEntriesFast();
553   Int_t      i           = n - (pedNotNoise ? kNoiseOffset : kPedestalOffset);
554   return static_cast<TH1F*>(secArray->At(i));
555 }
556
557 //_____________________________________________________________________
558 //
559 //EOF
560 //