Moving the FMD offline trigger to its own class, AliFMDOfflineTrigger
[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 <iostream>
33 #include <fstream>
34 #include "AliLog.h"
35 #include "TF1.h"
36 #include "TObject.h"
37 #include "TMath.h"
38 #include <TSystem.h>
39
40 //_____________________________________________________________________
41 ClassImp(AliFMDPedestalDA)
42
43 //_____________________________________________________________________
44 AliFMDPedestalDA::AliFMDPedestalDA() : AliFMDBaseDA(),
45   fCurrentChannel(1),
46   fPedSummary("PedestalSummary","pedestals",51200,0,51200),
47   fNoiseSummary("NoiseSummary","noise",51200,0,51200),
48   fZSfileFMD1(),
49   fZSfileFMD2(),
50   fZSfileFMD3(), 
51   fMinTimebin(3 * 4 * 3 * 16), // 3 ddls, 4 FECs, 3 Altros, 16 channels
52   fMaxTimebin(3 * 4 * 3 * 16)  // 3 ddls, 4 FECs, 3 Altros, 16 channels
53 {
54   // Default constructor 
55   fOutputFile.open("peds.csv");
56   fZSfileFMD1.open("ddl3072.csv");
57   fZSfileFMD2.open("ddl3073.csv");
58   fZSfileFMD3.open("ddl3074.csv");  
59 }
60
61 //_____________________________________________________________________
62 AliFMDPedestalDA::AliFMDPedestalDA(const AliFMDPedestalDA & pedDA) : 
63   AliFMDBaseDA(pedDA),
64   fCurrentChannel(1),
65   fPedSummary("PedestalSummary","pedestals",51200,0,51200),
66   fNoiseSummary("NoiseSummary","noise",51200,0,51200),
67   fZSfileFMD1(),
68   fZSfileFMD2(),
69   fZSfileFMD3(),
70   fMinTimebin(pedDA.fMinTimebin),
71   fMaxTimebin(pedDA.fMaxTimebin)
72 {
73   // Copy constructor 
74 }
75
76 //_____________________________________________________________________
77 AliFMDPedestalDA::~AliFMDPedestalDA() 
78
79   // Destructor.
80 }
81
82 //_____________________________________________________________________
83 void AliFMDPedestalDA::Init() 
84
85   // Initialise 
86   SetRequiredEvents(1000);
87   fMinTimebin.Reset(1024);
88   fMaxTimebin.Reset(-1);
89 }
90
91 //_____________________________________________________________________
92 void AliFMDPedestalDA::AddChannelContainer(TObjArray* sectorArray, 
93                                            UShort_t det, 
94                                            Char_t   ring, 
95                                            UShort_t sec, 
96                                            UShort_t strip) 
97 {
98   // Add a channel to the containers.
99   //
100   // Parameters:
101   //     sectorArray  Array of sectors
102   //     det          Detector 
103   //     ring         Ring
104   //     sec          Sector 
105   //     strip        Strip
106   AliFMDParameters* pars        = AliFMDParameters::Instance();
107   UInt_t            samples     = pars->GetSampleRate(det, ring, sec, strip);
108   TObjArray*        sampleArray = new TObjArray(samples);
109   sampleArray->SetOwner();
110   for (UInt_t sample = 0; sample < samples; sample++) {
111     TH1S* hSample = new TH1S(Form("FMD%d%c[%02d,03%d]_%d",
112                                   det,ring,sec,strip,sample),
113                              Form("FMD%d%c[%02d,%03%d]_%d",
114                                   det,ring,sec,strip),
115                              1024,-.5,1023.5);
116     hSample->SetXTitle("ADC");
117     hSample->SetYTitle("Events");
118     hSample->SetDirectory(0);
119     sampleArray->AddAt(hSample, sample);
120   }
121   sectorArray->AddAtAndExpand(sampleArray, strip);
122 }
123
124 //_____________________________________________________________________
125 void AliFMDPedestalDA::FillChannels(AliFMDDigit* digit) 
126 {
127   // Fill ADC values from a digit into the corresponding histogram.
128   //
129   // Parameters:
130   //    digit Digit to fill ADC values for.
131   UShort_t det   = digit->Detector();
132   Char_t   ring  = digit->Ring();
133   UShort_t sec   = digit->Sector();
134   UShort_t strip = digit->Strip();
135   
136   AliFMDParameters* pars     = AliFMDParameters::Instance();
137   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
138   for (UInt_t sample = 0; sample < samples; sample++) {
139     TH1S* hSample = GetChannel(det, ring, sec, strip, sample);
140     hSample->Fill(digit->Count(sample));
141   }
142   
143 }
144
145 //_____________________________________________________________________
146 void AliFMDPedestalDA::Analyse(UShort_t det, 
147                                Char_t   ring, 
148                                UShort_t sec, 
149                                UShort_t strip) 
150 {
151   // Analyse a strip.  That is, compute the mean and spread of the ADC
152   // spectra for all strips.  Also output on files the values. 
153   //
154   // Parameters:
155   //     det   Detector
156   //     ring  Ring
157   //     sec   Sector 
158   //     strip Strip.
159   AliFMDParameters* pars     = AliFMDParameters::Instance();
160   // Float_t           factor   = pars->GetPedestalFactor();
161   UInt_t            samples  = pars->GetSampleRate(det, ring, sec, strip);
162   for (UShort_t sample = 0; sample < samples; sample++) {
163   
164     TH1S* hChannel = GetChannel(det, ring, sec, strip,sample);
165     if(hChannel->GetEntries() == 0) {
166       //AliWarning(Form("No entries for FMD%d%c, sector %d, strip %d",
167       //            det,ring,sec,strip));
168       return;
169     }
170     
171     AliDebug(50, Form("Fitting FMD%d%c_%d_%d with %d entries",
172                       det,ring,sec,strip, hChannel->GetEntries()));
173     TF1 fitFunc("fitFunc","gausn",0,300);
174     fitFunc.SetParameters(100,100,1);
175     hChannel->Fit("fitFunc","Q0+","",10,200);
176     
177     Float_t mean = hChannel->GetMean();
178     Float_t rms  = hChannel->GetRMS();
179     
180     
181     
182     hChannel->GetXaxis()->SetRangeUser(mean-5*rms,mean+5*rms);
183     
184     mean = hChannel->GetMean();
185     rms  = hChannel->GetRMS();
186   
187     
188     UShort_t ddl, board, altro, channel;
189     UShort_t timebin;
190     
191     pars->Detector2Hardware(det,ring,sec,strip,sample,
192                             ddl,board,altro,channel,timebin);
193     Int_t idx = HWIndex(ddl, board, altro, channel);
194     if (idx >= 0) { 
195       fMinTimebin[idx] = TMath::Min(Short_t(timebin),   fMinTimebin[idx]);
196       fMaxTimebin[idx] = TMath::Max(Short_t(timebin+1), fMaxTimebin[idx]);
197     }
198     
199     std::ostream* zsFile = 0;
200     switch(det) {
201     case 1:  zsFile = &fZSfileFMD1; break;
202     case 2:  zsFile = &fZSfileFMD2; break;
203     case 3:  zsFile = &fZSfileFMD3; break;
204     default: AliWarning("Unknown sample!"); break;
205       
206     }
207     *zsFile << board   << ',' 
208             << altro   << ',' 
209             << channel << ',' 
210             << timebin << ','
211             << mean    << ',' 
212             << rms     << "\n";
213     
214     Float_t chi2ndf = 0;
215     
216     
217     if(fitFunc.GetNDF())
218       chi2ndf = fitFunc.GetChisquare() / fitFunc.GetNDF();
219     
220     
221     Int_t sampleToWrite = 2;
222     if      (samples == 2) sampleToWrite = 1;
223     else if (samples <  2) sampleToWrite = 0;
224     
225     if(sample != sampleToWrite) continue;
226     
227     
228     fOutputFile << det                         << ','
229                 << ring                        << ','
230                 << sec                         << ','
231                 << strip                       << ','
232                 << mean                        << ','
233                 << rms                         << ','
234                 << fitFunc.GetParameter(1)     << ','
235                 << fitFunc.GetParameter(2)     << ','
236                 << chi2ndf                     <<"\n";
237     
238     if(fSaveHistograms  ) {
239       gDirectory->cd(GetSectorPath(det, ring, sec, kTRUE));
240       TH1F* sumPed   = dynamic_cast<TH1F*>(gDirectory->Get("Pedestals"));
241       TH1F* sumNoise = dynamic_cast<TH1F*>(gDirectory->Get("Noise"));
242       Int_t nStr = (ring == 'I' ? 512 : 256);
243       if (!sumPed) {
244         sumPed = new TH1F("Pedestals", 
245                           Form("Summary of pedestals in FMD%d%c[%02d]", 
246                                det, ring, sec), 
247                           nStr, -.5, nStr-.5);
248         sumPed->SetXTitle("Strip");
249         sumPed->SetYTitle("Pedestal [ADC]");
250         sumPed->SetDirectory(gDirectory);
251       }
252       if (!sumNoise) { 
253         sumNoise = new TH1F("Noise", 
254                             Form("Summary of noise in FMD%d%c[%02d]", 
255                                  det, ring, sec), 
256                             nStr, -.5, nStr-.5);
257         sumNoise->SetXTitle("Strip");
258         sumNoise->SetYTitle("Noise [ADC]");
259         
260         sumNoise->SetDirectory(gDirectory);
261       }
262       sumPed->SetBinContent(strip+1, mean);
263       sumPed->SetBinError(strip+1, rms);
264       sumNoise->SetBinContent(strip+1, rms);
265       
266       if(sumNoise->GetEntries() == nStr)
267         sumNoise->Write(sumNoise->GetName(),TObject::kOverwrite);
268       if(sumPed->GetEntries() == nStr)
269         sumPed->Write(sumPed->GetName(),TObject::kOverwrite);
270       
271       fPedSummary.SetBinContent(fCurrentChannel,mean);
272       
273       fNoiseSummary.SetBinContent(fCurrentChannel,rms);
274       fCurrentChannel++;
275       
276       gDirectory->cd(GetStripPath(det, ring, sec, strip, kTRUE));
277       hChannel->GetXaxis()->SetRange(1,1024);
278       
279       hChannel->Write();
280     }
281   }
282 }
283
284 //_____________________________________________________________________
285 void AliFMDPedestalDA::Terminate(TFile* diagFile) 
286 {
287   // Called at the end of a job.  Fills in missing time-bins and
288   // closes output files
289   if(fSaveHistograms) {
290     diagFile->cd();
291     
292     fPedSummary.Write();
293     fNoiseSummary.Write();
294   }
295   AliFMDAltroMapping* map = AliFMDParameters::Instance()->GetAltroMap();
296   for (Int_t i = 0; i < 3; i++) { 
297     std::ofstream& out = (i == 0 ? fZSfileFMD1 : 
298                           i == 1 ? fZSfileFMD2 : 
299                           fZSfileFMD3);
300     if (out.is_open() && fSeenDetectors[i]) { 
301       FillinTimebins(out, map->Detector2DDL(i+1));
302     }
303     if (!fSeenDetectors[i]) {
304       TString n(Form("ddl%d.csv",3072+map->Detector2DDL(i+1)));
305       gSystem->Unlink(n.Data());
306     }
307   }
308   
309 }
310
311 //_____________________________________________________________________
312 void AliFMDPedestalDA::FillinTimebins(std::ofstream& out, UShort_t /*ddl*/)
313 {
314 #if 0
315   unsigned short  boards[] = { 0x0, 0x1, 0x10, 0x11, 0xFFFF };
316   unsigned short* board    = boards;
317   while ((*boards) != 0xFFFF) { 
318     for (UShort_t altro = 0; altro < 3; altro++) { 
319       for (UShort_t channel = 0; channel < 16; channel++) { 
320         Int_t idx = HWIndex(ddl, *board, altro, channel);
321         if (idx < 0) { 
322           AliWarning(Form("Invalid index for %4d/0x%02x/0x%x/0x%x: %d", 
323                           ddl, *board, altro, channel, idx));
324           continue;
325         }
326         Short_t min = fMinTimebin[idx];
327         Short_t max = fMaxTimebin[idx];
328         
329         // Channel not seen at all. 
330         if (min > 1023 || max < 0) continue;
331
332         out << "# Extra timebins for 0x" << std::hex 
333             << board << ',' << altro << ',' << channel 
334             << " got time-bins " << min << " to " << max-1 
335             << std::dec << std::endl;
336         
337         for (UShort_t t = 15; t < min; t++) 
338           // Write a phony line 
339           out << board << "," << altro << "," << channel << "," 
340               << t << "," << 1023 << "," << 0 << std::endl;
341
342         for (UShort_t t = max; t < 1024; t++) 
343           // Write a phony line 
344           out << board << "," << altro << "," << channel << "," 
345               << t << "," << 1023 << "," << 0 << std::endl;
346       } // channel loop
347     } // altro loop
348   } // board loop
349   // Write trailer, and close 
350 #endif
351   out.write("# EOF\n", 6);
352   out.close();
353 }
354
355 //_____________________________________________________________________
356 void AliFMDPedestalDA::WriteHeaderToFile() 
357 {
358   // Write headers to output files 
359   AliFMDParameters* pars       = AliFMDParameters::Instance();
360   fOutputFile.write(Form("# %s \n",pars->GetPedestalShuttleID()),13);
361   fOutputFile.write("# Detector, "
362                     "Ring, "
363                     "Sector, "
364                     "Strip, "
365                     "Pedestal, "
366                     "Noise, "
367                     "Mu, "
368                     "Sigma, "
369                     "Chi2/NDF \n", 71);
370
371   std::ostream* zss[] = { &fZSfileFMD1, &fZSfileFMD2, &fZSfileFMD3, 0 };
372   for (size_t i = 0; i < 3; i++) 
373     *(zss[i]) << "# FMD 1 pedestals \n"
374               << "# board, "
375               << "altro, "
376               << "channel, "
377               << "timebin, "
378               << "pedestal, "
379               << "noise\n";
380 }
381
382 //_____________________________________________________________________
383 TH1S* AliFMDPedestalDA::GetChannel(UShort_t det, 
384                                    Char_t   ring, 
385                                    UShort_t sec, 
386                                    UShort_t strip,
387                                    UInt_t   sample) 
388 {
389   // Get the histogram corresponding to a strip sample.
390   // 
391   // Parameters:
392   //     det    Detector
393   //     ring   Ring
394   //     sec    Sector
395   //     strip  Strip
396   //     sample Sample
397   //
398   // Return:
399   //     ADC spectra of a strip.
400   UShort_t   iring       = (ring == 'O' ? 0 : 1);
401   TObjArray* detArray    = static_cast<TObjArray*>(fDetectorArray.At(det));
402   TObjArray* ringArray   = static_cast<TObjArray*>(detArray->At(iring));
403   TObjArray* secArray    = static_cast<TObjArray*>(ringArray->At(sec));
404   TObjArray* sampleArray = static_cast<TObjArray*>(secArray->At(strip));
405   TH1S*      hSample     = static_cast<TH1S*>(sampleArray->At(sample));
406   return hSample;
407   
408 }
409
410 //_____________________________________________________________________
411 //
412 //EOF
413 //