AliRoot Core  edcc906 (edcc906)
TestRecPoints.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: TestRecPoints.C 23207 2007-12-20 09:59:20Z ivana $ */
17 
24 
25 #if !defined(__CINT__) || defined(__MAKECINT__)
26 // ROOT
27 #include "TString.h"
28 #include "TH1.h"
29 #include "TList.h"
30 #include "TFile.h"
31 #include "TDirectory.h"
32 
33 #include "TTree.h"
34 #include "TClonesArray.h"
35 
36 // STEER
37 #include "AliRun.h"
38 #include "AliRunLoader.h"
39 #include "AliLoader.h"
40 #include "AliCDBManager.h"
41 
42 // tracker
43 #include "AliGeomManager.h"
44 
45 // MUON
46 #include "AliMpConstants.h"
47 
48 // trigger
49 #include "AliMUONVTriggerStore.h"
50 #include "AliMUONDigitStoreV1.h"
51 #include "AliMUONDigitMaker.h"
52 #include "AliMUONLocalTrigger.h"
53 
54 // tracker
55 #include "AliMUONClusterStoreV2.h"
60 #include "AliMUONVDigitStore.h"
61 #include "AliMUONDigitStoreV2R.h"
62 #include "AliMUONVCluster.h"
63 #include "AliMUONRecoParam.h"
64 #include "AliMUONCDB.h"
65 
66 #endif
67 
70 
72 
74 
75 Int_t GetPlane(Int_t ch, Int_t cath){return kNcathodes * ch + cath;}
76 void ClusterSize(TList&, AliMUONVDigit*, Int_t&, Int_t);
77 
78 // Main Method
79 void TestRecPoints(TString baseDir=".", TString outDir=".", Float_t adcCut = 10., Int_t whatToTest=kTrackTrig, Int_t runNumber=0, TString cdbStorage="local://$ALICE_ROOT/OCDB")
80 {
81  const Int_t kNplanes = kNtrigChambers * kNcathodes;
82  const Int_t kNslats = 18;
83 
84  TList histoListTrig;
85 
86  TH1F* trigStripMult[kNplanes];
87  TH1F* trigClusterSize[kNplanes];
88  TH1F* trigClusterMult[kNplanes];
89  TString histoName, histoTitle;
90 
91  if(whatToTest!=kOnlyTracker){
92  for(Int_t ch=0; ch<kNtrigChambers; ch++){
93  for(Int_t cath=0; cath<kNcathodes; cath++){
94  Int_t iplane = GetPlane(ch, cath);
95  histoName = "stripMult";
96  histoName.Append(Form("Ch%iCath%i",ch,cath));
97  histoTitle = Form("Strip multiplicity Ch %i Cathode %i", ch, cath);
98  trigStripMult[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), kNslats, 0.-0.5, (Float_t)kNslats - 0.5);
99  trigStripMult[iplane]->SetXTitle("slat");
100  histoListTrig.Add(trigStripMult[iplane]);
101 
102  histoName = "clusterSize";
103  histoName.Append(Form("Ch%iCath%i",ch,cath));
104  histoTitle = Form("Cluster size Ch %i Cathode %i", ch, cath);
105  trigClusterSize[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), 10, 0.-0.5, (Float_t)10 - 0.5);
106  trigClusterSize[iplane]->SetXTitle("cluster size");
107  histoListTrig.Add(trigClusterSize[iplane]);
108 
109  histoName = "clusterMult";
110  histoName.Append(Form("Ch%iCath%i",ch,cath));
111  histoTitle = Form("Cluster multiplicity Ch %i Cathode %i", ch, cath);
112  trigClusterMult[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), kNslats, 0.-0.5, (Float_t)kNslats - 0.5);
113  trigClusterMult[iplane]->SetXTitle("slat");
114  histoListTrig.Add(trigClusterMult[iplane]);
115  }
116  }
117  }
118 
119  const Int_t kMaxDetElem = 26;
120  TList histoListTrack;
121  TH1F* trackClusterSize[kNtrackChambers];
122  TH1F* trackClusterMult[kNtrackChambers];
123  TH1F* trackClusterYield[kNtrackChambers];
124 
125  if(whatToTest!=kOnlyTrigger){
126  for(Int_t ch=0; ch<kNtrackChambers; ch++){
127  histoName = "clusterSize";
128  histoName.Append(Form("Ch%i",ch));
129  histoTitle = Form("Cluster size Ch %i", ch);
130  trackClusterSize[ch] = new TH1F(histoName.Data(), histoTitle.Data(), kMaxDetElem, 0.-0.5, (Float_t)kMaxDetElem - 0.5);
131  trackClusterSize[ch]->SetXTitle("cluster size");
132  histoListTrack.Add(trackClusterSize[ch]);
133 
134  histoName = "clusterMult";
135  histoName.Append(Form("Ch%i",ch));
136  histoTitle = Form("Cluster multiplicity vs. Charge Ch %i", ch);
137  trackClusterMult[ch] = new TH1F(histoName.Data(), histoTitle.Data(), 200, 0, 15000.);
138  trackClusterMult[ch]->SetXTitle("Charge (ADC)");
139  histoListTrack.Add(trackClusterMult[ch]);
140 
141  histoName = "clusterYield";
142  histoName.Append(Form("Ch%i",ch));
143  histoTitle = Form("Cluster yield vs. DetElem Ch %i", ch);
144  trackClusterYield[ch] = new TH1F(histoName.Data(), histoTitle.Data(), kMaxDetElem, 0.-0.5, (Float_t)kMaxDetElem - 0.5);
145  trackClusterYield[ch]->SetXTitle("Detector element #");
146  histoListTrack.Add(trackClusterYield[ch]);
147  }
148  }
149 
150 
151  // Creating Run Loader and opening RecPoints
152  TString filename = baseDir.Data();
153  filename.Append("/galice.root");
154  AliRunLoader * RunLoader = AliRunLoader::Open(filename.Data(),"MUONLoader","UPDATE");
155  if (RunLoader ==0x0) {
156  printf(">>> Error : Error Opening %s file \n",filename.Data());
157  return;
158  }
159 
160  // Loading MUON subsystem
161  AliLoader* MUONLoader = RunLoader->GetDetectorLoader("MUON");
162  if(whatToTest!=kOnlyTracker) MUONLoader->LoadRecPoints("READ");
163  if(whatToTest!=kOnlyTrigger) MUONLoader->LoadDigits("READ");
164 
165  TTree *treeR = 0x0, *treeD = 0x0;
166  AliMUONVTriggerStore* trigStore = 0x0;
167  AliMUONLocalTrigger* locTrg = 0x0;
168 
169  AliMUONVDigitStore* digitStoreTrack = 0x0;
170  AliMUONDigitStoreV2R digitStoreTrackCut;
171  AliMUONVCluster* cluster = 0x0;
172 
173  // Load mapping
174  AliCDBManager::Instance()->SetDefaultStorage(cdbStorage.Data());
176  if (!AliMUONCDB::LoadMapping()) return;
177 
178  AliMUONGeometryTransformer* transformer = 0x0;
179 
180  AliMUONRecoParam* recoParam = 0x0;
181 
182  AliMUONClusterStoreV2 clusterStore;
184 
185  AliMUONSimpleClusterServer* clusterServer = 0x0;
186 
187  Int_t firstChamber(0);
188  Int_t lastChamber(9);
189 
190  if(whatToTest!=kOnlyTrigger){
191  // Import TGeo geometry
192  TString geoFilename = baseDir.Data();
193  geoFilename.Append("/geometry.root");
194  if ( ! AliGeomManager::GetGeometry() ) {
195  AliGeomManager::LoadGeometry(geoFilename);
196  if (! AliGeomManager::GetGeometry() ) {
197  printf(">>> Error : getting geometry from file %s failed\n", geoFilename.Data());
198  return;
199  }
200  }
201  transformer = new AliMUONGeometryTransformer();
202  // Load geometry data
203  transformer->LoadGeometryData();
204  // Load reconstruction parameters
205  recoParam = AliMUONCDB::LoadRecoParam();
206  if (!recoParam) return;
207  recoParam->Print("FULL");
208  clusterFinder = new AliMUONClusterFinderMLEM(kFALSE,new AliMUONPreClusterFinder);
209  clusterServer = new AliMUONSimpleClusterServer(clusterFinder,*transformer);
210  }
211 
212  AliMUONDigitStoreV1 digitStore;
213  AliMUONVDigit* mDigit;
214 
215  Int_t clusterSize;
216 
217  AliMUONDigitMaker digitMaker;
218  TList digitsList[kNplanes];
219 
220  Int_t nevents = RunLoader->GetNumberOfEvents();
221  Int_t analysisFrac = nevents/10+1;
222 
223  printf("\nNumber of events = %i\n\n",nevents);
224 
225  for(Int_t ievent=0; ievent<nevents; ievent++){
226  if(ievent%analysisFrac==0) printf("Analysing event = %i\n",ievent);
227  RunLoader->GetEvent(ievent);
228  if(whatToTest!=kOnlyTracker){
229  digitStore.Clear();
230  treeR = MUONLoader->TreeR();
231  trigStore = AliMUONVTriggerStore::Create(*treeR);
232 
233  if ( trigStore == 0x0 ) continue;
234  trigStore->Clear();
235  trigStore->Connect(*treeR);
236  treeR->GetEvent(0);
237 
238  TIter nextLocal(trigStore->CreateLocalIterator());
239  while ( (locTrg = static_cast<AliMUONLocalTrigger*>( nextLocal() )) != NULL )
240  {
241  TArrayS xyPattern[2];
242  locTrg->GetXPattern(xyPattern[0]);
243  locTrg->GetYPattern(xyPattern[1]);
244 
245  Int_t nBoard = locTrg->LoCircuit();
246  digitMaker.TriggerDigits(nBoard, xyPattern, digitStore);
247  }
248 
249  TIter next(digitStore.CreateIterator());
250 
251  while ( ( mDigit = static_cast<AliMUONVDigit*>(next()) ) )
252  {
253  Int_t detElemId = mDigit->DetElemId();
254  Int_t ch = detElemId/100 - 11;
255  Int_t cathode = mDigit->Cathode();
256  Int_t slat = detElemId%100;
257  Int_t iplane = GetPlane(ch, cathode);
258  trigStripMult[iplane]->Fill(slat);
259  digitsList[iplane].Add(mDigit);
260  } // loop on digits
261 
262  for(Int_t iplane=0; iplane<kNplanes; iplane++){
263  while(digitsList[iplane].GetEntries()){
264  clusterSize=1;
265  mDigit = (AliMUONVDigit*)digitsList[iplane].At(0);
266  digitsList[iplane].Remove(mDigit);
267  ClusterSize(digitsList[iplane],mDigit,clusterSize,1);
268  //if(clusterSize>1) printf("Cluster size = %i\n\n",clusterSize);
269  trigClusterSize[iplane]->Fill(clusterSize);
270 
271  Int_t detElemId = mDigit->DetElemId();
272  Int_t slat = detElemId%100;
273  trigClusterMult[iplane]->Fill(slat);
274  } // loop o sorted digits
275  } // loop on planes
276  } // trigger part
277 
278  if(whatToTest!=kOnlyTrigger){
279  clusterStore.Clear();
280  treeD = MUONLoader->TreeD();
281  digitStoreTrack = AliMUONVDigitStore::Create(*treeD);
282 
283  if ( digitStoreTrack == 0x0 ) continue;
284  digitStoreTrack->Clear();
285  digitStoreTrackCut.Clear();
286  digitStoreTrack->Connect(*treeD);
287  treeD->GetEvent(0);
288 
289  // Cut low charge channels (pedestal subtraction)
290  TIter nextDigitTrack(digitStoreTrack->CreateIterator());
291  while ( ( mDigit = static_cast<AliMUONVDigit*>(nextDigitTrack()) ) )
292  {
293  //printf("Digit class = %s",mDigit->ClassName());
294  Float_t charge = mDigit->Charge();
295  if(charge<adcCut) continue;
296  digitStoreTrackCut.Add(*mDigit, AliMUONVDigitStore::kDeny);
297  } // loop on digits
298 
299  TIter nextDigitTrackCut(digitStoreTrackCut.CreateIterator());
300  clusterServer->UseDigits(nextDigitTrackCut,&digitStoreTrackCut);
301 
302  for (Int_t ch = firstChamber; ch <= lastChamber; ++ch ){
303  clusterServer->Clusterize(ch, clusterStore, AliMpArea(),recoParam);
304  }
305 
306  TIter next(clusterStore.CreateIterator());
307  while ( ( cluster = static_cast<AliMUONVCluster*>(next()) ) )
308  {
309  Float_t charge = cluster->GetCharge();
310  if(charge==0.) continue;
311  Int_t ch = cluster->GetChamberId();
312  Int_t npads = cluster->GetNDigits();
313  Int_t detElemId = cluster->GetDetElemId()%100;
314  //printf("charge = %f pads = %i detElemId = %i\n",charge,npads,detElemId);
315 
316  trackClusterSize[ch]->Fill(npads);
317  trackClusterMult[ch]->Fill(charge);
318  trackClusterYield[ch]->Fill(detElemId);
319  } // loop on clusters
320  } // tracker part
321  } // loop on events
322 
323  MUONLoader->UnloadRecPoints();
324  MUONLoader->UnloadDigits();
325  RunLoader->UnloadAll();
326  delete RunLoader;
327  TString outFileName = outDir.Data();
328  outFileName.Append("/outTestRecPoints.root");
329  TFile* outFile = new TFile(outFileName.Data(), "RECREATE");
330  TDirectory* dir = 0x0;
331  if(whatToTest!=kOnlyTracker){
332  outFile->cd();
333  dir = outFile->mkdir("Trigger");
334  dir->cd();
335  histoListTrig.Write();
336  }
337  if(whatToTest!=kOnlyTrigger){
338  outFile->cd();
339  dir = outFile->mkdir("Tracker");
340  dir->cd();
341  histoListTrack.Write();
342  }
343  outFile->Close();
344  printf("\nSee results in %s\n",outFileName.Data());
345 }
346 
347 void ClusterSize(TList& list, AliMUONVDigit* refDigit, Int_t& clusterSize, Int_t depthLevel)
348 {
349  AliMUONVDigit* mDigit = 0x0;
350  for(Int_t idigit=0; idigit<list.GetEntries(); idigit++){
351  mDigit = (AliMUONVDigit*) list.At(idigit);
352  if(mDigit->DetElemId() != refDigit->DetElemId()) continue;
353  Int_t diffX = TMath::Abs(mDigit->PadX() - refDigit->PadX());
354  Int_t diffY = TMath::Abs(mDigit->PadY() - refDigit->PadY());
355  if(diffX>1) continue;
356  if(diffY>1) continue;
357  if(diffX*diffY != 0) continue;
358  clusterSize++;
359  list.Remove(mDigit);
360  //printf("DetElemId = %i Level = %i Ref. (%2i,%2i) Found (%2i,%2i)\n",mDigit->DetElemId(),depthLevel,refDigit->PadX(),refDigit->PadY(),mDigit->PadX(),mDigit->PadY());
361  ClusterSize(list, mDigit, clusterSize, depthLevel+1);
362  Int_t val = idigit + depthLevel - clusterSize;
363  //printf("Level = %i Current digit = %i\t",depthLevel,idigit);
364  idigit = TMath::Max(-1,val);
365  //printf("New digit = %i\n",idigit);
366  }
367 }
void GetYPattern(TArrayS &array) const
return Y pattern array
const Int_t kNcathodes
Definition: TestRecPoints.C:69
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
Reading Raw data class for trigger and tracker chambers.
void UnloadAll(Option_t *detectors="all")
Interface of a cluster finder.
virtual TIterator * CreateLocalIterator() const =0
Create iterator on local trigger.
virtual Int_t PadY() const =0
The y-index of this digit (>=0)
void UnloadRecPoints() const
Definition: AliLoader.h:132
Interface for a digit container.
Bool_t LoadGeometryData(const TString &fileName)
Concrete implementation of AliMUONVDigitStore for real digits.
Top container class for geometry transformations.
Int_t TriggerDigits(Int_t nBoard, const TArrayS *xyPattern, AliMUONVDigitStore &digitStore, Bool_t warn=kTRUE) const
void GetXPattern(TArrayS &array) const
return X pattern array
TTree * TreeD() const
Definition: AliLoader.h:87
AliLoader * GetDetectorLoader(const char *detname)
A rectangle area positioned in plane..
Definition: AliMpArea.h:20
virtual void Clear(Option_t *opt="")
Clear ourselves (i.e. Reset)
Class with MUON reconstruction parameters.
Int_t LoadRecPoints(Option_t *opt="")
Definition: AliLoader.h:113
Bool_t LoadMapping(Bool_t segmentationOnly)
Definition: AliMUONCDB.cxx:502
virtual Float_t Charge() const =0
The charge of this digit, calibrated or not depending on IsCalibrated()
Implementation of VClusterStore.
(Legacy) implementation of AliMUONVDigitStore
virtual Int_t DetElemId() const =0
The detection element this digit belongs to.
virtual TIterator * CreateIterator() const
Return an iterator to loop over the whole store.
Base class of a trigger information store.
virtual TIterator * CreateIterator() const =0
Create an iterator to loop over all our digits.
const Int_t kNtrigChambers
Definition: TestRecPoints.C:68
static AliRunLoader * Open(const char *filename="galice.root", const char *eventfoldername=AliConfig::GetDefaultEventFolderName(), Option_t *option="READ")
virtual AliMUONVDigit * Add(const AliMUONVDigit &digit, EReplacePolicy replace)
TTree * TreeR() const
Definition: AliLoader.h:89
return clusterFinder
abstract base class for clusters
Reconstructed Local Trigger object.
void UnloadDigits() const
Definition: AliLoader.h:131
virtual AliMUONVDigitStore * Create() const =0
Create an (empty) object of the same concrete class as *this.
virtual AliMUONVStore * Create() const =0
Create an empty copy of this.
Int_t Clusterize(Int_t chamberId, AliMUONVClusterStore &clusterStore, const AliMpArea &area, const AliMUONRecoParam *recoParam=0x0)
Find and add clusters from a given region of a given chamber to the store.
static Int_t NofCathodes()
Return number of cathodes.
Implementation of AliMUONVClusterServer interface.
virtual Int_t Cathode() const =0
Cathode number this digit is on (0 or 1)
Int_t GetPlane(Int_t ch, Int_t cath)
Definition: TestRecPoints.C:75
const Int_t kNtrackChambers
Definition: TestRecPoints.C:71
virtual TIterator * CreateIterator() const
Create an iterator to loop over all our digits.
virtual Int_t GetNDigits() const =0
Return number of associated digits.
virtual TIterator * CreateIterator() const
Create an iterator to loop over all our digits.
virtual void Clear(Option_t *opt="")
Clear ourselves (i.e. Reset)
Int_t GetNumberOfEvents()
TString baseDir
Definition: UnitTest.C:38
virtual Int_t PadX() const =0
The x-index of this digit (>=0)
void SetRun(Int_t run)
void SetDefaultStorage(const char *dbString)
static Int_t GetChamberId(UInt_t uniqueID)
Return chamber id (0..), part of the uniqueID.
Int_t GetEvent(Int_t evno)
virtual void Clear(Option_t *opt="")
Clear container.
void TestRecPoints(TString baseDir=".", TString outDir=".", Float_t adcCut=10., Int_t whatToTest=kTrackTrig, Int_t runNumber=0, TString cdbStorage="local://$ALICE_ROOT/OCDB")
Definition: TestRecPoints.C:79
void UseDigits(TIter &next, AliMUONVDigitStore *digitStore=0x0)
Specify an iterator to loop over the digits needed to perform our job.
ABC of a MUON digit.
Definition: AliMUONVDigit.h:18
static Int_t runNumber
Definition: pdc06_config.C:126
static Int_t GetDetElemId(UInt_t uniqueID)
Return detection element id, part of the uniqueID.
virtual void Clear(Option_t *opt="")=0
Clear ourselves (i.e. Reset)
A basic pre-cluster finder.
void ClusterSize(TList &, AliMUONVDigit *, Int_t &, Int_t)
Int_t LoadDigits(Option_t *opt="")
Definition: AliLoader.h:106
static AliCDBManager * Instance(TMap *entryCache=NULL, Int_t run=-1)
virtual Double_t GetCharge() const =0
Set the cluster charge.
static TGeoManager * GetGeometry()
virtual Bool_t Connect(TTree &tree, Bool_t alone=kTRUE) const
Connect us to a TTree (only valid if CanConnect()==kTRUE)
static void LoadGeometry(const char *geomFileName=NULL)
static Int_t NofTrackingChambers()
Return number of tracking chambers.
Int_t LoCircuit() const
Return Circuit number.
Cluster finder in MUON arm of ALICE.
AliMUONRecoParam * LoadRecoParam()
Definition: AliMUONCDB.cxx:535
static Int_t NofTriggerChambers()
virtual void Print(Option_t *option="") const