AliRoot Core  a565103 (a565103)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MUONTriggerChamberEfficiency.C
Go to the documentation of this file.
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 #if !defined(__CINT__) || defined(__MAKECINT__)
19 
20 #include "Riostream.h"
21 
22 // ROOT includes
23 #include "TGrid.h"
24 #include "TString.h"
25 #include "TFile.h"
26 #include "TH1.h"
27 
28 // MUON includes
29 #include "AliMUONCDB.h"
30 #include "AliMUONCalibrationData.h"
33 #include "AliCDBManager.h"
34 #include "AliCDBRunRange.h"
35 
36 #endif
37 
46 
47 void MUONTriggerChamberEfficiency(TString inputFile = "./MUON.TriggerEfficiencyMap.root",
48  TString outputCDB = "",
49  Int_t firstRun=0, Int_t lastRun = AliCDBRunRange::Infinity()
50 )
51 {
61 
62  AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
63 
64  AliMUONTriggerEfficiencyCells* effMap = new AliMUONTriggerEfficiencyCells(inputFile.Data());
65 
66  if ( outputCDB.IsNull() ){
67  // Draw the efficiency and exit
68  AliCDBManager::Instance()->SetRun(firstRun);
70 
71  trigChEff->DisplayEfficiency(kFALSE,kFALSE);
72  return;
73  }
74 
75 
76  // Write efficiency on OCDB
77 
78  AliCDBManager::Instance()->SetSpecificStorage("MUON/Calib/TriggerEfficiency", outputCDB.Data());
79 
80  AliMUONCDB::WriteToCDB(effMap, "MUON/Calib/TriggerEfficiency", firstRun, lastRun, "Measured efficiencies");
81 }
82 
83 //____________________________________________________________
84 void ShowOCDBmap(Int_t runNumber = 0, TString specificCDB="", TString ocdbPath = "local://$ALICE_ROOT/OCDB", TString runType="Full")
85 {
92  if ( ocdbPath.BeginsWith("alien://") || ocdbPath.BeginsWith("raw://"))
93  TGrid::Connect("alien://");
94 
95  if (!ocdbPath.CompareTo("MC"))
96  AliCDBManager::Instance()->SetDefaultStorage(ocdbPath.Data(),runType.Data());
97  else
98  AliCDBManager::Instance()->SetDefaultStorage(ocdbPath.Data());
99 
100  if ( !specificCDB.IsNull() )
101  AliCDBManager::Instance()->SetSpecificStorage("MUON/Calib/TriggerEfficiency", specificCDB.Data());
102  AliCDBManager::Instance()->SetRun(runNumber);
103  AliMUONCalibrationData calib(runNumber);
104 
106  trigChEff->DisplayEfficiency(kFALSE,kFALSE);
107 }
108 
109 
110 //____________________________________________________________
111 void FillHisto(TH1* histo, Int_t nevents)
112 {
114  for ( Int_t ibin=1; ibin<=histo->GetXaxis()->GetNbins(); ++ibin ) {
115  Double_t binCenter = histo->GetXaxis()->GetBinCenter(ibin);
116  for ( Int_t ievent=0; ievent<nevents; ++ievent ) {
117  histo->Fill(binCenter);
118  }
119  }
120 }
121 
122 //____________________________________________________________
123 void BuildDefaultMap(TString outFilename="/tmp/defTrigChEff.root", Double_t globalValue = 1., Int_t nevents = 100000)
124 {
126 
127  // Create histograms
128  enum { kBendingEff, kNonBendingEff, kBothPlanesEff, kAllTracks, kNcounts};
129  TString countTypeName[kNcounts] = {"bendPlane", "nonBendPlane","bothPlanes", "allTracks"};
130 
131  const Char_t* yAxisTitle = "counts";
132 
133  const Int_t kNboards = 234; //AliMpConstants::NofLocalBoards();
134  const Int_t kFirstTrigCh = 11;//AliMpConstants::NofTrackingChambers()+1;
135  const Int_t kNchambers = 4;
136  const Int_t kNslats = 18;
137 
138  Int_t chamberBins = kNchambers;
139  Float_t chamberLow = kFirstTrigCh-0.5, chamberHigh = kFirstTrigCh+kNchambers-0.5;
140  const Char_t* chamberName = "chamber";
141 
142  Int_t slatBins = kNslats;
143  Float_t slatLow = 0-0.5, slatHigh = kNslats-0.5;
144  const Char_t* slatName = "slat";
145 
146  Int_t boardBins = kNboards;
147  Float_t boardLow = 1-0.5, boardHigh = kNboards+1.-0.5;
148  const Char_t* boardName = "board";
149 
150  TString baseName, histoName, histoTitle;
151  TList* histoList = new TList();
152  histoList->SetOwner();
153 
154  TH1F* histo;
155 
156  for(Int_t icount=0; icount<kNcounts; icount++){
157  histoName = Form("%sCountChamber", countTypeName[icount].Data());
158  histo = new TH1F(histoName, histoName,
159  chamberBins, chamberLow, chamberHigh);
160  histo->GetXaxis()->SetTitle(chamberName);
161  histo->GetYaxis()->SetTitle(yAxisTitle);
162  Double_t nfills = ( icount == kAllTracks ) ? nevents : globalValue * (Double_t)nevents;
163  FillHisto(histo, (Int_t)nfills);
164  histoList->AddLast(histo);
165  } // loop on counts
166 
167  for(Int_t icount=0; icount<kNcounts; icount++){
168  for(Int_t ch=0; ch<kNchambers; ch++){
169  histoName = Form("%sCountSlatCh%i", countTypeName[icount].Data(), kFirstTrigCh+ch);
170  histo = new TH1F(histoName, histoName,
171  slatBins, slatLow, slatHigh);
172  histo->GetXaxis()->SetTitle(slatName);
173  histo->GetYaxis()->SetTitle(yAxisTitle);
174  Double_t nfills = ( icount == kAllTracks ) ? nevents : globalValue * (Double_t)nevents;
175  FillHisto(histo, (Int_t)nfills);
176  histoList->AddLast(histo);
177  } // loop on chamber
178  } // loop on counts
179 
180  for(Int_t icount=0; icount<kNcounts; icount++){
181  for(Int_t ch=0; ch<kNchambers; ch++){
182  histoName = Form("%sCountBoardCh%i", countTypeName[icount].Data(), kFirstTrigCh+ch);
183  histo = new TH1F(histoName, histoName,
184  boardBins, boardLow, boardHigh);
185  histo->GetXaxis()->SetTitle(boardName);
186  histo->GetYaxis()->SetTitle(yAxisTitle);
187  Double_t nfills = ( icount == kAllTracks ) ? nevents : globalValue * (Double_t)nevents;
188  FillHisto(histo, (Int_t)nfills);
189  histoList->AddLast(histo);
190  } // loop on chamber
191  } // loop on counts
192 
193  TFile* outFile = TFile::Open(outFilename,"create");
194  histoList->Write("triggerChamberEff",TObject::kSingleKey);
195  outFile->Close();
196 }
197 
198 
199 //____________________________________________________________
200 void CompleteEfficiency(TString effFileWithHoles, TString effFileCompatible, TString outFilename)
201 {
202  //
205  //
206 
207  TList* histoList[2] = {0x0, 0x0};
208  TString filenames[2] = {effFileWithHoles, effFileCompatible};
209  for ( Int_t ifile=0; ifile<2; ifile++ ) {
210  TFile* file = TFile::Open(filenames[ifile].Data());
211  if ( ! file ) {
212  printf("Fatal: cannot find %s\n", filenames[ifile].Data());
213  return;
214  }
215  histoList[ifile] = static_cast<TList*> (file->FindObjectAny("triggerChamberEff"));
216  if ( ! histoList[ifile] ) {
217  printf("Cannot find histo list in file %s\n", filenames[ifile].Data());
218  return;
219  }
220  }
221 
222  TString detElemName[2] = {"Slat", "Board"};
223  enum { kBendingEff, kNonBendingEff, kBothPlanesEff, kAllTracks, kNcounts};
224  TString countTypeName[kNcounts] = {"bendPlane", "nonBendPlane","bothPlanes", "allTracks"};
225 
226  Bool_t isChanged = kFALSE;
227  TString histoName = "";
228  for ( Int_t idet=0; idet<2; idet++ ) {
229  for ( Int_t ich=11; ich<=14; ich++ ) {
230  histoName = Form("%sCount%sCh%i", countTypeName[kAllTracks].Data(), detElemName[idet].Data(), ich);
231  TH1* allTracksHisto = static_cast<TH1*> (histoList[0]->FindObject(histoName.Data()));
232  for ( Int_t ibin=1; ibin<=allTracksHisto->GetXaxis()->GetNbins(); ibin++ ) {
233  if ( allTracksHisto->GetBinContent(ibin) > 0. ) continue;
234  isChanged = kTRUE;
235  printf("Modifying info for Ch %i %s %3i\n", ich, detElemName[idet].Data(), (Int_t)allTracksHisto->GetXaxis()->GetBinCenter(ibin));
236  // If allTracks has no entries, it means that efficiency could not be calculated for this bin:
237  // fill information from the compatible histogram
238 
239  // Check the statistics collected by the switched off detection element
240  Double_t nTracks = 0;
241  for ( Int_t jch=11; jch<=14; jch++ ) {
242  histoName = Form("%sCount%sCh%i", countTypeName[kAllTracks].Data(), detElemName[idet].Data(), jch);
243  TH1* allTracksOtherCh = static_cast<TH1*> (histoList[0]->FindObject(histoName.Data()));
244  nTracks = allTracksOtherCh->GetBinContent(ibin);
245  if ( nTracks > 0. ) {
246  //printf("Statistics for %s : %g\n", histoName.Data(), nTracks); // REMEMBER TO CUT
247  break;
248  }
249  }
250 
251  histoName = Form("%sCount%sCh%i", countTypeName[kAllTracks].Data(), detElemName[idet].Data(), ich);
252  TH1* allTracksHistoAux = static_cast<TH1*> (histoList[1]->FindObject(histoName.Data()));
253  Double_t nTracksNew = allTracksHistoAux->GetBinContent(ibin);
254  if ( nTracksNew == 0.) {
255  printf("Warning: new histogram has no entries for Ch %i %s %3i\n", ich, detElemName[idet].Data(), (Int_t)allTracksHisto->GetXaxis()->GetBinCenter(ibin));
256  continue;
257  }
258  Double_t scaleFactor = TMath::Min(nTracksNew, nTracks) / nTracksNew;
259  //printf("Statistics ineff %g new %g scaleFactor %g\n", nTracks, nTracksNew, scaleFactor); // REMEMBER TO CUT
260 
261  for ( Int_t icount=0; icount<kNcounts; icount++ ) {
262  histoName = Form("%sCount%sCh%i", countTypeName[icount].Data(), detElemName[idet].Data(), ich);
263  TH1* auxHisto = static_cast<TH1*> (histoList[1]->FindObject(histoName.Data()));
264  TH1* histo = static_cast<TH1*> (histoList[0]->FindObject(histoName.Data()));
265  histo->SetBinContent(ibin, auxHisto->GetBinContent(ibin) * scaleFactor);
266  if ( histo->GetSumw2N() > 0 ) histo->SetBinError(ibin, auxHisto->GetBinError(ibin)*scaleFactor);
267  } // loop on cont types
268  } // loop on histogram bins
269  } // loop on chamber
270  } // loop on detection element (slat or board)
271 
272  if ( ! isChanged ) {
273  printf("Input histograms not modified\n");
274  return;
275  }
276  TFile* outFile = TFile::Open(outFilename,"create");
277  histoList[0]->Write("triggerChamberEff",TObject::kSingleKey);
278  outFile->Close();
279 }
280 
281 
void ShowOCDBmap(Int_t runNumber=0, TString specificCDB="", TString ocdbPath="local://$ALICE_ROOT/OCDB", TString runType="Full")
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
TFile * Open(const char *filename, Long64_t &nevents)
void FillHisto(TH1 *histo, Int_t nevents)
Store and give access to the trigger chamber efficiency.
void MUONTriggerChamberEfficiency(TString inputFile="./MUON.TriggerEfficiencyMap.root", TString outputCDB="", Int_t firstRun=0, Int_t lastRun=AliCDBRunRange::Infinity())
void CompleteEfficiency(TString effFileWithHoles, TString effFileCompatible, TString outFilename)
AliMUONTriggerEfficiencyCells * TriggerEfficiency() const
Get the trigger efficiency map.
void DisplayEfficiency(Bool_t perSlat=kFALSE, Bool_t show2Dhisto=kTRUE)
Single entry point to access MUON calibration data.
void BuildDefaultMap(TString outFilename="/tmp/defTrigChEff.root", Double_t globalValue=1., Int_t nevents=100000)
Calculate, apply and possibly draw trigger chamber efficiency.
void WriteToCDB(const char *calibpath, TObject *object, Int_t startRun, Int_t endRun, Bool_t defaultValues)