Non-implemented private copy constructor and assignment operator
[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 }
77
78 //_____________________________________________________________________
79 AliFMDPedestalDA::AliFMDPedestalDA(const AliFMDPedestalDA & pedDA) : 
80   AliFMDBaseDA(pedDA),
81   fCurrentChannel(1),
82   fPedSummary("PedestalSummary","pedestals",51200,0,51200),
83   fNoiseSummary("NoiseSummary","noise",51200,0,51200),
84   fZSfileFMD1(),
85   fZSfileFMD2(),
86   fZSfileFMD3(),
87   fMinTimebin(pedDA.fMinTimebin),
88   fMaxTimebin(pedDA.fMaxTimebin),
89   fSummaryFMD1i(pedDA.fSummaryFMD1i),
90   fSummaryFMD2i(pedDA.fSummaryFMD2i),
91   fSummaryFMD2o(pedDA.fSummaryFMD2o),
92   fSummaryFMD3i(pedDA.fSummaryFMD3i),
93   fSummaryFMD3o(pedDA.fSummaryFMD3o)
94 {
95   // Copy constructor 
96 }
97
98 //_____________________________________________________________________
99 AliFMDPedestalDA::~AliFMDPedestalDA() 
100
101   // Destructor.
102 }
103
104 //_____________________________________________________________________
105 void AliFMDPedestalDA::Init() 
106
107   // Initialise 
108   SetRequiredEvents(1000);
109   fMinTimebin.Reset(1024);
110   fMaxTimebin.Reset(-1);
111
112
113 }
114
115 //_____________________________________________________________________
116 void AliFMDPedestalDA::AddChannelContainer(TObjArray* sectorArray, 
117                                            UShort_t det, 
118                                            Char_t   ring, 
119                                            UShort_t sec, 
120                                            UShort_t strip) 
121 {
122   // Add a channel to the containers.
123   //
124   // Parameters:
125   //     sectorArray  Array of sectors
126   //     det          Detector 
127   //     ring         Ring
128   //     sec          Sector 
129   //     strip        Strip
130   AliFMDParameters* pars        = AliFMDParameters::Instance();
131   UInt_t            samples     = pars->GetSampleRate(det, ring, sec, strip);
132   TObjArray*        sampleArray = new TObjArray(samples);
133   sampleArray->SetOwner();
134   for (UInt_t sample = 0; sample < samples; sample++) {
135     TString name(Form("FMD%d%c[%02d,%03d]_%d", det,ring,sec,strip,sample));
136     TH1S* hSample = new TH1S(name.Data(),name.Data(), 1024,-.5,1023.5);
137     hSample->SetXTitle("ADC");
138     hSample->SetYTitle("Events");
139     hSample->SetDirectory(0);
140     sampleArray->AddAt(hSample, sample);
141   }
142   sectorArray->AddAtAndExpand(sampleArray, strip);
143 }
144
145 //_____________________________________________________________________
146 void AliFMDPedestalDA::FillChannels(AliFMDDigit* digit) 
147 {
148   // Fill ADC values from a digit into the corresponding histogram.
149   //
150   // Parameters:
151   //    digit Digit to fill ADC values for.
152   UShort_t det   = digit->Detector();
153   Char_t   ring  = digit->Ring();
154   UShort_t sec   = digit->Sector();
155   UShort_t strip = digit->Strip();
156   
157   AliFMDParameters* pars     = AliFMDParameters::Instance();
158   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
159   for (UInt_t sample = 0; sample < samples; sample++) {
160     TH1S* hSample = GetChannel(det, ring, sec, strip, sample);
161     hSample->Fill(digit->Count(sample));
162   }
163   
164 }
165
166 //_____________________________________________________________________
167 void AliFMDPedestalDA::MakeSummary(UShort_t det, Char_t ring)
168 {
169   //Create summary hists for FMD pedestals
170   std::cout << "Making summary for FMD" << det << ring << " ..." << std::endl;
171   switch (det) { 
172   case 1: 
173     fSummaryFMD1i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
174     break;
175   case 2:
176     switch (ring) { 
177     case 'I': case 'i':
178       fSummaryFMD2i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
179       break;
180     case 'O': case 'o':
181       fSummaryFMD2o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
182       break;
183     }
184     break;
185   case 3:
186     switch (ring) { 
187     case 'I': case 'i':
188       fSummaryFMD3i = MakeSummaryHistogram("ped", "Pedestals", det, ring);
189       break;
190     case 'O': case 'o':
191       fSummaryFMD3o = MakeSummaryHistogram("ped", "Pedestals", det, ring);
192       break;
193     }
194     break;
195   }
196 }
197
198 //_____________________________________________________________________
199 void AliFMDPedestalDA::Analyse(UShort_t det, 
200                                Char_t   ring, 
201                                UShort_t sec, 
202                                UShort_t strip) 
203 {
204   // Analyse a strip.  That is, compute the mean and spread of the ADC
205   // spectra for all strips.  Also output on files the values. 
206   //
207   // Parameters:
208   //     det   Detector
209   //     ring  Ring
210   //     sec   Sector 
211   //     strip Strip.
212   AliFMDParameters* pars     = AliFMDParameters::Instance();
213   TH2* summary = 0;
214   switch (det) { 
215   case 1: summary = fSummaryFMD1i; break;
216   case 2: 
217     switch (ring) { 
218     case 'I':  summary = fSummaryFMD2i; break;
219     case 'O':  summary = fSummaryFMD2o; break;
220     }
221     break;
222   case 3:
223     switch (ring) { 
224     case 'I':  summary = fSummaryFMD3i; break;
225     case 'O':  summary = fSummaryFMD3o; break;
226     }
227     break;
228   }
229   static bool first = true;
230   if (summary && first) { 
231     std::cout << "Filling summary " << summary->GetName() << std::endl;
232     first = false;
233   }
234
235   // Float_t           factor   = pars->GetPedestalFactor();
236   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
237   for (UShort_t sample = 0; sample < samples; sample++) {
238   
239     TH1S* hChannel = GetChannel(det, ring, sec, strip,sample);
240     if(hChannel->GetEntries() == 0) {
241       //AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
242       //            det,ring,sec,strip));
243       return;
244     }
245     
246     AliDebug(50, Form("Fitting FMD%d%c[%02d,%03d] with %d entries",
247                       det,ring,sec,strip, int(hChannel->GetEntries())));
248     TF1 fitFunc("fitFunc","gausn",0,300);
249     fitFunc.SetParameters(100,100,1);
250     hChannel->Fit("fitFunc","Q0+","",10,200);
251     
252     Float_t mean = hChannel->GetMean();
253     Float_t rms  = hChannel->GetRMS();
254     
255     
256     
257     hChannel->GetXaxis()->SetRangeUser(mean-5*rms,mean+5*rms);
258     
259     mean = hChannel->GetMean();
260     rms  = hChannel->GetRMS();
261   
262     
263     UShort_t ddl, board, altro, channel;
264     UShort_t timebin;
265     
266     pars->Detector2Hardware(det,ring,sec,strip,sample,
267                             ddl,board,altro,channel,timebin);
268     Int_t idx = HWIndex(ddl, board, altro, channel);
269     if (idx >= 0) { 
270       fMinTimebin[idx] = TMath::Min(Short_t(timebin),   fMinTimebin[idx]);
271       fMaxTimebin[idx] = TMath::Max(Short_t(timebin+1), fMaxTimebin[idx]);
272     }
273     
274     std::ostream* zsFile = 0;
275     switch(det) {
276     case 1:  zsFile = &fZSfileFMD1; break;
277     case 2:  zsFile = &fZSfileFMD2; break;
278     case 3:  zsFile = &fZSfileFMD3; break;
279     default: AliWarning("Unknown sample!"); break;
280       
281     }
282     *zsFile << board   << ',' 
283             << altro   << ',' 
284             << channel << ',' 
285             << timebin << ','
286             << mean    << ',' 
287             << rms     << "\n";
288     
289     Float_t chi2ndf = 0;
290     
291     
292     if(fitFunc.GetNDF())
293       chi2ndf = fitFunc.GetChisquare() / fitFunc.GetNDF();
294     
295     
296     Int_t sampleToWrite = 2;
297     if      (samples == 2) sampleToWrite = 1;
298     else if (samples <  2) sampleToWrite = 0;
299     
300     if(sample != sampleToWrite) continue;
301     
302     
303     fOutputFile << det                         << ','
304                 << ring                        << ','
305                 << sec                         << ','
306                 << strip                       << ','
307                 << mean                        << ','
308                 << rms                         << ','
309                 << fitFunc.GetParameter(1)     << ','
310                 << fitFunc.GetParameter(2)     << ','
311                 << chi2ndf                     <<"\n";
312
313     if (summary) { 
314       Int_t bin = summary->FindBin(sec, strip);
315       summary->SetBinContent(bin, mean);
316       summary->SetBinError(bin, rms);
317     }
318
319     if(fSaveHistograms  ) {
320       gDirectory->cd(GetSectorPath(det, ring, sec, kTRUE));
321       TH1F* sumPed   = dynamic_cast<TH1F*>(gDirectory->Get("Pedestals"));
322       TH1F* sumNoise = dynamic_cast<TH1F*>(gDirectory->Get("Noise"));
323       Int_t nStr = (ring == 'I' ? 512 : 256);
324       if (!sumPed) {
325         sumPed = new TH1F("Pedestals", 
326                           Form("Summary of pedestals in FMD%d%c[%02d]", 
327                                det, ring, sec), 
328                           nStr, -.5, nStr-.5);
329         sumPed->SetXTitle("Strip");
330         sumPed->SetYTitle("Pedestal [ADC]");
331         sumPed->SetDirectory(gDirectory);
332       }
333       if (!sumNoise) { 
334         sumNoise = new TH1F("Noise", 
335                             Form("Summary of noise in FMD%d%c[%02d]", 
336                                  det, ring, sec), 
337                             nStr, -.5, nStr-.5);
338         sumNoise->SetXTitle("Strip");
339         sumNoise->SetYTitle("Noise [ADC]");
340         
341         sumNoise->SetDirectory(gDirectory);
342       }
343       sumPed->SetBinContent(strip+1, mean);
344       sumPed->SetBinError(strip+1, rms);
345       sumNoise->SetBinContent(strip+1, rms);
346       
347       if(sumNoise->GetEntries() == nStr)
348         sumNoise->Write(sumNoise->GetName(),TObject::kOverwrite);
349       if(sumPed->GetEntries() == nStr)
350         sumPed->Write(sumPed->GetName(),TObject::kOverwrite);
351       
352       fPedSummary.SetBinContent(fCurrentChannel,mean);
353       
354       fNoiseSummary.SetBinContent(fCurrentChannel,rms);
355       fCurrentChannel++;
356       
357       gDirectory->cd(GetStripPath(det, ring, sec, strip, kTRUE));
358       hChannel->GetXaxis()->SetRange(1,1024);
359       
360       hChannel->Write();
361     }
362   }
363 }
364
365 //_____________________________________________________________________
366 void AliFMDPedestalDA::Terminate(TFile* diagFile) 
367 {
368   // Called at the end of a job.  Fills in missing time-bins and
369   // closes output files
370   if(fSaveHistograms) {
371     diagFile->cd();
372     
373     fPedSummary.Write();
374     fNoiseSummary.Write();
375   }
376   AliFMDAltroMapping* map = AliFMDParameters::Instance()->GetAltroMap();
377   for (Int_t i = 0; i < 3; i++) { 
378     std::ofstream& out = (i == 0 ? fZSfileFMD1 : 
379                           i == 1 ? fZSfileFMD2 : 
380                           fZSfileFMD3);
381     if (out.is_open() && fSeenDetectors[i]) { 
382       FillinTimebins(out, map->Detector2DDL(i+1));
383     }
384     if (!fSeenDetectors[i]) {
385       TString n(Form("ddl%d.csv",3072+map->Detector2DDL(i+1)));
386       gSystem->Unlink(n.Data());
387     }
388   }
389   
390 }
391
392 //_____________________________________________________________________
393 void AliFMDPedestalDA::FillinTimebins(std::ofstream& out, UShort_t /*ddl*/)
394 {
395   // 
396   // Fill missing timebins
397   //
398 #if 0
399   unsigned short  boards[] = { 0x0, 0x1, 0x10, 0x11, 0xFFFF };
400   unsigned short* board    = boards;
401   while ((*boards) != 0xFFFF) { 
402     for (UShort_t altro = 0; altro < 3; altro++) { 
403       for (UShort_t channel = 0; channel < 16; channel++) { 
404         Int_t idx = HWIndex(ddl, *board, altro, channel);
405         if (idx < 0) { 
406           AliWarning(Form("Invalid index for %4d/0x%02x/0x%x/0x%x: %d", 
407                           ddl, *board, altro, channel, idx));
408           continue;
409         }
410         Short_t min = fMinTimebin[idx];
411         Short_t max = fMaxTimebin[idx];
412         
413         // Channel not seen at all. 
414         if (min > 1023 || max < 0) continue;
415
416         out << "# Extra timebins for 0x" << std::hex 
417             << board << ',' << altro << ',' << channel 
418             << " got time-bins " << min << " to " << max-1 
419             << std::dec << std::endl;
420         
421         for (UShort_t t = 15; t < min; t++) 
422           // Write a phony line 
423           out << board << "," << altro << "," << channel << "," 
424               << t << "," << 1023 << "," << 0 << std::endl;
425
426         for (UShort_t t = max; t < 1024; t++) 
427           // Write a phony line 
428           out << board << "," << altro << "," << channel << "," 
429               << t << "," << 1023 << "," << 0 << std::endl;
430       } // channel loop
431     } // altro loop
432   } // board loop
433   // Write trailer, and close 
434 #endif
435   out.write("# EOF\n", 6);
436   out.close();
437 }
438
439 //_____________________________________________________________________
440 void AliFMDPedestalDA::WriteHeaderToFile() 
441 {
442   //
443   // Write headers to output files 
444   //
445   AliFMDParameters* pars       = AliFMDParameters::Instance();
446   fOutputFile.write(Form("# %s \n",pars->GetPedestalShuttleID()),13);
447   TDatime now;
448   fOutputFile << "# This file created from run # " << fRunno 
449               << " @ " << now.AsString() << std::endl;
450   fOutputFile.write("# Detector, "
451                     "Ring, "
452                     "Sector, "
453                     "Strip, "
454                     "Pedestal, "
455                     "Noise, "
456                     "Mu, "
457                     "Sigma, "
458                     "Chi2/NDF \n", 71);
459
460   std::ostream* zss[] = { &fZSfileFMD1, &fZSfileFMD2, &fZSfileFMD3, 0 };
461   for (size_t i = 0; i < 3; i++)  {
462     *(zss[i]) << "# FMD " << (i+1) << " pedestals \n"
463               << "# board, "
464               << "altro, "
465               << "channel, "
466               << "timebin, "
467               << "pedestal, "
468               << "noise\n";
469     *(zss[i]) << "# This file created from run # " << fRunno 
470               << " @ " << now.AsString() << std::endl;
471   }
472 }
473
474 //_____________________________________________________________________
475 TH1S* AliFMDPedestalDA::GetChannel(UShort_t det, 
476                                    Char_t   ring, 
477                                    UShort_t sec, 
478                                    UShort_t strip,
479                                    UInt_t   sample) 
480 {
481   // Get the histogram corresponding to a strip sample.
482   // 
483   // Parameters:
484   //     det    Detector
485   //     ring   Ring
486   //     sec    Sector
487   //     strip  Strip
488   //     sample Sample
489   //
490   // Return:
491   //     ADC spectra of a strip.
492   UShort_t   iring       = (ring == 'O' ? 0 : 1);
493   TObjArray* detArray    = static_cast<TObjArray*>(fDetectorArray.At(det));
494   TObjArray* ringArray   = static_cast<TObjArray*>(detArray->At(iring));
495   TObjArray* secArray    = static_cast<TObjArray*>(ringArray->At(sec));
496   TObjArray* sampleArray = static_cast<TObjArray*>(secArray->At(strip));
497   TH1S*      hSample     = static_cast<TH1S*>(sampleArray->At(sample));
498   return hSample;
499   
500 }
501
502 //_____________________________________________________________________
503 //
504 //EOF
505 //