Moving from AliMUONCalibParam1I to AliMUONCalibParamNI
[u/mrichter/AliRoot.git] / MUON / AliMUONPadStatusMaker.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 /// \class AliMUONPadStatusMaker
19 ///
20 /// Make a 2DStore of pad statuses, using different sources of information,
21 /// like pedestal values, gain values, and HV values.
22 ///
23 /// \author Laurent Aphecetche
24
25 #include "AliMUONPadStatusMaker.h"
26
27 #include "AliMUON2DMap.h"
28 #include "AliMUON2DStoreValidator.h"
29 #include "AliMUONCalibParamNI.h"
30 #include "AliMUONCalibrationData.h"
31 #include "AliMUONHVNamer.h"
32 #include "AliMUONObjectPair.h"
33 #include "AliMUONVCalibParam.h"
34 #include "AliMUONVDataIterator.h"
35 #include "AliMpArea.h"
36 #include "AliMpDEIterator.h"
37 #include "AliMpDEManager.h"
38 #include "AliMpIntPair.h"
39 #include "AliMpManuList.h"
40 #include "AliMpMotifMap.h"
41 #include "AliMpMotifPosition.h"
42 #include "AliMpPCB.h"
43 #include "AliMpPad.h"
44 #include "AliMpSector.h"
45 #include "AliMpSectorSegmentation.h"
46 #include "AliMpSegmentation.h"
47 #include "AliMpSlat.h"
48 #include "AliMpSlatSegmentation.h"
49 #include "AliMpStationType.h"
50 #include "AliMpVPadIterator.h"
51 #include "AliCDBEntry.h"
52 #include "AliCDBManager.h"
53 #include "AliDCSValue.h"
54 #include "AliLog.h"
55 #include <Riostream.h>
56 #include <TMap.h>
57 #include <TStopwatch.h>
58 #include <TString.h>
59
60
61 /// \cond CLASSIMP
62 ClassImp(AliMUONPadStatusMaker)
63 /// \endcond
64
65 //_____________________________________________________________________________
66 AliMUONPadStatusMaker::AliMUONPadStatusMaker(const AliMUONCalibrationData& calibData)
67 : fCalibrationData(calibData),
68   fPedMeanLimits(0,4095),
69   fPedSigmaLimits(0,4095),
70   fHVSt12Limits(0,5000),
71   fHVSt345Limits(0,5000)
72 {
73    /// ctor
74 }
75
76 //_____________________________________________________________________________
77 AliMUONPadStatusMaker::~AliMUONPadStatusMaker()
78 {
79   /// dtor.
80 }
81
82 //_____________________________________________________________________________
83 AliMUONV2DStore*
84 AliMUONPadStatusMaker::Combine(const AliMUONV2DStore& store1,
85                                const AliMUONV2DStore& store2,
86                                Int_t binShift) const
87 {
88   /// Combine two status containers into one, shifting store2 status bits
89   /// to the left by binShift before making an OR with store1.
90   
91   TStopwatch timer;
92   timer.Start(kTRUE);
93   
94   AliMUONV2DStore* combined = static_cast<AliMUONV2DStore*>(store1.Clone());
95   
96   AliMUONVDataIterator* it = store1.Iterator();
97   AliMUONObjectPair* pair;
98   
99   while ( ( pair = static_cast<AliMUONObjectPair*>(it->Next()) ) )
100   {
101     AliMpIntPair* ip = static_cast<AliMpIntPair*>(pair->First());
102     Int_t detElemId = ip->GetFirst();
103     Int_t manuId = ip->GetSecond();
104     AliMUONVCalibParam* param1 = static_cast<AliMUONVCalibParam*>(store1.Get(detElemId,manuId));
105     AliMUONVCalibParam* param2 = static_cast<AliMUONVCalibParam*>(store2.Get(detElemId,manuId));
106     if (!param2)
107     {
108       AliWarning(Form("Could not get statuses for store2 for DE %d ManuId %d. Marking as missing.",
109                     detElemId,manuId));
110       param2 = static_cast<AliMUONVCalibParam*>(param1->Clone());
111       for ( Int_t manuChannel = 0; manuChannel < param2->Size(); ++manuChannel )
112       {
113         param2->SetValueAsInt(manuChannel,0,kMissing);
114       }      
115     }
116     AliMUONVCalibParam* paramCombined = static_cast<AliMUONVCalibParam*>(combined->Get(detElemId,manuId));
117     if (!paramCombined)
118     {
119       paramCombined = static_cast<AliMUONVCalibParam*>(param2->Clone());
120       combined->Set(detElemId,manuId,paramCombined,kFALSE);
121     }
122     
123     for ( Int_t manuChannel = 0; manuChannel < param1->Size(); ++manuChannel )
124     {
125       if ( AliMpManuList::DoesChannelExist(detElemId, manuId, manuChannel) )
126       {
127         Int_t status1(param1->ValueAsInt(manuChannel));
128         Int_t status2(param2->ValueAsInt(manuChannel));
129         
130         Int_t status = status1 | (status2 << binShift);
131         
132         paramCombined->SetValueAsInt(manuChannel,0,status);
133       }
134     }
135   }
136   
137   delete it;
138   
139   AliInfo("Timer:");
140   StdoutToAliInfo(timer.Print(););
141   
142   return combined;
143 }
144
145 //_____________________________________________________________________________
146 AliMUONV2DStore* 
147 AliMUONPadStatusMaker::GeneratePadStatus(Int_t value)
148 {
149   /// Generate a "fake" store, with all (detElemId,manuId) present,
150   /// and containing all the same value
151   
152   AliMUONCalibParamNI param(1,64,value);
153   
154   return AliMUON2DMap::Generate(param);
155 }
156
157 //_____________________________________________________________________________
158 Bool_t 
159 AliMUONPadStatusMaker::GetSt12Status(const TMap& hvMap,
160                                      Int_t detElemId, Int_t sector,
161                                      Bool_t& hvChannelTooLow,
162                                      Bool_t& hvChannelTooHigh,
163                                      Bool_t& hvChannelON) const
164 {
165   /// Get HV status for one HV sector of St12
166   
167   /// For a given PCB in a given DE, get the HV status (both the channel
168   /// and the switch).
169   /// Returns false if hv switch changed during the run.
170   
171   Bool_t error = kFALSE;
172   hvChannelTooLow = kFALSE;
173   hvChannelTooHigh = kFALSE;
174   hvChannelON = kTRUE;
175   
176   AliMUONHVNamer hvNamer;
177   
178   TString hvChannel(hvNamer.DCSHVChannelName(detElemId,sector));
179   
180   TPair* hvPair = static_cast<TPair*>(hvMap.FindObject(hvChannel.Data()));
181   if (!hvPair)
182   {
183     AliError(Form("Did not find expected alias (%s) for DE %d",
184                   hvChannel.Data(),detElemId));  
185     error = kTRUE;
186   }
187   else
188   {
189     TObjArray* values = static_cast<TObjArray*>(hvPair->Value());
190     if (!values)
191     {
192       AliError(Form("Could not get values for alias %s",hvChannel.Data()));
193       error = kTRUE;
194     }
195     else
196     {
197       // find out min and max value, and makes a cut
198       Float_t hvMin(1E9);
199       Float_t hvMax(0);
200       TIter next(values);
201       AliDCSValue* val;
202       
203       while ( ( val = static_cast<AliDCSValue*>(next()) ) )
204       {
205         Float_t hv = val->GetFloat();
206         hvMin = TMath::Min(hv,hvMin);
207         hvMax = TMath::Max(hv,hvMax);
208       }
209       
210       float lowThreshold = fHVSt12Limits.X();
211       float highThreshold = fHVSt12Limits.Y();
212             
213       if ( hvMin < lowThreshold ) hvChannelTooLow = kTRUE;
214       if ( hvMax > highThreshold ) hvChannelTooHigh = kTRUE;
215       if ( hvMin < 1 ) hvChannelON = kFALSE;
216     }
217   }
218   
219   return error;
220 }
221
222 //_____________________________________________________________________________
223 Bool_t 
224 AliMUONPadStatusMaker::GetSt345Status(const TMap& hvMap,
225                                       Int_t detElemId, Int_t pcbIndex,
226                                       Bool_t& hvChannelTooLow,
227                                       Bool_t& hvChannelTooHigh,
228                                       Bool_t& hvChannelON,
229                                       Bool_t& hvSwitchON) const
230 {
231   /// For a given PCB in a given DE, get the HV status (both the channel
232   /// and the switch).
233   /// Returns false if something goes wrong (in particular if 
234   /// hv switch changed during the run).
235   
236   Bool_t error = kFALSE;
237   hvChannelTooLow = kFALSE;
238   hvChannelTooHigh = kFALSE;
239   hvSwitchON = kTRUE;
240   hvChannelON = kTRUE;
241   
242   AliMUONHVNamer hvNamer;
243   
244   TString hvChannel(hvNamer.DCSHVChannelName(detElemId));
245   
246   TPair* hvPair = static_cast<TPair*>(hvMap.FindObject(hvChannel.Data()));
247   if (!hvPair)
248   {
249     AliError(Form("Did not find expected alias (%s) for DE %d",
250                   hvChannel.Data(),detElemId));  
251     error = kTRUE;
252   }
253   else
254   {
255     TObjArray* values = static_cast<TObjArray*>(hvPair->Value());
256     if (!values)
257     {
258       AliError(Form("Could not get values for alias %s",hvChannel.Data()));
259       error = kTRUE;
260     }
261     else
262     {
263       // find out min and max value, and makes a cut
264       Float_t hvMin(1E9);
265       Float_t hvMax(0);
266       TIter next(values);
267       AliDCSValue* val;
268       
269       while ( ( val = static_cast<AliDCSValue*>(next()) ) )
270       {
271         Float_t hv = val->GetFloat();
272         hvMin = TMath::Min(hv,hvMin);
273         hvMax = TMath::Max(hv,hvMax);
274       }
275
276       float lowThreshold = fHVSt345Limits.X();
277       float highThreshold = fHVSt345Limits.Y();
278
279       if ( hvMin < lowThreshold ) hvChannelTooLow = kTRUE;
280       if ( hvMax > highThreshold ) hvChannelTooHigh = kTRUE;
281       if ( hvMin < 1 ) hvChannelON = kFALSE;
282     }
283   }
284   
285   TString hvSwitch(hvNamer.DCSHVSwitchName(detElemId,pcbIndex));
286   TPair* switchPair = static_cast<TPair*>(hvMap.FindObject(hvSwitch.Data()));
287   if (!switchPair)
288   {
289     AliError(Form("Did not find expected alias (%s) for DE %d PCB %d",
290                   hvSwitch.Data(),detElemId,pcbIndex));
291     error = kTRUE;
292   }
293   else
294   {
295     TObjArray* values = static_cast<TObjArray*>(switchPair->Value());
296     if (!values)
297     {    
298       AliError(Form("Could not get values for alias %s",hvSwitch.Data()));
299       error = kTRUE;
300     }
301     else
302     {
303       // we'll count the number of ON/OFF for this pad, to insure
304       // consistency (i.e. if status changed during the run, we should
305       // at least notify this fact ;-) and hope it's not the norm)
306       Int_t nTrue(0);
307       Int_t nFalse(0);
308       TIter next(values);
309       AliDCSValue* val;
310       
311       while ( ( val = static_cast<AliDCSValue*>(next()) ) )
312       {
313         if ( val->GetBool() )
314         {
315           ++nTrue;
316         }
317         else
318         {
319           ++nFalse;
320         }
321       }
322       
323       if ( (nTrue>0 && nFalse>0) )
324       {
325         AliWarning(Form("Status of HV Switch %s changed during this run nTrue=%d nFalse=%d! Will consider it OFF",
326                         hvSwitch.Data(),nTrue,nFalse));
327         error = kTRUE;
328       }
329       
330       if ( nFalse ) hvSwitchON = kFALSE;
331     }
332   }
333   return error;
334 }
335
336 //_____________________________________________________________________________
337 AliMUONV2DStore* 
338 AliMUONPadStatusMaker::MakeGainStatus(const AliMUONV2DStore& /*gainValues*/) const
339 {
340   /// FIXME: to be implemented
341   AliWarning("Not implemented yet");
342   return 0x0;
343 }
344
345 //_____________________________________________________________________________
346 AliMUONV2DStore* 
347 AliMUONPadStatusMaker::MakeHVStatus(const TMap& hvValues) const
348 {
349   /// Scrutinize HV values and deduce an HV status for each pad
350   
351   TStopwatch timerSt12;
352   TStopwatch timerSt345;
353   
354   timerSt12.Start(kTRUE);
355   timerSt12.Stop();
356   timerSt345.Start(kTRUE);
357   timerSt345.Stop();
358   
359   AliMUONHVNamer hvNamer;
360   
361   AliMpDEIterator deIt;
362   
363   deIt.First();
364   
365   AliMUONV2DStore* hv = new AliMUON2DMap(kTRUE);
366   
367   while ( !deIt.IsDone() )
368   {
369     Int_t detElemId = deIt.CurrentDEId();
370     
371     switch ( AliMpDEManager::GetStationType(detElemId) )
372     {
373       case AliMp::kStation1:
374       case AliMp::kStation2:
375         timerSt12.Start(kFALSE);
376         for ( int sector = 0; sector < 3; ++sector)
377         {
378           AliDebug(1,Form("detElemId %5d sector %d",detElemId,sector));
379
380           Bool_t hvChannelTooLow, hvChannelTooHigh, hvChannelON;
381           Bool_t error = GetSt12Status(hvValues,
382                                        detElemId,sector,
383                                        hvChannelTooLow,hvChannelTooHigh,
384                                        hvChannelON);
385           Int_t status = 0;
386           if ( error ) status |= kHVError;
387           if ( hvChannelTooLow ) status |= kHVTooLow;
388           if ( hvChannelTooHigh ) status |= kHVTooHigh; 
389           if ( !hvChannelON ) status |= kHVChannelOFF;
390           SetStatusSt12(*hv,detElemId,sector,status);
391           
392         }
393           timerSt12.Stop();
394         break;
395       case AliMp::kStation345:
396       {
397         timerSt345.Start(kFALSE);
398         for ( Int_t pcbIndex = 0; pcbIndex < hvNamer.NumberOfPCBs(detElemId); ++pcbIndex)
399         {
400           AliDebug(1,Form("detElemId %5d pcbIndex %d",detElemId,pcbIndex));
401           Bool_t hvChannelTooLow, hvChannelTooHigh, hvChannelON,hvSwitchON;
402           Bool_t error = GetSt345Status(hvValues,
403                                         detElemId,pcbIndex,
404                                         hvChannelTooLow,hvChannelTooHigh,
405                                         hvChannelON,hvSwitchON);
406           Int_t status = 0;
407           if ( error ) status |= kHVError;
408           if ( hvChannelTooLow ) status |= kHVTooLow;
409           if ( hvChannelTooHigh ) status |= kHVTooHigh; 
410           if ( !hvSwitchON ) status |= kHVSwitchOFF; 
411           if ( !hvChannelON) status |= kHVChannelOFF;
412           SetStatusSt345(*hv,detElemId,pcbIndex,status);
413         }
414         timerSt345.Stop();
415       }
416         break;
417       default:
418         break;
419     }
420     deIt.Next();
421   }
422   
423   AliInfo("St12 timer:");
424   StdoutToAliInfo(timerSt12.Print(););
425   AliInfo("St345 timer:");
426   StdoutToAliInfo(timerSt345.Print(););
427   
428   return hv;
429 }
430
431 //_____________________________________________________________________________
432 AliMUONV2DStore* 
433 AliMUONPadStatusMaker::MakePedestalStatus(const AliMUONV2DStore& pedValues) const
434 {
435   /// Assign a pedestal status to each pad
436   
437   TStopwatch timer;
438   
439   timer.Start(kTRUE);
440   
441   AliMUONV2DStore* pedStatuses = new AliMUON2DMap(kTRUE);
442   
443   AliMUONVDataIterator* it = pedValues.Iterator();
444   AliMUONObjectPair* pair;
445   Int_t nofManus(0);
446   
447   while ( ( pair = static_cast<AliMUONObjectPair*>(it->Next() ) ) )
448   {
449     AliMpIntPair* ip = static_cast<AliMpIntPair*>(pair->First());
450     Int_t detElemId = ip->GetFirst();
451     Int_t manuId = ip->GetSecond();
452     AliMUONVCalibParam* pedestals = static_cast<AliMUONVCalibParam*>(pair->Second());
453     ++nofManus;
454     for ( Int_t manuChannel = 0; manuChannel < pedestals->Size(); ++manuChannel )
455     {
456       Int_t status(0);
457       if ( AliMpManuList::DoesChannelExist(detElemId, manuId, manuChannel) )
458       {
459         Float_t pedMean = pedestals->ValueAsFloat(manuChannel,0);
460         Float_t pedSigma = pedestals->ValueAsFloat(manuChannel,1);
461         if ( pedMean < fPedMeanLimits.X() ) status |= kPedMeanTooLow;
462         if ( pedMean > fPedMeanLimits.Y() ) status |= kPedMeanTooHigh;
463         if ( pedSigma < fPedSigmaLimits.X() ) status |= kPedSigmaTooLow;
464         if ( pedSigma > fPedSigmaLimits.Y() ) status |= kPedSigmaTooHigh;
465         if ( pedMean == 0 ) status |= kPedMeanZero;
466         
467         AliMUONVCalibParam* vStatus = 
468           static_cast<AliMUONVCalibParam*>(pedStatuses->Get(detElemId,manuId));
469         if ( !vStatus ) 
470         {
471           vStatus = new AliMUONCalibParamNI(1,64,0);
472           pedStatuses->Set(detElemId,manuId,vStatus,false);
473         }
474         vStatus->SetValueAsInt(manuChannel,0,status);
475       }
476     }
477   }
478   
479   AliInfo(Form("%d manus checked in :",nofManus));
480   StdoutToAliInfo(timer.Print(););
481   return pedStatuses;  
482 }
483
484 //_____________________________________________________________________________
485 AliMUONV2DStore* 
486 AliMUONPadStatusMaker::MakeStatus() const
487 {
488   /// Read ped, gains and hv values from CDB, apply some Q&A and produces
489   /// a combined status for each pad.
490
491   TMap* hvValues = fCalibrationData.HV();
492   AliMUONV2DStore* hvStatus(0x0);
493   
494   if (!hvValues)
495   {
496     AliError("Could not get HV values from CDB. Will create dummy ones and mark those as missing");
497     AliMUONCalibParamNI param(1,64,kHVMissing);
498     hvStatus = AliMUON2DMap::Generate(param);
499   }
500   else
501   {
502     hvStatus = MakeHVStatus(*hvValues);
503   }
504   
505   AliMUONV2DStore* pedValues = fCalibrationData.Pedestals();
506   AliMUONV2DStore* pedStatus(0x0);
507
508   if (!pedValues)
509   {
510     AliError("Could not get pedestals values from CDB. Will create dummy ones and mark those as missing");
511     AliMUONCalibParamNI param(1,64,kPedMissing);
512     pedStatus = AliMUON2DMap::Generate(param);
513   }
514   else
515   {
516     pedStatus = MakePedestalStatus(*pedValues);
517   }
518   
519   // FIXME: should do the same for gains as for hv and ped.    
520   
521   AliMUONV2DStore* status = Combine(*hvStatus,*pedStatus,8);
522   
523   delete hvStatus;
524   delete pedStatus;
525   
526   // Insure we get all channels there (some or even all can be bad, but they
527   // must be there somehow).
528   
529   AliMUON2DStoreValidator validator;
530         
531   TObjArray* a = validator.Validate(*status);
532     
533   if (a) 
534   {
535     // this should not happen.
536     AliError("Status store not complete. Crash to follow soon...");
537     StdoutToAliError(a->Print(););
538     AliFatal("this should not happen at all!");
539     delete status;
540     status = 0x0;
541   }
542     
543   return status;
544 }
545
546 //_____________________________________________________________________________
547 void
548 AliMUONPadStatusMaker::SetStatusSt12(AliMUONV2DStore& hvStatus,
549                                      Int_t detElemId, 
550                                      Int_t isector,
551                                      Int_t status) const
552 {
553   /// Flag all pads of detElemId (for St12) as bad.
554   
555   // FIXME: need a way to iterator on pads over a given HV sector for St12... 
556   // we currently suppose that one sector is about a third of the chamber...
557   // FIXME !! This has to be checked very carefully...
558   
559   const AliMp::CathodType kCathodes[] = { AliMp::kCath0, AliMp::kCath1 };
560   
561   for ( Int_t icathode = 0; icathode < 2; ++icathode )
562   {
563     const AliMpSectorSegmentation* seg = 
564     static_cast<const AliMpSectorSegmentation*>(AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,kCathodes[icathode]));
565     const AliMpSector* sector = seg->GetSector();
566     AliMpMotifMap* mMap = sector->GetMotifMap();
567     TArrayI a;
568     
569     mMap->GetAllMotifPositionsIDs(a);
570     
571     TVector2 dim = seg->Dimensions();
572     Double_t x = dim.X()*2;
573     Double_t xmin = isector*x/3.0;
574     Double_t xmax = xmin + x/3.0;   
575     
576     for ( Int_t i = 0; i < a.GetSize(); ++i ) 
577     {
578       AliMpMotifPosition* pos = mMap->FindMotifPosition(a[i]);
579       Int_t manuId = pos->GetID();
580       TVector2 position = pos->Position();
581       if ( position.X() >= xmin && position.X() <= xmax) 
582       {
583         AliMUONVCalibParam* dead =
584         static_cast<AliMUONVCalibParam*>(hvStatus.Get(detElemId,manuId));
585         if (!dead)
586         {
587           dead = new AliMUONCalibParamNI(1,64,status);
588           hvStatus.Set(detElemId,manuId,dead,false);
589         }        
590         else
591         {
592           // FIXME: this should really not happen, if we'd know really the
593           // relationship between manuId and HV sector...
594           // For the time being, let's leave it like that, for testing
595           // purposes only. For production, this will have to be fixed.
596           AliWarning("Please fixme.");
597         }
598       }
599     }
600   }  
601 }
602
603 //_____________________________________________________________________________
604 void
605 AliMUONPadStatusMaker::SetStatusSt345(AliMUONV2DStore& hvStatus,
606                                       Int_t detElemId, Int_t pcbIndex,
607                                       Int_t status) const
608 {
609   /// Flag all pads of pcbIndex-th PCB of detElemId (for St345) as bad.
610   
611   const AliMp::CathodType kCathodes[] = { AliMp::kCath0, AliMp::kCath1 };
612   
613   for ( Int_t icathode = 0; icathode < 2; ++icathode )
614   {
615     const AliMpSlatSegmentation* seg = static_cast<const AliMpSlatSegmentation*>
616     (AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,kCathodes[icathode]));
617     const AliMpSlat* slat = seg->Slat();
618     const AliMpPCB* pcb = slat->GetPCB(pcbIndex);
619
620     for ( Int_t i = 0; i < pcb->GetSize(); ++i ) 
621     {
622       AliMpMotifPosition* pos = pcb->GetMotifPosition(i);
623       Int_t manuId = pos->GetID();
624       AliMUONVCalibParam* dead = 
625         static_cast<AliMUONVCalibParam*>(hvStatus.Get(detElemId,manuId));
626       if (dead)
627       {
628         AliError(Form("dead is not null as expected from DE %d manuId %d",
629                       detElemId,manuId));
630       }
631       if (!dead)
632       {
633         dead = new AliMUONCalibParamNI(1,64,status);
634         hvStatus.Set(detElemId,manuId,dead,false);
635       }
636     }    
637   }
638 }
639
640
641