92d8f7bf3a66628c02fd7744309d28a93fcbfd59
[u/mrichter/AliRoot.git] / MUON / MUONCDB.C
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 /// By Laurent Aphecetche
19
20 #if !defined(__CINT__) || defined(__MAKECINT__)
21
22 #include "MUONCDB.h"
23
24 #include "AliCDBEntry.h"
25 #include "AliCDBManager.h"
26 #include "AliMUON1DArray.h"
27 #include "AliMUON2DMap.h"
28 #include "AliMUONCalibParam1I.h"
29 #include "AliMUONCalibParam2F.h"
30 #include "AliMUONConstants.h"
31 #include "AliMUONObjectPair.h"
32 #include "AliMUONTriggerEfficiencyCells.h"
33 #include "AliMUONTriggerLut.h"
34 #include "AliMUONV2DStore.h"
35 #include "AliMUONVCalibParam.h"
36 #include "AliMUONVDataIterator.h"
37 #include "AliMpDEIterator.h"
38 #include "AliMpDEManager.h"
39 #include "AliMpManuList.h"
40 #include "AliMpSegFactory.h"
41 #include "AliMpStationType.h"
42 #include "AliMpVSegmentation.h"
43 #include "Riostream.h"
44 #include "TH1F.h"
45 #include "TList.h"
46 #include "TRandom.h"
47 #include "TStopwatch.h"
48 #include "TSystem.h"
49 #endif
50
51 //_____________________________________________________________________________
52 Int_t countChannels(AliMpVSegmentation& seg)
53 {
54   Int_t n(0);
55   
56   for ( Int_t ix = 0; ix < seg.MaxPadIndexX(); ++ix )
57   {
58     for ( Int_t iy = 0; iy < seg.MaxPadIndexY(); ++iy )
59     {
60       if ( seg.HasPad(AliMpIntPair(ix,iy)) ) ++n;
61     }
62   }
63   return n;
64 }
65
66 //_____________________________________________________________________________
67 AliMpSegFactory* segFactory()
68 {
69   static AliMpSegFactory* sf = new AliMpSegFactory();
70   return sf;
71 }
72
73 //_____________________________________________________________________________
74 void countChannels()
75 {
76   AliMpDEIterator it;
77   Int_t ntotal(0);
78   Int_t ntracker(0);
79   Int_t ntrigger(0);
80   
81   for ( Int_t station = 0; station < AliMUONConstants::NCh(); ++station )
82   {
83     Int_t n(0);
84     it.First(station);
85     while (!it.IsDone())
86     {
87       Int_t de = it.CurrentDE();
88       for ( Int_t cathode = 0; cathode < 2; ++cathode )
89       {
90         AliMpVSegmentation* seg = segFactory()->CreateMpSegmentation(de,cathode);
91         n += countChannels(*seg);
92       }
93       it.Next();
94     }
95     cout << "Station " << station << " has " << n << " channels" << endl;
96     if ( station < AliMUONConstants::NTrackingCh() )
97     {
98       ntracker += n;
99     }
100     else
101     {
102       ntrigger += n;
103     }
104     ntotal += n;
105   }
106   cout << "Tracker channels = " << ntracker << endl;
107   cout << "Trigger channels = " << ntrigger << endl;
108   cout << "Total channels =" << ntotal << endl;
109 }
110
111 //_____________________________________________________________________________
112 AliMUONV2DStore* read2D(const char* calibType, Int_t runNumber)
113 {
114   AliCDBManager* man = AliCDBManager::Instance();
115   man->SetDefaultStorage(CDBPath);
116
117   AliCDBEntry* entry = man->Get(calibType,runNumber);
118
119   if (entry)
120     {
121       return (AliMUONV2DStore*)entry->GetObject();
122     }
123   return 0;
124 }
125
126 //_____________________________________________________________________________
127 AliMUONV1DStore* read1D(const char* calibType, Int_t runNumber)
128 {
129   AliCDBManager* man = AliCDBManager::Instance();
130   man->SetDefaultStorage(CDBPath);
131   
132   AliCDBEntry* entry = man->Get(calibType,runNumber);
133   
134   if (entry)
135   {
136     return (AliMUONV1DStore*)entry->GetObject();
137   }
138   return 0;
139 }
140
141 //_____________________________________________________________________________
142 void checkCDB(const char* calibType)
143 {
144   TString c(calibType);
145   Float_t refValue(0);
146   
147   if ( c == "MUON/Calib/DeadChannels" )
148   {
149     refValue=5;
150   }
151    
152   AliMUONV2DStore* store = read2D(calibType);
153   if (!store) return;
154   
155   TList* list = AliMpManuList::ManuList();
156   TIter next(list);
157   AliMpIntPair* p;
158   
159   while ( ( p = (AliMpIntPair*)next() ) )
160   {
161     Int_t detElemId = p->GetFirst();
162     Int_t manuId = p->GetSecond();
163     
164     AliMUONVCalibParam* value = 
165       dynamic_cast<AliMUONVCalibParam*>(store->Get(detElemId,manuId));
166     
167     if (value)
168     {
169       for ( Int_t manuChannel = 0; manuChannel < value->Size(); ++manuChannel )
170       {
171         Float_t testValue = value->ValueAsFloat(manuChannel,0);
172         if ( testValue && testValue != refValue )
173         {
174           cout << "Got a strange value for DE=" << detElemId << " manuId="
175           << manuId << " manuChannel=" << manuChannel << " was expecting "
176           << refValue << " and I got " << testValue << endl;
177         }
178       }
179     }
180     else
181     {
182       cout << "Got a null value for DE=" << detElemId << " manuId="
183       << manuId << endl;
184     }
185   }
186   
187   delete list;
188   delete store;
189 }
190
191 //_____________________________________________________________________________
192 //void testDump(AliMUONV2DStore& store, int n)
193 //{  
194 //  AliMUONObjectPair* p;
195 //  
196 //  Int_t c(0);
197 //  
198 //  for ( Int_t i = 0; i < n; ++i )
199 //  {
200 //    AliMUONVDataIterator* it = store.Iterator();
201 //    
202 //    while ( ( p = dynamic_cast<AliMUONObjectPair*>(it->Next() ) ) )
203 //    {
204 //      AliMpIntPair* dm = dynamic_cast<AliMpIntPair*>(p->Key());
205 //      if (dm)
206 //      {
207 //        Int_t a(dm->GetFirst()+dm->GetSecond());
208 //        a=2;
209 //        ++c;
210 //        AliMUONVCalibParam* param = dynamic_cast<AliMUONVCalibParam*>(p->Value());
211 //      }
212 //      delete p;
213 //    }
214 //    delete it;
215 //  }
216 //    
217 //  cout << c << endl;
218 //}
219
220 //_____________________________________________________________________________
221 AliMUONV2DStore* diff(AliMUONV2DStore& store1, AliMUONV2DStore& store2, 
222                       const char* /* opt */)
223 {
224   // creates a store which contains store1-store2
225   // if opt="abs" the difference is absolute one,
226   // if opt="rel" then what is stored is (store1-store2)/store1
227   
228   AliMUONV2DStore* d = static_cast<AliMUONV2DStore*>(store1.Clone());
229
230   AliMUONVDataIterator* it = d->Iterator();
231   
232   AliMUONObjectPair* p;
233   
234   while ( ( p = dynamic_cast<AliMUONObjectPair*>(it->Next() ) ) )
235   {
236     AliMpIntPair* dm = dynamic_cast<AliMpIntPair*>(p->Key());
237     //FIXME: this might happen (if a full manu is missing, for instance)
238     //handle it.
239     assert(dm!=0);
240     AliMUONVCalibParam* param = dynamic_cast<AliMUONVCalibParam*>(p->Value());
241     //FIXMENOT: this should *not* happen
242     assert(param!=0);
243
244     AliMUONVCalibParam* param2 = dynamic_cast<AliMUONVCalibParam*>(store2.Get(dm->GetFirst(),dm->GetSecond()));
245     //FIXME: this might happen. Handle it.
246     assert(param2!=0);
247     
248     for ( Int_t i = 0; i < param->Size(); ++i )
249     {
250       for ( Int_t j = 0; j < param->Dimension(); ++j )
251       {
252         param->SetValueAsFloat(i,j,param->ValueAsFloat(i,j) - param2->ValueAsFloat(i,j));
253       }      
254     }
255   }
256   return d;
257 }
258
259 //_____________________________________________________________________________
260 void getBoundaries(const AliMUONV2DStore& store,
261                    Float_t& x0min, Float_t& x0max,
262                    Float_t& x1min, Float_t& x1max)
263 {
264   x0min=1E30;
265   x0max=-1E30;
266   x1min=1E30;
267   x1max=-1E30;
268   
269   TList* list = AliMpManuList::ManuList();
270   TIter next(list);
271   AliMpIntPair* p;
272   
273   while ( ( p = (AliMpIntPair*)next() ) )
274   {
275     Int_t detElemId = p->GetFirst();
276     Int_t manuId = p->GetSecond();
277         
278     AliMpVSegmentation* seg = 
279       segFactory()->CreateMpSegmentationByElectronics(detElemId,manuId);
280           
281     AliMUONVCalibParam* value = 
282       dynamic_cast<AliMUONVCalibParam*>(store.Get(detElemId,manuId));
283     
284     if (!value) continue;
285
286     for ( Int_t manuChannel = 0; manuChannel < value->Size(); ++manuChannel )
287     {
288       AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,manuChannel),kFALSE);
289       if (!pad.IsValid()) continue;
290
291       Float_t x0 = value->ValueAsFloat(manuChannel,0);
292       
293       x0min = TMath::Min(x0min,x0);
294       x0max = TMath::Max(x0max,x0);
295       if ( value->Dimension()>1 )
296       {
297         Float_t x1 = value->ValueAsFloat(manuChannel,1);
298         x1min = TMath::Min(x1min,x1);
299         x1max = TMath::Max(x1max,x1);
300       }
301     }
302   }  
303   delete list;
304 }
305
306 //_____________________________________________________________________________
307 void plot(const AliMUONV2DStore& store, const char* name, Int_t nbins)
308 {
309   Float_t x0min, x0max, x1min, x1max;
310   
311   getBoundaries(store,x0min,x0max,x1min,x1max);
312   
313   if ( x0min > x0max ) 
314   {
315     cerr << Form("Something is wrong with boundaries : x0(min,max)=%e,%e",
316                  x0min,x0max) << endl;
317     return;
318   }
319
320   if ( TMath::Abs(x0min-x0max) < 1E-3 ) 
321   {
322     x0min -= 1;
323     x0max += 1;
324   }
325   
326   TH1* h0 = new TH1F(Form("%s_0",name),Form("%s_0",name),
327                     nbins,x0min,x0max);
328   
329   TH1* h1(0);
330   
331   if ( x1max > x1min )
332   {
333     h1 = new TH1F(Form("%s_1",name),Form("%s_1",name),
334                   nbins,x1min,x1max);
335   }
336   
337   TList* list = AliMpManuList::ManuList();
338   TIter next(list);
339   AliMpIntPair* p;
340   Int_t n(0);
341   Int_t nPerStation[7];
342   
343   for ( Int_t i = 0; i < 7; ++i ) nPerStation[i]=0;
344   
345   while ( ( p = (AliMpIntPair*)next() ) )
346   {
347     Int_t detElemId = p->GetFirst();
348     Int_t manuId = p->GetSecond();
349     Int_t station = AliMpDEManager::GetChamberId(detElemId);
350     
351     AliMpVSegmentation* seg = 
352       segFactory()->CreateMpSegmentationByElectronics(detElemId,manuId);
353     
354     AliMUONVCalibParam* value = 
355       dynamic_cast<AliMUONVCalibParam*>(store.Get(detElemId,manuId));
356     
357     if (value)
358     {
359       for ( Int_t manuChannel = 0; manuChannel < value->Size(); ++manuChannel )
360       {
361         AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,manuChannel),kFALSE);
362         if (!pad.IsValid()) continue;
363
364         ++n;
365         ++nPerStation[station];
366         Float_t x = value->ValueAsFloat(manuChannel,0);
367         h0->Fill(x);
368         if (h1)
369         {
370           h1->Fill(value->ValueAsFloat(manuChannel,1));
371         }
372       }
373     }
374     else
375     {
376       cout << "Got a null value for DE=" << detElemId << " manuId="
377       << manuId << endl;
378     }
379   }
380   
381   cout << "Number of channels = " << n << endl;
382   for ( Int_t i = 0; i < 7; ++i )
383   {
384     cout << "Station " << (i+1) << " " << nPerStation[i] << endl;
385   }
386   
387   delete list;
388 }
389
390 //_____________________________________________________________________________
391 void plotCDB(const char* calibType, Int_t runNumber)
392 {
393   AliMUONV2DStore* store = read2D(calibType,runNumber);
394   if (!store) return;
395
396   TString c(calibType);
397   c.ReplaceAll("/","_");
398   
399   plot(*store,c.Data());
400   
401   delete store;
402 }
403
404 //_____________________________________________________________________________
405 void testReadStore(const AliMUONV2DStore& store, Int_t n)
406 {
407   TList* list = AliMpManuList::ManuList();
408   TIter next(list);
409   AliMpIntPair* p;
410   
411   while ( ( p = (AliMpIntPair*)next() ) )
412   {
413     for ( Int_t i = 0; i < n; ++i )
414     {
415       store.Get(p->GetFirst(),p->GetSecond());
416     }
417   }
418   delete list;
419
420
421 //_____________________________________________________________________________
422 Int_t makePedestalStore(AliMUONV2DStore& pedestalStore, Bool_t defaultValues)
423 {
424   TList* list = AliMpManuList::ManuList();
425   TIter next(list);
426   
427   AliMpIntPair* p;
428   
429   Int_t nchannels(0);
430   Int_t nmanus(0);
431   
432   Bool_t replace = kFALSE;
433   
434   const Int_t nChannels(64);
435   const Float_t kPedestalMeanMean(150);
436   const Float_t kPedestalMeanSigma(10);
437   const Float_t kPedestalSigmaMean(1.0);
438   const Float_t kPedestalSigmaSigma(0.2);
439   
440   while ( ( p = (AliMpIntPair*)next() ) )
441   {
442     ++nmanus;
443     AliMUONVCalibParam* ped = new AliMUONCalibParam2F(nChannels,AliMUONVCalibParam::InvalidFloatValue());
444     
445     Int_t detElemId = p->GetFirst();
446         
447     Int_t manuId = p->GetSecond();
448     
449     AliMpVSegmentation* seg = 
450       segFactory()->CreateMpSegmentationByElectronics(detElemId,manuId);
451     
452     for ( Int_t manuChannel = 0; manuChannel < nChannels; ++manuChannel )
453     {
454       AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,manuChannel),kFALSE);
455       if (!pad.IsValid()) continue;
456       
457       ++nchannels;
458       
459       Float_t meanPedestal;
460       Float_t sigmaPedestal;
461       
462       if ( defaultValues ) 
463       {
464         meanPedestal = 0.0;
465         sigmaPedestal = 1.0;
466       }
467       else
468       {
469         meanPedestal = -1;
470         while ( meanPedestal < 0 )
471         {
472           meanPedestal = gRandom->Gaus(kPedestalMeanMean,kPedestalMeanSigma);
473         }
474         sigmaPedestal = -1;
475         while ( sigmaPedestal < 0 )
476         {
477           sigmaPedestal = gRandom->Gaus(kPedestalSigmaMean,kPedestalSigmaSigma);
478         }
479       }
480       ped->SetValueAsFloat(manuChannel,0,meanPedestal);
481       ped->SetValueAsFloat(manuChannel,1,sigmaPedestal);
482       
483     }
484     Bool_t ok = pedestalStore.Set(detElemId,manuId,ped,replace);
485     if (!ok)
486     {
487       cout << "Could not set DetElemId=" << detElemId << " manuId="
488       << manuId << endl;
489     }
490   }
491   
492   delete list;
493   cout << nmanus << " Manus and " << nchannels << " channels." << endl;
494   return nchannels;
495   
496 }
497
498 //_____________________________________________________________________________
499 Int_t makeGainStore(AliMUONV2DStore& gainStore, Bool_t defaultValues)
500 {  
501   TList* list = AliMpManuList::ManuList();
502   TIter next(list);
503   
504   AliMpIntPair* p;
505   
506   Int_t nchannels(0);
507   Int_t nmanus(0);
508   
509   Bool_t replace = kFALSE;
510   
511   const Int_t nChannels(64);
512   const Double_t kSaturation(3000);
513   const Double_t kGainMean(1.0);
514   const Double_t kGainSigma(0.05);
515   
516   while ( ( p = (AliMpIntPair*)next() ) )
517   {
518     ++nmanus;
519     AliMUONVCalibParam* gain = new AliMUONCalibParam2F(nChannels,AliMUONVCalibParam::InvalidFloatValue());
520
521     Int_t detElemId = p->GetFirst();
522     Int_t manuId = p->GetSecond();
523
524     AliMpVSegmentation* seg = 
525       segFactory()->CreateMpSegmentationByElectronics(detElemId,manuId);
526
527     for ( Int_t manuChannel = 0; manuChannel < nChannels; ++manuChannel )
528     {
529       AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,manuChannel),kFALSE);
530       if (!pad.IsValid()) continue;
531       
532       ++nchannels;
533       
534       Float_t meanGain;
535       Float_t saturation(kSaturation);
536     
537       if ( defaultValues ) 
538       {
539         meanGain = 1.0;
540       }
541       else
542       {
543         meanGain = -1;
544         while ( meanGain < 0 )
545         {
546           meanGain = gRandom->Gaus(kGainMean,kGainSigma);
547         }
548       }
549       gain->SetValueAsFloat(manuChannel,0,meanGain);
550       gain->SetValueAsFloat(manuChannel,1,saturation);
551       
552     }
553     Bool_t ok = gainStore.Set(detElemId,manuId,gain,replace);
554     if (!ok)
555     {
556       cout << "Could not set DetElemId=" << detElemId << " manuId="
557         << manuId << endl;
558     }
559   }
560   
561   delete list;
562   cout << nmanus << " Manus and " << nchannels << " channels." << endl;
563   return nchannels;
564 }
565
566 //_____________________________________________________________________________
567 Int_t makeDeadStore(AliMUONV2DStore& deadStore, Bool_t defaultValues)
568 {  
569   TList* list = AliMpManuList::ManuList();
570   TIter next(list);
571   
572   AliMpIntPair* p;
573   
574   Int_t nchannels(0);
575   Int_t nmanus(0);
576   
577   Bool_t replace = kFALSE;
578   
579   const Int_t nChannels(64);
580   const Double_t deadProba = 1.0; // 1%
581   
582   while ( ( p = (AliMpIntPair*)next() ) )
583   {
584     ++nmanus;
585     AliMUONVCalibParam* dead = new AliMUONCalibParam1I(nChannels,-9999);
586     
587     Int_t detElemId = p->GetFirst();
588     Int_t manuId = p->GetSecond();
589     
590     AliMpVSegmentation* seg = 
591       segFactory()->CreateMpSegmentationByElectronics(detElemId,manuId);
592     
593     for ( Int_t manuChannel = 0; manuChannel < nChannels; ++manuChannel )
594     {
595       AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,manuChannel),kFALSE);
596       if (!pad.IsValid()) continue;
597       
598       ++nchannels;
599             
600       if (!defaultValues)
601       {
602         // probability that this channel is dead ~ 1%
603         if ( gRandom->Uniform(100.0) < deadProba ) 
604         {
605           Int_t reason = 5; // that value could be used to distinguish
606                             // why the channel is dead or how it was flagged bad (online,
607                             // offline, by chance...). 5 is of course a fake number.
608           dead->SetValueAsInt(manuChannel,0,reason);
609         }
610       }
611     }
612     Bool_t ok = deadStore.Set(detElemId,manuId,dead,replace);
613     if (!ok)
614     {
615       cout << "Could not set DetElemId=" << detElemId << " manuId="
616       << manuId << endl;
617     }
618   }
619   
620   delete list;
621   cout << nmanus << " Manus and " << nchannels << " channels." << endl;
622   return nchannels;
623 }
624
625 //_____________________________________________________________________________
626 void testMakeStores(Int_t readLoop)
627 {
628   AliMUONV2DStore* pedestalStore = new AliMUON2DMap;
629   AliMUONV2DStore* gainStore = new AliMUON2DMap;
630   AliMUONV2DStore* deadStore = new AliMUON2DMap;
631   
632   TStopwatch timer;
633   
634   cout << "Creating" << endl;
635   
636   Bool_t defaultValues = kTRUE;
637   
638   timer.Start(kTRUE);
639   makePedestalStore(*pedestalStore,defaultValues);
640   makeGainStore(*gainStore,defaultValues);
641   makeDeadStore(*deadStore,defaultValues);
642   timer.Print();
643   
644   cout << "Reading..." << endl;
645   timer.Start(kTRUE);
646   testReadStore(*pedestalStore,readLoop);
647   testReadStore(*gainStore,readLoop);
648   testReadStore(*deadStore,readLoop);
649   cout << timer.CpuTime()/readLoop << " CPUs (mean of " << readLoop 
650     <<" samples." << endl;
651   
652   delete pedestalStore;
653   delete gainStore;
654   delete deadStore;
655 }
656
657 //_____________________________________________________________________________
658 void generateTrigger(const char* cdbpath)
659 {
660   // 
661   // Generate trigger related conditions :
662   //
663   // - trigger masks for board (locals, regionals, global)
664   // - trigger lut
665   // - trigger efficiency
666   // - trigger switches (to be implemented FIXME)
667   //
668   const Int_t nlboards = 234;
669   AliMUONV1DStore* localBoardMasks = new AliMUON1DArray(nlboards+1);
670   
671   // Generate fake mask values for 234 localboards and put that into
672   // one single container (localBoardMasks)
673   for ( Int_t i = 1; i <= nlboards; ++i )
674   {
675     AliMUONVCalibParam* localBoard = new AliMUONCalibParam1I(8);
676     for ( Int_t x = 0; x < 2; ++x )
677     {
678       for ( Int_t y = 0; y < 4; ++y )
679       {
680         Int_t index = x*4+y;
681         localBoard->SetValueAsInt(index,0,0xFFFF);
682       }
683     }
684     localBoardMasks->Set(i,localBoard,kFALSE);
685   }
686   
687   // Generate values for regional boards
688   const Int_t nrboards = 16;
689   AliMUONV1DStore* regionalBoardMasks = new AliMUON1DArray(16);
690   
691   for ( Int_t i = 0; i < nrboards; ++i )
692   {
693     AliMUONVCalibParam* regionalBoard = new AliMUONCalibParam1I(16);
694     for ( Int_t j = 0; j < 16; ++j )
695     {
696       regionalBoard->SetValueAsInt(j,0,0x3F);
697     }
698     regionalBoardMasks->Set(i,regionalBoard,kFALSE);
699   }
700   
701   // Generate values for global board
702   AliMUONVCalibParam* globalBoardMasks = new AliMUONCalibParam1I(16);
703   
704   for ( Int_t j = 0; j < 16; ++j )
705   {
706     globalBoardMasks->SetValueAsInt(j,0,0xFFF);
707   }
708     
709   AliMUONTriggerLut lut;
710   lut.ReadFromFile("$(ALICE_ROOT)/MUON/data/lutAptLpt1Hpt1p7.root");
711   
712   AliMUONTriggerEfficiencyCells cells("$ALICE_ROOT/MUON/data/efficiencyCells.dat");
713   
714   //--------------------------------------------
715   // Store the resulting containers into the CDB
716   Int_t ever = 99999999;
717   
718   AliCDBId id("MUON/Calib/LocalTriggerBoardMasks",0,ever);
719   
720   AliCDBMetaData md;
721   md.SetBeamPeriod(1);
722   md.SetAliRootVersion(gROOT->GetVersion());
723   md.SetComment("Test with default values");
724   md.SetResponsible("Rachid Guernane");
725   
726   AliCDBManager* man = AliCDBManager::Instance();
727   man->SetDefaultStorage(cdbpath);
728   man->Put(localBoardMasks,id,&md);
729   
730   id.SetPath("MUON/Calib/RegionalTriggerBoardMasks");
731   
732   man->Put(regionalBoardMasks,id,(AliCDBMetaData*)md.Clone());
733   
734   id.SetPath("MUON/Calib/GlobalTriggerBoardMasks");
735   
736   man->Put(globalBoardMasks,id,(AliCDBMetaData*)md.Clone());
737
738   id.SetPath("MUON/Calib/TriggerLut");
739   
740   man->Put(&lut,id,(AliCDBMetaData*)md.Clone());
741   
742   id.SetPath("MUON/Calib/TriggerEfficiency");
743   md.SetResponsible("Diego Stocco");
744   
745   man->Put(&cells,id,(AliCDBMetaData*)md.Clone());
746   
747   delete localBoardMasks;
748   delete regionalBoardMasks;
749   delete globalBoardMasks;
750 }
751
752
753 //_____________________________________________________________________________
754 void writeToCDB(const char* cdbpath, const char* calibpath, TObject* object, 
755                 Int_t startRun, Int_t endRun, Bool_t defaultValues)
756 {
757   AliCDBId id(calibpath,startRun,endRun);
758   AliCDBMetaData md;
759   md.SetBeamPeriod(1);
760   md.SetAliRootVersion(gROOT->GetVersion());
761   if ( defaultValues )
762   {
763     md.SetComment("Test with default values");
764   }
765   else
766   {
767     md.SetComment("Test with random values");
768   }
769   md.SetResponsible("Laurent Aphecetche");
770   
771   AliCDBManager* man = AliCDBManager::Instance();
772   man->SetDefaultStorage(cdbpath);
773   man->Put(object,id,&md);
774 }
775
776 //_____________________________________________________________________________
777 void writePedestals(const char* cdbpath, Bool_t defaultValues,
778                     Int_t startRun, Int_t endRun)
779 {
780   AliMUONV2DStore* pedestalStore = new AliMUON2DMap;
781   Int_t ngenerated = makePedestalStore(*pedestalStore,defaultValues);
782   cout << "Ngenerated = " << ngenerated << endl;
783   
784   writeToCDB(cdbpath,"MUON/Calib/Pedestals",pedestalStore,startRun,endRun,defaultValues);
785 }
786
787 //_____________________________________________________________________________
788 void writeGains(const char* cdbpath, Bool_t defaultValues,
789                     Int_t startRun, Int_t endRun)
790 {
791   AliMUONV2DStore* gainStore = new AliMUON2DMap;
792   Int_t ngenerated = makeGainStore(*gainStore,defaultValues);
793   cout << "Ngenerated = " << ngenerated << endl;
794   
795   writeToCDB(cdbpath,"MUON/Calib/Gains",gainStore,startRun,endRun,defaultValues);
796 }
797
798 //_____________________________________________________________________________
799 void writeDeadChannels(const char* cdbpath, Bool_t defaultValues,
800                        Int_t startRun, Int_t endRun)
801 {
802   AliMUONV2DStore* deadStore = new AliMUON2DMap;
803   Int_t ngenerated = makeDeadStore(*deadStore,defaultValues);
804   cout << "Ngenerated = " << ngenerated << endl;
805   
806   writeToCDB(cdbpath,"MUON/Calib/DeadChannels",deadStore,startRun,endRun,defaultValues);
807 }
808
809