- Channels (or pads) with bad addressing of manus or occupancy rate > 1 are not calib...
[u/mrichter/AliRoot.git] / MUON / AliMUONPedestal.cxx
1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 // $Id$
17
18 #include "AliMUONPedestal.h"
19 #include "AliMUONErrorCounter.h"
20 #include "AliMUONVStore.h"
21 #include "AliMUON2DMap.h"
22 #include "AliMUONCalibParamND.h"
23
24 #include <TString.h>
25 #include <THashTable.h>
26 #include <THashList.h>
27 #include <TTimeStamp.h>
28 #include <TMath.h>
29 #include <TTree.h>
30 #include <TFile.h>
31 #include <TH1F.h>
32 #include <Riostream.h>
33
34 #include <sstream>
35 //
36 //AMORE
37 //
38 #ifdef ALI_AMORE
39 #include <AmoreDA.h>
40 #include "TSystem.h"
41 #endif
42
43 //-----------------------------------------------------------------------------
44 /// \class AliMUONPedestal
45 ///
46 /// Implementation of the pedestal computing
47 ///
48 /// add
49 /// 
50 ///
51 /// \author Alberto Baldisseri, JL Charvet (05/05/2009)
52 //-----------------------------------------------------------------------------
53
54 /// \cond CLASSIMP
55 ClassImp(AliMUONPedestal)
56 /// \endcond
57
58 //______________________________________________________________________________
59 AliMUONPedestal::AliMUONPedestal()
60 : TObject(),
61 fN(0),
62 fNEvents(0),
63 fRunNumber(0),
64 fNChannel(0),
65 fNManu(0),
66 fErrorBuspatchTable(new THashTable(100,2)),
67 fManuBuspatchTable(new THashTable(100,2)),
68 fDate(new TTimeStamp()),
69 fFilcout(0),
70 fPedestalStore(new AliMUON2DMap(kFALSE)),
71 fIndex(-1)
72 {
73 /// Default constructor
74   sprintf(fHistoFileName," ");
75   sprintf(fprefixDA," "); 
76 }
77 //  AliMUONPedestal& operator=(const AliMUONPedestal& other); Copy ctor
78
79 //______________________________________________________________________________
80 AliMUONPedestal::~AliMUONPedestal()
81 {
82 /// Destructor
83 }
84
85 //______________________________________________________________________________
86 void AliMUONPedestal::MakePed(Int_t busPatchId, Int_t manuId, Int_t channelId, Int_t charge)
87 {
88   /// Compute pedestals values
89   AliMUONVCalibParam* ped = 
90     static_cast<AliMUONVCalibParam*>(fPedestalStore ->FindObject(busPatchId, manuId));
91
92   if (!ped) {
93     fNManu++;
94     ped = new AliMUONCalibParamND(2, kNChannels,busPatchId, manuId, -1.); // put default wise -1, not connected channel
95     fPedestalStore ->Add(ped);  
96   }
97
98   // Initialization for the first value
99   if (ped->ValueAsDouble(channelId, 0) == -1) ped->SetValueAsDouble(channelId, 0, 0.);
100   if (ped->ValueAsDouble(channelId, 1) == -1) ped->SetValueAsDouble(channelId, 1, 0.);
101
102   Double_t pedMean  = ped->ValueAsDouble(channelId, 0) + (Double_t) charge;
103   Double_t pedSigma = ped->ValueAsDouble(channelId, 1) + (Double_t) charge*charge;
104
105   ped->SetValueAsDouble(channelId, 0, pedMean);
106   ped->SetValueAsDouble(channelId, 1, pedSigma);
107
108   char bpmanuname[256];
109   AliMUONErrorCounter* manuCounter;
110   sprintf(bpmanuname,"bp%dmanu%d",busPatchId,manuId);           
111   if (!(manuCounter = (AliMUONErrorCounter*)GetManuBuspatchTable()->FindObject(bpmanuname)))
112     {
113       // New (buspatch,manu)
114       manuCounter = new AliMUONErrorCounter(busPatchId,manuId);
115       manuCounter->SetName(bpmanuname);
116       GetManuBuspatchTable()->Add(manuCounter);
117     }
118   else
119     {
120       // Existing buspatch
121       manuCounter->Increment();
122     }   
123 }
124
125 //______________________________________________________________________________
126 TString AliMUONPedestal::WritePedHeader(void) 
127 {
128 ///
129
130   ostringstream stream;
131   stream<<"//===========================================================================" << endl;
132   stream<<"//                       Pedestal file calculated by MUONTRKda"<<endl;
133   stream<<"//===========================================================================" << endl;
134   stream<<"//       * Run           : " << fRunNumber << endl; 
135   stream<<"//       * Date          : " << fDate->AsString("l") <<endl;
136   stream<<"//       * Statictics    : " << fNEvents << endl;
137   stream<<"//       * # of MANUS    : " << fNManu << endl;
138   stream<<"//       * # of channels : " << fNChannel << endl;
139   if (fErrorBuspatchTable->GetSize())
140   {
141     stream<<"//"<<endl;
142     stream<<"//       * Buspatches with less statistics (due to parity errors)"<<endl;
143     TIterator* iter = fErrorBuspatchTable->MakeIterator();
144     AliMUONErrorCounter* parityerror;
145     while((parityerror = (AliMUONErrorCounter*) iter->Next()))
146     {
147       stream<<"//         bp "<<parityerror->BusPatch()<<" events used "<<fNEvents-parityerror->Events()<<endl;
148     }
149   }  
150   stream<<"//"<<endl;
151   stream<<"//---------------------------------------------------------------------------" << endl;
152   stream<<"//---------------------------------------------------------------------------" << endl;
153   stream<<"//      BP     MANU     CH.      MEAN    SIGMA"<<endl;
154   stream<<"//---------------------------------------------------------------------------" << endl;
155
156   return TString(stream.str().c_str());
157 }
158
159 //______________________________________________________________________________
160 TString AliMUONPedestal::WritePedData(Int_t BP, Int_t Manu, Int_t ch, Double_t pedMean, Double_t pedSigma) 
161 {
162 ///
163
164   ostringstream stream("");
165   stream << "\t" << BP << "\t" << Manu <<"\t"<< ch << "\t"
166          << pedMean <<"\t"<< pedSigma << endl;
167   return TString(stream.str().c_str());
168
169 }
170
171 //______________________________________________________________________________
172 void AliMUONPedestal::MakePedStore(TString shuttleFile_1 = "") 
173 {
174
175   /// Store pedestals in ASCII files
176   Double_t pedMean;
177   Double_t pedSigma;
178   ofstream fileout;
179 #ifdef ALI_AMORE
180   ostringstream stringout; // String to be sent to AMORE_DB
181 #endif
182   TString tempstring;  
183   Int_t busPatchId;
184   Int_t manuId;
185   Int_t channelId;
186
187 // histo
188   TFile*  histoFile = 0;
189   TTree* tree = 0;
190   TH1F* pedMeanHisto = 0;
191   TH1F* pedSigmaHisto = 0;
192
193   if (fIndex<0) // Pedestal run (fIndex=-1)
194   {
195     sprintf(fHistoFileName,"%s.root",fprefixDA);
196     histoFile = new TFile(fHistoFileName,"RECREATE","MUON Tracking pedestals");
197
198     Char_t name[255];
199     Char_t title[255];
200     sprintf(name,"pedmean_allch");
201     sprintf(title,"Pedestal mean all channels");
202     Int_t nx = kADCMax+1;
203     Int_t xmin = 0;
204     Int_t xmax = kADCMax; 
205     pedMeanHisto = new TH1F(name,title,nx,xmin,xmax);
206     pedMeanHisto->SetDirectory(histoFile);
207
208     sprintf(name,"pedsigma_allch");
209     sprintf(title,"Pedestal sigma all channels");
210     nx = 201;
211     xmin = 0;
212     xmax = 200; 
213     pedSigmaHisto = new TH1F(name,title,nx,xmin,xmax);
214     pedSigmaHisto->SetDirectory(histoFile);
215
216     tree = new TTree("t","Pedestal tree");
217     tree->Branch("bp",&busPatchId,"bp/I");
218     tree->Branch("manu",&manuId,",manu/I");
219     tree->Branch("channel",&channelId,",channel/I");
220     tree->Branch("pedMean",&pedMean,",pedMean/D");
221     tree->Branch("pedSigma",&pedSigma,",pedSigma/D");
222   }
223
224   if (!shuttleFile_1.IsNull()) {
225     fileout.open(shuttleFile_1.Data());
226     tempstring = WritePedHeader();
227     fileout << tempstring;
228 #ifdef ALI_AMORE
229     stringout << tempstring;
230 #endif
231   }
232   // print in logfile
233   if (fErrorBuspatchTable->GetSize())
234   {
235     cout<<"\n* Buspatches with less statistics (due to parity errors)"<<endl;
236     (*fFilcout)<<"\n* Buspatches with less statistics (due to parity errors)"<<endl;
237     TIterator* iter = fErrorBuspatchTable->MakeIterator();
238     AliMUONErrorCounter* parityerror;
239     while((parityerror = (AliMUONErrorCounter*) iter->Next()))
240     {
241       cout<<"  bp "<<parityerror->BusPatch()<<": events used = "<<fNEvents-parityerror->Events()<<endl;
242       (*fFilcout)<<"  bp "<<parityerror->BusPatch()<<": events used = "<<fNEvents-parityerror->Events()<<endl;
243     }
244
245   }
246
247
248 // iterator over pedestal
249   TIter next(fPedestalStore ->CreateIterator());
250   AliMUONVCalibParam* ped;
251
252   while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) )
253   {
254     busPatchId              = ped->ID0();
255     manuId                  = ped->ID1();
256     if(manuId==0)
257       {
258         cout << " !!! BIG WARNING: ManuId = " << manuId << " !!! in  BP = " << busPatchId << endl;
259         (*fFilcout) << " !!! BIG WARNING: ManuId = " << manuId << " !!! in  BP = " << busPatchId << endl;
260       }
261     Int_t eventCounter;
262     // Correct the number of events for buspatch with errors
263     char bpname[256];
264     AliMUONErrorCounter* errorCounter;
265     sprintf(bpname,"bp%d",busPatchId);
266     if ((errorCounter = (AliMUONErrorCounter*)fErrorBuspatchTable->FindObject(bpname)))
267     {
268       eventCounter = fNEvents - errorCounter->Events();
269     }
270     else
271     {
272       eventCounter = fNEvents;
273     }
274
275     Int_t occupancy;
276     // value of (buspatch, manu) occupancy
277     char bpmanuname[256];
278     AliMUONErrorCounter* manuCounter;
279     sprintf(bpmanuname,"bp%dmanu%d",busPatchId,manuId);
280     manuCounter = (AliMUONErrorCounter*)fManuBuspatchTable->FindObject(bpmanuname);
281     occupancy = manuCounter->Events()/64/eventCounter;
282     if(occupancy>1)
283       {
284         cout << " !!! BIG WARNING: ManuId = " << manuId << " !!! in  BP = " << busPatchId << " occupancy (>1) = " << occupancy << endl;
285         (*fFilcout) << " !!! BIG WARNING: ManuId = " << manuId << " !!! in  BP = " << busPatchId << " occupancy (>1) = " << occupancy <<endl;
286      }
287
288     for (channelId = 0; channelId < ped->Size() ; ++channelId) {
289       pedMean  = ped->ValueAsDouble(channelId, 0);
290
291       if (pedMean > 0) // connected channels
292         {
293           ped->SetValueAsDouble(channelId, 0, pedMean/(Double_t)eventCounter);
294           pedMean  = ped->ValueAsDouble(channelId, 0);
295           pedSigma = ped->ValueAsDouble(channelId, 1);
296           ped->SetValueAsDouble(channelId, 1, TMath::Sqrt(TMath::Abs(pedSigma/(Double_t)eventCounter - pedMean*pedMean)));
297           if(manuId == 0 || occupancy>1)
298             {
299               ped->SetValueAsDouble(channelId, 0, kADCMax);
300               ped->SetValueAsDouble(channelId, 1, kADCMax);
301             }
302         }
303       else
304         {
305           ped->SetValueAsDouble(channelId, 0, kADCMax);
306           ped->SetValueAsDouble(channelId, 1, kADCMax);
307         }
308
309       pedMean  = ped->ValueAsDouble(channelId, 0);
310       pedSigma = ped->ValueAsDouble(channelId, 1);
311
312         if (!shuttleFile_1.IsNull()) {
313           tempstring = WritePedData(busPatchId,manuId,channelId,pedMean,pedSigma);
314           fileout << tempstring;
315 #ifdef ALI_AMORE
316           stringout << tempstring;
317 #endif
318         }
319         if(fIndex<0) // Pedestal Run
320         {
321           pedMeanHisto->Fill(pedMean);
322           pedSigmaHisto->Fill(pedSigma);
323           tree->Fill();
324         }
325     }
326   }
327
328 // file outputs
329   if (!shuttleFile_1.IsNull())  fileout.close();
330
331 // Outputs to root file and AMORE DB
332 #ifdef ALI_AMORE
333   //
334   //Send objects to the AMORE DB
335   //
336     amore::da::AmoreDA amoreDA(amore::da::AmoreDA::kSender);
337     TObjString peddata(stringout.str().c_str());
338     Int_t status =0;
339     status = amoreDA.Send("Pedestals",&peddata);
340     if ( status )
341       cout << "Warning: Failed to write Pedestals in the AMORE database : " << status << endl;
342     else 
343       cout << "amoreDA.Send(Pedestals) ok" << endl;  
344 #else
345   cout << "Warning: MCH DA not compiled with AMORE support" << endl;
346 #endif
347   if(fIndex<0) // Pedestal Run
348   { 
349     histoFile->Write();  
350     histoFile->Close(); 
351     delete fPedestalStore ;
352   }
353 }