Overlaps corrected, new shape of sectors
[u/mrichter/AliRoot.git] / FMD / FMDutil / AliFMDBaseDA.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    AliFMDBaseDA.cxx
17     @author  Hans Hjersing Dalsgaard <canute@nbi.dk>
18     @date    Wed Mar 26 11:30:45 2008
19     @brief   Base class for detector algorithms.
20 */
21 //
22 // This is the implementation of the (virtual) base class for the FMD
23 // detector algorithms(DA). It implements the creation of the relevant
24 // containers and handles the loop over the raw data. The derived
25 // classes can control the parameters and action to be taken making
26 // this the base class for the Pedestal, Gain and Physics DA.
27 //
28
29 #include "AliFMDBaseDA.h"
30 #include "AliRawReader.h"
31 #include "AliFMDDigit.h"
32 #include "AliFMDParameters.h"
33 #include "AliFMDRawReader.h"
34 #include "AliFMDCalibSampleRate.h"
35 #include "AliFMDCalibStripRange.h"
36 #include "AliLog.h"
37 #include "AliRawEventHeaderBase.h"
38 #include "AliFMDDigit.h"
39 #include <TClonesArray.h>
40 #include <TFile.h>
41 #include <TDatime.h>
42 #include <TSystem.h>
43 #include <TH2F.h>
44 #include <TStopwatch.h>
45 #include <TROOT.h>
46 #include <TPluginManager.h>
47 #include <iostream>
48 #include <fstream>
49
50 //_____________________________________________________________________
51 ClassImp(AliFMDBaseDA)
52 #if 0 
53 ; // Do not delete  - to let Emacs for mat the code
54 #endif
55
56 //_____________________________________________________________________
57 TString
58 AliFMDBaseDA::GetStripPath(UShort_t det, 
59                            Char_t   ring, 
60                            UShort_t sec, 
61                            UShort_t str, 
62                            Bool_t   full) const
63 {
64   // Get the strip path 
65   // 
66   // Parameters 
67   //     det      Detector number
68   //     ring     Ring identifier 
69   //     sec      Sector number 
70   //     str      Strip number
71   //     full     If true, return full path 
72   // 
73   // Return 
74   //     The path
75   return TString::Format("%s%sFMD%d%c[%02d,%03d]", 
76                          (full ? GetSectorPath(det, ring, sec, full).Data() : ""), 
77                          (full ? "/" : ""), det, ring, sec, str);
78 }
79 //_____________________________________________________________________
80 TString
81 AliFMDBaseDA::GetSectorPath(UShort_t det, 
82                             Char_t   ring, 
83                             UShort_t sec, 
84                             Bool_t   full) const
85 {
86   // Get the strip path 
87   // 
88   // Parameters 
89   //     det      Detector number
90   //     ring     Ring identifier 
91   //     sec      Sector number 
92   //     str      Strip number
93   //     full     If true, return full path 
94   // 
95   // Return 
96   //     The path
97   return TString::Format("%s%sFMD%d%c[%02d]", 
98                          (full ? GetRingPath(det, ring, full).Data() : ""), 
99                          (full ? "/" : ""), det, ring, sec);
100 }
101 //_____________________________________________________________________
102 TString
103 AliFMDBaseDA::GetRingPath(UShort_t det, 
104                           Char_t   ring, 
105                           Bool_t   full) const
106 {
107   // Get the strip path 
108   // 
109   // Parameters 
110   //     det      Detector number
111   //     ring     Ring identifier 
112   //     sec      Sector number 
113   //     str      Strip number
114   //     full     If true, return full path 
115   // 
116   // Return 
117   //     The path
118   return TString::Format("%s%sFMD%d%c", 
119                          (full ? GetDetectorPath(det, full).Data() : ""), 
120                          (full ? "/" : ""), det, ring);
121 }
122 //_____________________________________________________________________
123 TString
124 AliFMDBaseDA::GetDetectorPath(UShort_t det, 
125                               Bool_t   full) const
126 {
127   // Get the strip path 
128   // 
129   // Parameters 
130   //     det      Detector number
131   //     ring     Rg identifier 
132   //     sec      Sector number 
133   //     str      Strip number
134   //     full     If true, return full path 
135   // 
136   // Return 
137   //     The path
138   return TString::Format("%s%sFMD%d", 
139                          (full ? fDiagnosticsFilename.Data() : ""), 
140                          (full ? ":/" : ""), det);
141 }
142
143 //_____________________________________________________________________
144 AliFMDBaseDA::AliFMDBaseDA() : 
145   TNamed(),
146   fDiagnosticsFilename("diagnosticsHistograms.root"),
147   fOutputFile(),
148   fConditionsFile(),
149   fSaveHistograms(kFALSE),
150   fMakeSummaries(kFALSE),
151   fDetectorArray(),
152   fPulseSize(10),
153   fPulseLength(10),
154   fRequiredEvents(0),
155   fCurrentEvent(0), 
156   fRunno(0),
157   fSummaries(0),
158   fAll(false)
159 {
160   //Constructor
161   for(Int_t i = 0; i< 3;i++) {
162     fSeenDetectors[i] = false;
163     fNEventsPerDetector[i] = 0;
164   }
165 }
166
167
168 //_____________________________________________________________________
169 AliFMDBaseDA::AliFMDBaseDA(const AliFMDBaseDA & baseDA) : 
170   TNamed(baseDA),
171   fDiagnosticsFilename(baseDA.fDiagnosticsFilename.Data()),
172   fOutputFile(),
173   fConditionsFile(),
174   fSaveHistograms(baseDA.fSaveHistograms),
175   fMakeSummaries(baseDA.fMakeSummaries),
176   fDetectorArray(baseDA.fDetectorArray),
177   fPulseSize(baseDA.fPulseSize),
178   fPulseLength(baseDA.fPulseLength),
179   fRequiredEvents(baseDA.fRequiredEvents),
180   fCurrentEvent(baseDA.fCurrentEvent),
181   fRunno(baseDA.fRunno),
182   fSummaries(0),
183   fAll(baseDA.fAll)
184 {
185   //Copy constructor
186   for(Int_t i = 0; i< 3;i++) {
187     fSeenDetectors[i] = baseDA.fSeenDetectors[0];
188     fNEventsPerDetector[i] = baseDA.fNEventsPerDetector[i];
189   }
190
191   fDetectorArray.SetOwner();  
192   fDetectorArray.ResetBit(TObject::kMustCleanup);
193   fSummaries.ResetBit(TObject::kMustCleanup);
194   fDetectorArray.ResetBit(TObject::kIsOnHeap);
195   fSummaries.ResetBit(TObject::kIsOnHeap);
196 }
197
198
199 //_____________________________________________________________________
200 AliFMDBaseDA::~AliFMDBaseDA() 
201 {
202   //destructor
203   // fDetectorArray.ls();
204   // fSummaries.ls();
205   // fDetectorArray.SetOwner(false);
206   // fSummaries.SetOwner(false);
207
208 }
209
210 //_____________________________________________________________________
211 Bool_t
212 AliFMDBaseDA::OpenFiles(Bool_t appendRun)
213 {
214   fDetectorArray.SetOwner();
215   if (!appendRun || fRunno == 0) {
216     Rotate("conditions.csv", 3);
217     fConditionsFile.open("conditions.csv");
218   }
219   else {
220     fConditionsFile.open(Form("conditions_%09d.csv", 
221                               fRunno));
222   }
223   if (!fConditionsFile) { 
224     Error("OpenFiles", "Failed to open conditions file");
225     return false;
226   }
227   return true;
228 }
229
230 //_____________________________________________________________________
231 Bool_t AliFMDBaseDA::HaveEnough(Int_t nEvents) const
232 {
233   // if (!fAll) return nEvents > GetRequiredEvents();
234   if (nEvents <= 1) return false;
235
236   Bool_t ret = true; // Assume we have it 
237   for (Int_t i = 0; i < 3; i++) { 
238     if (!fSeenDetectors[i]) continue;
239     if (Int_t(fNEventsPerDetector[i]) <= GetRequiredEvents()) ret = false;
240   }
241   return ret;
242 }
243 //_____________________________________________________________________
244 UShort_t AliFMDBaseDA::GetProgress(Int_t nEvents) const
245 {
246   // if (!fAll) 
247   //  return UShort_t((nEvents *100)/ GetRequiredEvents());
248
249   if (nEvents <= 1) return 0;
250
251   Int_t got = 0;
252   Int_t total = 0;
253   for (Int_t i = 0; i < 3; i++) {
254     if (!fSeenDetectors[i]) continue;
255     got   += fNEventsPerDetector[i];
256     total += GetRequiredEvents();
257   }
258   return UShort_t(total > 0 ? (got * 100.) / total : 0);
259 }
260
261 //_____________________________________________________________________
262 Bool_t
263 AliFMDBaseDA::Run(AliRawReader* reader, Bool_t appendRun, Bool_t isBase) 
264 {
265   //Run the FMD DA
266   TFile* diagFile = 0;
267   // if (fSaveHistograms)
268   //  diagFile = TFile::Open(fDiagnosticsFilename.Data(),"RECREATE");
269
270   reader->Reset();
271   fRunno = reader->GetRunNumber();
272   if (fRunno <= 0) { 
273     TString dateRunNumber(gSystem->Getenv("DATE_RUN_NUMBER"));
274     if (!dateRunNumber.IsNull()) fRunno = dateRunNumber.Atoll();
275   }
276
277   // Open our output files
278   if (!OpenFiles(appendRun)) return false;
279  
280   // Create our reader 
281   AliFMDRawReader* fmdReader  = new AliFMDRawReader(reader,0);
282   TClonesArray*    digitArray = new TClonesArray("AliFMDDigit",0);
283   
284   // Look for SOD  
285   Bool_t sodread = kFALSE;
286   for(Int_t i=0;i<3;i++) {
287     if (!reader->NextEvent()) { 
288       // Read Start-of-Run / Start-of-Files event
289       AliWarning(Form("Failed to read the %d%s event",
290                       i+1, (i == 0 ? "st" : (i == 1 ? "nd" : "rd"))));
291       break;
292     }
293     
294     UInt_t eventType = reader->GetType();
295     if(eventType == AliRawEventHeaderBase::kStartOfData || 
296        eventType == AliRawEventHeaderBase::kFormatError) { 
297       
298       WriteConditionsData(fmdReader);
299       sodread = kTRUE;
300       break;
301     }
302   }
303   if (isBase) 
304     // IF this is the base DA, do thing more 
305     return true;
306
307   // Initialize everything 
308   Init();
309   // Initialize our containers 
310   InitContainer(diagFile);
311   if (AliLog::GetDebugLevel("FMD","") >= 3) { 
312     fDetectorArray.ls();
313   }
314   
315   if(!sodread) 
316     AliWarning("No SOD event detected!");
317   
318   int lastProgress = 0;
319
320   // Reset event counters 
321   for(Int_t i = 0; i< 3;i++) fNEventsPerDetector[i] = 0;
322
323   // Loop until we have enough 
324   for(Int_t n = 1; !HaveEnough(n); n++) {
325     AliDebugF(1,"Get the next event %d", n);
326     if(!reader->NextEvent()) { n--; continue; }
327     UInt_t eventType = reader->GetType();
328     AliDebugF(3, "Event type is %d", eventType);
329     if(eventType != AliRawEventHeaderBase::kPhysicsEvent) { n--; continue; }
330
331     Bool_t seen[] = { false, false, false };
332     SetCurrentEvent(n);
333     digitArray->Clear();
334     fmdReader->ReadAdcs(digitArray);
335   
336     for(Int_t i = 0; i<digitArray->GetEntriesFast();i++) {
337       AliFMDDigit* digit = static_cast<AliFMDDigit*>(digitArray->At(i));
338       UShort_t det = digit->Detector();
339       fSeenDetectors[det-1] = true;
340       seen[det-1]           = true;
341
342       // Only fill if we do not have enough for this detector 
343       if (Int_t(fNEventsPerDetector[det-1]) < GetRequiredEvents()) 
344         FillChannels(digit);
345     }
346     for(Int_t i = 0; i< 3;i++) 
347       if (seen[i]) (fNEventsPerDetector[i])++;
348       
349     FinishEvent();
350     
351     Int_t nReq = GetRequiredEvents();
352     AliDebugF(5, "%9d: %6d/%6d %6d/%6d %6d/%6d", n, 
353               fNEventsPerDetector[0], nReq,
354               fNEventsPerDetector[1], nReq,
355               fNEventsPerDetector[2], nReq);
356
357     int progress = GetProgress(n);
358     if (progress <= lastProgress) continue;
359     lastProgress = progress;
360     std::cout << "Progress: " << lastProgress << " / 100 " << std::endl;
361   }
362
363   // Now at end of processing 
364   AliInfoF("Looped over %d events (%d,%d,%d)",GetCurrentEvent(),
365            fNEventsPerDetector[0], 
366            fNEventsPerDetector[1], 
367            fNEventsPerDetector[2]);
368   WriteHeaderToFile();
369   
370   // Analyse the data and make summaries 
371   for(UShort_t det=1;det<=3;det++) {
372     if (!fSeenDetectors[det-1]) continue;
373     std::cout << "FMD" << det << std::endl;
374     UShort_t firstRing = (det == 1 ? 1 : 0);
375     for (UShort_t ir = firstRing; ir < 2; ir++) {
376       Char_t   ring = (ir == 0 ? 'O' : 'I');
377       UShort_t nsec = (ir == 0 ? 40  : 20);
378       UShort_t nstr = (ir == 0 ? 256 : 512);
379
380       if (fMakeSummaries) MakeSummary(det, ring);
381
382       std::cout << " Ring " << ring << ": " << std::flush;
383       for(UShort_t sec =0; sec < nsec;  sec++)  {
384         for(UShort_t strip = 0; strip < nstr; strip++) {
385           Analyse(det,ring,sec,strip);
386         }
387         std::cout << '.' << std::flush;
388       }
389       // if(fSaveHistograms)
390       // diagFile->Flush();
391       std::cout << "done" << std::endl;
392     }
393   }
394   
395   // If we have an output file, close it  
396   if(fOutputFile.is_open()) {
397     fOutputFile.write("# EOF\n",6);
398     fOutputFile.close();
399   }
400   
401   // Do final stuff on the diagnostics file 
402   Terminate(diagFile);
403     
404   // Possibly write histograms to diagnostics file 
405   if(fSaveHistograms ) {
406     diagFile = TFile::Open(fDiagnosticsFilename.Data(),"RECREATE");
407     fDetectorArray.Write("FMD", TObject::kSingleKey);
408     fSummaries.Write();
409     AliInfo("Closing diagnostics file - please wait ...");
410     // diagFile->Write();
411     diagFile->Close();
412     AliInfo("done");
413   }
414
415   // All is good 
416   return true;
417 }
418 //_____________________________________________________________________
419
420 void AliFMDBaseDA::InitContainer(TDirectory* diagFile)
421 {
422   //Prepare container for diagnostics    
423   TDirectory* savDir   = gDirectory;
424   Bool_t owners = true;
425
426   for(UShort_t det=1;det<=3;det++) {
427     Array* detArray = new Array(det == 1 ? 1 : 2);
428     detArray->SetOwner(owners);
429     detArray->SetName(GetDetectorPath(det, false));
430     fDetectorArray.AddAtAndExpand(detArray,det);
431
432     TDirectory* detDir = 0;
433     if (diagFile) {
434       diagFile->cd();
435       detDir = diagFile->mkdir(GetDetectorPath(det, kFALSE));
436     }
437
438     UShort_t FirstRing = (det == 1 ? 1 : 0);
439     for (UShort_t ir = FirstRing; ir < 2; ir++) {
440       Char_t   ring = (ir == 0 ? 'O' : 'I');
441       UShort_t nsec = (ir == 0 ? 40  : 20);
442       UShort_t nstr = (ir == 0 ? 256 : 512);
443       Array* ringArray = new Array(nsec);
444       ringArray->SetOwner(owners);
445       ringArray->SetName(GetRingPath(det, ring, false));
446       ringArray->ResetBit(TObject::kMustCleanup);
447       detArray->AddAtAndExpand(ringArray,ir);
448
449
450       TDirectory* ringDir = 0;
451       if (detDir) { 
452         detDir->cd();
453         ringDir = detDir->mkdir(GetRingPath(det,ring, kFALSE));
454       }
455       
456
457       for(UShort_t sec =0; sec < nsec;  sec++)  {
458         Array* sectorArray = new Array(nstr);
459         sectorArray->SetOwner(owners);
460         sectorArray->SetName(GetSectorPath(det, ring, sec, false));
461         sectorArray->ResetBit(TObject::kMustCleanup);
462         ringArray->AddAtAndExpand(sectorArray,sec);
463
464
465         TDirectory* secDir = 0;
466         if (ringDir) { 
467           ringDir->cd();
468           secDir = ringDir->mkdir(GetSectorPath(det, ring, sec, kFALSE));
469         }
470         
471         for(UShort_t strip = 0; strip < nstr; strip++) {
472           if (secDir) { 
473             secDir->cd();
474             secDir->mkdir(GetStripPath(det, ring, sec, strip, kFALSE));
475           }
476           Array* stripArray = new Array(0);
477           stripArray->SetOwner(owners);
478           stripArray->SetName(GetStripPath(det, ring, sec, strip, false));
479           stripArray->ResetBit(TObject::kMustCleanup);
480           sectorArray->AddAtAndExpand(stripArray, strip);
481           AddChannelContainer(stripArray, det, ring, sec, strip);
482         }
483         AddSectorSummary(sectorArray, det, ring, sec, nstr);
484       }
485     }
486   }
487   savDir->cd();
488 }
489
490 //_____________________________________________________________________ 
491 void AliFMDBaseDA::WriteConditionsData(AliFMDRawReader* fmdReader) 
492 {
493   //Write the conditions data to file
494   AliFMDParameters* pars       = AliFMDParameters::Instance();
495   fConditionsFile.write(Form("# %s \n",pars->GetConditionsShuttleID()),14);
496   TDatime now;
497   fConditionsFile << "# This file created from run number " << fRunno 
498                   << " at " << now.AsString() << std::endl;
499   
500   AliFMDCalibSampleRate* sampleRate = new AliFMDCalibSampleRate();
501   AliFMDCalibStripRange* stripRange = new AliFMDCalibStripRange();
502   
503   fmdReader->ReadSODevent(sampleRate,stripRange,fPulseSize,fPulseLength,
504                           fSeenDetectors);
505
506   sampleRate->WriteToFile(fConditionsFile, fSeenDetectors);
507   stripRange->WriteToFile(fConditionsFile, fSeenDetectors);
508
509  
510   // Zero Suppresion
511   
512   // Strip Range
513   
514   fConditionsFile.write("# Gain Events \n",15);
515   
516   for(UShort_t det=1; det<=3;det++) {
517     if (!fSeenDetectors[det-1]) { 
518       continue;
519     }
520     UShort_t firstring = (det == 1 ? 1 : 0);
521     for(UShort_t iring = firstring; iring <=1;iring++) {
522       Char_t ring = (iring == 1 ? 'I' : 'O');
523       for(UShort_t board =0 ; board <=1; board++) {
524         
525         Int_t idx = GetHalfringIndex(det,ring,board);
526         
527         fConditionsFile << det                     << ','
528                         << ring                    << ','
529                         << board                   << ','
530                         << fPulseLength.At(idx)    << "\n";
531         
532       }
533     }
534   }
535   
536   fConditionsFile.write("# Gain Pulse \n",14);
537   
538   for(UShort_t det=1; det<=3;det++) {
539     if (!fSeenDetectors[det-1]) { 
540       continue;
541     }
542     UShort_t firstring = (det == 1 ? 1 : 0);
543     for(UShort_t iring = firstring; iring <=1;iring++) {
544       Char_t ring = (iring == 1 ? 'I' : 'O');
545       for(UShort_t board =0 ; board <=1; board++) {
546         
547         Int_t idx = GetHalfringIndex(det,ring,board);
548         
549         fConditionsFile << det                     << ','
550                         << ring                    << ','
551                         << board                   << ','
552                         << fPulseSize.At(idx)      << "\n";
553         
554       }
555     }
556   }
557   // sampleRate->WriteToFile(std::cout, fSeenDetectors);
558   // stripRange->WriteToFile(std::cout, fSeenDetectors);
559
560   if(fConditionsFile.is_open()) {
561     
562     fConditionsFile.write("# EOF\n",6);
563     fConditionsFile.close();
564     
565   }
566   
567 }
568 //_____________________________________________________________________ 
569 Int_t AliFMDBaseDA::GetHalfringIndex(UShort_t det, Char_t ring, 
570                                      UShort_t board) const 
571 {
572   // Get the index corresponding to a half-ring 
573   // 
574   // Parameters: 
575   //   det    Detector number 
576   //   ring   Ring identifier 
577   //   board  Board number 
578   //
579   // Return 
580   //   Internal index of the board 
581   UShort_t iring  =  (ring == 'I' ? 1 : 0);
582   
583   Int_t index = (((det-1) << 2) | (iring << 1) | (board << 0));
584   
585   return index-2;
586   
587 }
588 //_____________________________________________________________________ 
589 void AliFMDBaseDA::Rotate(const char* base, int max) const
590 {
591   // 
592   // Rotate a set of files.   base is the basic name of the files.
593   // If the file base.max exists it is removed. 
594   // If the file base.n exists (where n < max) it is renamed to
595   // base.(n-1).  
596   // If the file base exists, it is renamed to base.1 
597   //
598   // Parameters:
599   //   base Base name of the files
600   //   max  Maximum number to keep (minus one for the current).
601
602   // Note:  TSystem::AccessPathName returns false if the condition is
603   // fulfilled! 
604
605   // Check if we have base.max, and if so, remove it. 
606   TString testName(Form("%s.%d", base, max));
607   if (!gSystem->AccessPathName(testName.Data())) 
608     gSystem->Unlink(testName.Data());
609     
610   // Loop down from max-1 to 1 and move files 
611   for (int i = max-1; i >= 1; i--) { 
612     testName = Form("%s.%d", base, i);
613     if (!gSystem->AccessPathName(testName.Data())) {
614       TString newName(Form("%s.%d", base, i+1));
615       gSystem->Rename(testName.Data(), newName.Data());
616     }
617   }
618
619   // If we have the file base, rename it to base.1 
620   testName = Form("%s", base);
621   if (!gSystem->AccessPathName(testName.Data())){
622     TString newName(Form("%s.%d", base, 1));
623     gSystem->Rename(testName.Data(), newName.Data());
624   }
625 }
626
627 //_____________________________________________________________________ 
628 TH2*
629 AliFMDBaseDA::MakeSummaryHistogram(const char* prefix, const char* title, 
630                                    UShort_t d, Char_t r) 
631 {
632   // 
633   // Utility function for defining summary histograms 
634   // 
635   // Parameters:
636   //    det    Detector 
637   //    ring   Ring identifier 
638   //    prefix Histogram prefix 
639   //    title  Histogram title 
640   //
641   Int_t nX = ((d == 1 || r == 'I' || r == 'i') ?  20 :  40);
642   Int_t nY = ((d == 1 || r == 'I' || r == 'i') ? 512 : 256);
643   
644   TString n = TString::Format("%sFMD%d%c", prefix, d, r);
645   TString t = TString::Format("%s for FMD%d%c", title, d, r);
646   TH2* ret = new TH2F(n, t, nX, -0.5, nX-0.5, nY, -0.5, nY-0.5);
647   ret->SetXTitle("Sector #");
648   ret->SetYTitle("Strip #");
649   ret->SetDirectory(0);
650   // if (!fSummaries) fSummaries = new Array;
651   fSummaries.Add(ret);
652   return ret;
653 }
654
655 //_____________________________________________________________________ 
656 AliFMDBaseDA::Array*
657 AliFMDBaseDA::GetDetectorArray(UShort_t det)
658 {
659   if (det < 1 || det > 3) { 
660     AliErrorF("Detector index %d out of bounds", det);
661     return 0;
662   }
663   return static_cast<Array*>(fDetectorArray.At(det));
664 }
665 //_____________________________________________________________________ 
666 AliFMDBaseDA::Array*
667 AliFMDBaseDA::GetRingArray(UShort_t det, Char_t ring)
668 {
669   Int_t idx = (ring == 'O' || ring == 'o' ? 0 : 1);
670   Array* array = GetDetectorArray(det);
671   if (!array) return 0;
672   array = static_cast<Array*>(array->At(idx));
673   if (!array) AliErrorF("No ring array for FMD%d%c (%d)", det, ring, idx);
674   return array;
675 }
676 //_____________________________________________________________________ 
677 AliFMDBaseDA::Array*
678 AliFMDBaseDA::GetSectorArray(UShort_t det, Char_t ring, UShort_t sector)
679 {
680   Array* array = GetRingArray(det,ring);
681   if (!array) return 0;
682   array = static_cast<Array*>(array->At(sector));
683   if (!array) AliErrorF("No sector array for FMD%d%c[%02d]", det, ring, sector);
684   return array;
685 }
686 //_____________________________________________________________________ 
687 AliFMDBaseDA::Array*
688 AliFMDBaseDA::GetStripArray(UShort_t det, Char_t ring, 
689                             UShort_t sector, UShort_t strip)
690 {
691   Array* array = GetSectorArray(det,ring,sector);
692   if (!array) return 0;
693   array = static_cast<Array*>(array->At(strip));
694   if (!array) AliErrorF("No strip array for FMD%d%c[%02d,%03d]", 
695                         det, ring, sector, strip);
696   return array;
697 }
698
699 //=====================================================================
700 AliFMDBaseDA::Runner::Runner()
701   : fReader(0),
702     fSource(""),
703     fDiagFile(""), 
704     fDiag(false),
705     fAll(false),
706     fFast(true),
707     fUpload(true),
708     fAppendRun(false),
709     fOwnUpload(false)
710 {}
711
712 //_____________________________________________________________________ 
713 void
714 AliFMDBaseDA::Runner::AddHandlers()
715 {
716   gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
717                                         "*",
718                                         "TStreamerInfo",
719                                         "RIO",
720                                         "TStreamerInfo()");
721   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", "Minuit", 
722                                         "TMinuitMinimizer",
723                                         "Minuit", 
724                                         "TMinuitMinimizer(const char *)");
725   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", 
726                                         "GSLMultiMin", 
727                                         "ROOT::Math::GSLMinimizer",
728                                         "MathMore", 
729                                         "GSLMinimizer(const char *)");
730   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", 
731                                         "GSLMultiFit", 
732                                         "ROOT::Math::GSLNLSMinimizer",
733                                         "MathMore", "GSLNLSMinimizer(int)");
734   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", 
735                                         "GSLSimAn", 
736                                         "ROOT::Math::GSLSimAnMinimizer",
737                                         "MathMore", 
738                                         "GSLSimAnMinimizer(int)");
739   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", 
740                                         "Linear", 
741                                         "TLinearMinimizer",
742                                         "Minuit", 
743                                         "TLinearMinimizer(const char *)");
744   gROOT->GetPluginManager()->AddHandler("ROOT::Math::Minimizer", 
745                                         "Fumili", 
746                                         "TFumiliMinimizer",
747                                         "Fumili", 
748                                         "TFumiliMinimizer(int)");
749 }
750 //_____________________________________________________________________ 
751 void
752 AliFMDBaseDA::Runner::ShowUsage(std::ostream& o, const char* progname)
753 {
754   o << "Usage: " << progname << " SOURCE [OPTIONS]\n\n"
755     << "Options:\n"
756     << "   -h,--help                Show this help\n"
757     << "   -d,--diagnostics[=FILE]  Write diagnostics to file\n"
758     << "   -D,--debug=LEVEL         Set the debug level\n"
759     << "   -A,--all                 Try to get data from all detectors (" << fAll << ")\n"
760     << "   -F,--fast                Fast exit (" <<fFast <<")\n"
761     << "   -U,--upload              Upload to FXS (" <<fUpload <<")\n"
762     << "   -R,--append              Append run # to filename (" <<fAppendRun <<")\n"
763     << "   -Z,--own                 Use custom upload (" << fOwnUpload << ")\n"
764     << "\n"
765     << "SOURCE is one of\n"
766     << " * FILE.raw                Raw data file\n"
767     << " * FILE.root               ROOT'ified raw data file\n"
768     << " * collection://FILE.root  Event list in a ROOT file\n"
769     << " * collection://FILE       File containing list of ROOT files\n"
770     << " * ^FMD                    Monitor source\n"
771     << "There are other options too.  Check the AliRawReader docs\n"
772     << std::endl;
773 }
774  
775 //_____________________________________________________________________ 
776 namespace {
777   Bool_t ExtractValue(const TString& arg, TString& val)
778   {
779     val = "";
780     Int_t eq = arg.Index("=");
781     if (eq == kNPOS) return false;
782     
783     val = arg(eq+1, arg.Length()-eq-1);
784     return true;
785   }
786 }
787       
788 //_____________________________________________________________________ 
789 Int_t
790 AliFMDBaseDA::Runner::Init(int argc, char** argv, Bool_t reader)
791 {
792   AddHandlers();
793
794   // --- Process the command line ------------------------------------
795   TString source;
796   Int_t   debugLevel  = 0;
797   Bool_t  help        = false;
798
799   for (int i = 1; i < argc; i++) { 
800     TString arg(argv[i]);
801     Bool_t  badOption   = false;
802     
803     if (arg[0] == '-') { // It's an option 
804       if (arg[1] == '-') { // It's a long option 
805         TString val;
806         if      (arg.EqualTo("--help"))     help = true; 
807         else if (arg.BeginsWith("--debug")) {
808           if (ExtractValue(arg, val))
809             debugLevel = val.Atoi();
810         }
811         else if (arg.BeginsWith("--diagnostics")) { 
812           fDiag = !fDiag;
813           if (ExtractValue(arg, val)) 
814             fDiagFile = val;
815         }
816         else if (arg.EqualTo("--all"))    fAll  = !fAll;
817         else if (arg.EqualTo("--fast"))   fFast = !fFast;
818         else if (arg.EqualTo("--upload")) fUpload = !fUpload;
819         else if (arg.EqualTo("--append")) fAppendRun = !fAppendRun;
820         else if (arg.EqualTo("--own"))    fOwnUpload = !fOwnUpload;
821         else badOption = true;
822       }
823       else { // Short option 
824         TString next(i < argc-1 ? argv[i+1] : "");
825         switch (arg[1]) { 
826         case 'h': help = true; break;
827         case 'd': fDiag = !fDiag; 
828           if (!next.IsNull() && next[0] != '-') {
829             fDiagFile = next;
830             i++;
831           }
832           break;
833         case 'D': 
834           if (!next.IsNull() && next[0] != '-') {
835             debugLevel = next.Atoi();
836             i++;
837           }
838           break;
839         case 'A': fAll       = !fAll ; break ;
840         case 'F': fFast      = !fFast ; break ;
841         case 'U': fUpload    = !fUpload ; break ;
842         case 'R': fAppendRun = !fAppendRun ; break ;
843         case 'Z': fOwnUpload = !fOwnUpload ; break ;
844         default: badOption = true;
845         }
846       } // End of options
847       if (badOption) { 
848         std::cerr << argv[0] << ": Unknown option " << argv[i] 
849                   << std::endl;
850         return -1;
851       }
852     }
853     else { // source or compatibility debug level 
854       source = arg;
855     }
856   }
857   
858   // --- Check if help was requested ---------------------------------
859   if (help) { 
860     ShowUsage(std::cout, argv[0]);
861     return 1;
862   }
863
864   // --- Check if we have a source -----------------------------------
865   if (source.IsNull()) { 
866     std::cerr << "No source given" << std::endl;
867     return -2;
868   }
869   fSource = source;
870
871   // --- Initialize various things -----------------------------------
872   AliFMDParameters::Instance()->Init(kFALSE,0);
873
874   //This will only work for FDR 1 data. When newer data becomes
875   //available the ! must be removed!
876   Bool_t old = kTRUE;
877   AliFMDParameters::Instance()->UseCompleteHeader(old);
878   
879   AliLog::EnableDebug(debugLevel > 0);
880   AliLog::SetModuleDebugLevel("FMD", debugLevel);
881
882   // --- Make our reader ---------------------------------------------
883   if (!reader) return 0;
884   fReader = AliRawReader::Create(source);
885   if (!fReader) { 
886     std::cerr << "Failed to make raw reader for source \"" << source 
887               << "\"" << std::endl;
888     return -3;
889   }
890   return 0;
891 }
892
893 //_____________________________________________________________________ 
894 Int_t
895 AliFMDBaseDA::Runner::RunNumber() const
896
897   Int_t run = -1;
898   if (fReader) 
899     run = fReader->GetRunNumber(); 
900   if (run <= 0) {
901     TString dateRunNumber(gSystem->Getenv("DATE_RUN_NUMBER"));
902     if (!dateRunNumber.IsNull()) 
903       run = dateRunNumber.Atoll();
904     else 
905       run = -1;
906   }
907   return run;
908 }
909
910 //_____________________________________________________________________ 
911 Bool_t
912 AliFMDBaseDA::Runner::Exec(AliFMDBaseDA& da)
913 {
914   TStopwatch timer;
915   timer.Start();
916
917   da.SetSaveDiagnostics(fDiag || !fDiagFile.IsNull());
918   da.SetTryAll(fAll);
919   if (!fDiagFile.IsNull()) da.SetDiagnosticsFilename(fDiagFile);
920
921   Bool_t ret = da.Run(fReader, fAppendRun);
922
923   timer.Stop();
924   timer.Print();
925
926   return ret;
927 }
928
929 #if 0
930 AliFMDBaseDA::_Array::~_Array()
931 {
932   // Printf("Deleting %s (%p)", this->GetName(), this);
933   // SetOwner(false);
934   // Clear();
935 }
936 #endif
937   
938   
939 //_____________________________________________________________________ 
940
941 //_____________________________________________________________________ 
942 //
943 // EOF
944 //
945