]>
Commit | Line | Data |
---|---|---|
cc727848 | 1 | /************************************************************************** |
6f6a3f3b | 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 | **************************************************************************/ | |
cc727848 | 15 | |
16 | /* $Id: TestRecPoints.C 23207 2007-12-20 09:59:20Z ivana $ */ | |
17 | ||
18 | /// \ingroup macros | |
19 | /// \file TestRecPoints.C | |
20 | /// \brief This macro is to be used to check the trigger and tracking chamber behaviour | |
21 | /// through the RecPoints analysis | |
22 | /// | |
23 | /// \author D. Stocco (Torino) | |
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" | |
35be7ed7 | 41 | #include "AliCDBEntry.h" |
42 | #include "AliCDBPath.h" | |
cc727848 | 43 | |
44 | // tracker | |
45 | #include "AliGeomManager.h" | |
46 | ||
47 | // MUON | |
48 | #include "AliMpCDB.h" | |
49 | #include "AliMpDDL.h" | |
50 | #include "AliMpDDLStore.h" | |
51 | #include "AliMpConstants.h" | |
52 | ||
53 | // trigger | |
54 | #include "AliMUONVTriggerStore.h" | |
55 | #include "AliMUONDigitStoreV1.h" | |
56 | #include "AliMUONDigitMaker.h" | |
57 | #include "AliMUONLocalTrigger.h" | |
58 | ||
59 | // tracker | |
60 | #include "AliMUONClusterStoreV2.h" | |
61 | #include "AliMUONClusterFinderMLEM.h" | |
62 | #include "AliMUONPreClusterFinder.h" | |
63 | #include "AliMUONSimpleClusterServer.h" | |
64 | #include "AliMUONGeometryTransformer.h" | |
65 | #include "AliMUONVDigitStore.h" | |
6f6a3f3b | 66 | #include "AliMUONDigitStoreV2R.h" |
cc727848 | 67 | #include "AliMUONVCluster.h" |
35be7ed7 | 68 | #include "AliMUONRecoParam.h" |
cc727848 | 69 | |
70 | #endif | |
71 | ||
72 | const Int_t kNtrigChambers = AliMpConstants::NofTriggerChambers(); | |
73 | const Int_t kNcathodes = AliMpConstants::NofCathodes(); | |
74 | ||
75 | const Int_t kNtrackChambers = AliMpConstants::NofTrackingChambers(); | |
76 | ||
77 | enum {kOnlyTracker, kOnlyTrigger, kTrackTrig}; | |
78 | ||
79 | Int_t GetPlane(Int_t ch, Int_t cath){return kNcathodes * ch + cath;} | |
80 | void ClusterSize(TList&, AliMUONVDigit*, Int_t&, Int_t); | |
81 | ||
82 | // Main Method | |
162637e4 | 83 | void TestRecPoints(TString baseDir=".", TString outDir=".", Float_t adcCut = 10., Int_t whatToTest=kTrackTrig, Int_t runNumber=0, TString cdbStorage="local://$ALICE_ROOT/OCDB") |
cc727848 | 84 | { |
85 | const Int_t kNplanes = kNtrigChambers * kNcathodes; | |
86 | const Int_t kNslats = 18; | |
6f6a3f3b | 87 | |
cc727848 | 88 | TList histoListTrig; |
6f6a3f3b | 89 | |
cc727848 | 90 | TH1F* trigStripMult[kNplanes]; |
91 | TH1F* trigClusterSize[kNplanes]; | |
92 | TH1F* trigClusterMult[kNplanes]; | |
93 | TString histoName, histoTitle; | |
6f6a3f3b | 94 | |
cc727848 | 95 | if(whatToTest!=kOnlyTracker){ |
96 | for(Int_t ch=0; ch<kNtrigChambers; ch++){ | |
97 | for(Int_t cath=0; cath<kNcathodes; cath++){ | |
6f6a3f3b | 98 | Int_t iplane = GetPlane(ch, cath); |
99 | histoName = "stripMult"; | |
100 | histoName.Append(Form("Ch%iCath%i",ch,cath)); | |
101 | histoTitle = Form("Strip multiplicity Ch %i Cathode %i", ch, cath); | |
102 | trigStripMult[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), kNslats, 0.-0.5, (Float_t)kNslats - 0.5); | |
103 | trigStripMult[iplane]->SetXTitle("slat"); | |
104 | histoListTrig.Add(trigStripMult[iplane]); | |
105 | ||
106 | histoName = "clusterSize"; | |
107 | histoName.Append(Form("Ch%iCath%i",ch,cath)); | |
108 | histoTitle = Form("Cluster size Ch %i Cathode %i", ch, cath); | |
109 | trigClusterSize[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), 10, 0.-0.5, (Float_t)10 - 0.5); | |
110 | trigClusterSize[iplane]->SetXTitle("cluster size"); | |
111 | histoListTrig.Add(trigClusterSize[iplane]); | |
112 | ||
113 | histoName = "clusterMult"; | |
114 | histoName.Append(Form("Ch%iCath%i",ch,cath)); | |
115 | histoTitle = Form("Cluster multiplicity Ch %i Cathode %i", ch, cath); | |
116 | trigClusterMult[iplane] = new TH1F(histoName.Data(), histoTitle.Data(), kNslats, 0.-0.5, (Float_t)kNslats - 0.5); | |
117 | trigClusterMult[iplane]->SetXTitle("slat"); | |
118 | histoListTrig.Add(trigClusterMult[iplane]); | |
cc727848 | 119 | } |
120 | } | |
121 | } | |
6f6a3f3b | 122 | |
cc727848 | 123 | const Int_t kMaxDetElem = 26; |
124 | TList histoListTrack; | |
125 | TH1F* trackClusterSize[kNtrackChambers]; | |
126 | TH1F* trackClusterMult[kNtrackChambers]; | |
127 | TH1F* trackClusterYield[kNtrackChambers]; | |
6f6a3f3b | 128 | |
cc727848 | 129 | if(whatToTest!=kOnlyTrigger){ |
130 | for(Int_t ch=0; ch<kNtrackChambers; ch++){ | |
131 | histoName = "clusterSize"; | |
132 | histoName.Append(Form("Ch%i",ch)); | |
133 | histoTitle = Form("Cluster size Ch %i", ch); | |
134 | trackClusterSize[ch] = new TH1F(histoName.Data(), histoTitle.Data(), kMaxDetElem, 0.-0.5, (Float_t)kMaxDetElem - 0.5); | |
135 | trackClusterSize[ch]->SetXTitle("cluster size"); | |
136 | histoListTrack.Add(trackClusterSize[ch]); | |
6f6a3f3b | 137 | |
cc727848 | 138 | histoName = "clusterMult"; |
139 | histoName.Append(Form("Ch%i",ch)); | |
140 | histoTitle = Form("Cluster multiplicity vs. Charge Ch %i", ch); | |
141 | trackClusterMult[ch] = new TH1F(histoName.Data(), histoTitle.Data(), 200, 0, 15000.); | |
6f6a3f3b | 142 | trackClusterMult[ch]->SetXTitle("Charge (ADC)"); |
cc727848 | 143 | histoListTrack.Add(trackClusterMult[ch]); |
6f6a3f3b | 144 | |
cc727848 | 145 | histoName = "clusterYield"; |
146 | histoName.Append(Form("Ch%i",ch)); | |
147 | histoTitle = Form("Cluster yield vs. DetElem Ch %i", ch); | |
148 | trackClusterYield[ch] = new TH1F(histoName.Data(), histoTitle.Data(), kMaxDetElem, 0.-0.5, (Float_t)kMaxDetElem - 0.5); | |
149 | trackClusterYield[ch]->SetXTitle("Detector element #"); | |
150 | histoListTrack.Add(trackClusterYield[ch]); | |
151 | } | |
152 | } | |
153 | ||
6f6a3f3b | 154 | |
cc727848 | 155 | // Creating Run Loader and opening RecPoints |
156 | TString filename = baseDir.Data(); | |
157 | filename.Append("/galice.root"); | |
158 | AliRunLoader * RunLoader = AliRunLoader::Open(filename.Data(),"MUONLoader","UPDATE"); | |
159 | if (RunLoader ==0x0) { | |
160 | printf(">>> Error : Error Opening %s file \n",filename.Data()); | |
161 | return; | |
162 | } | |
6f6a3f3b | 163 | |
cc727848 | 164 | // Loading MUON subsystem |
165 | AliLoader* MUONLoader = RunLoader->GetDetectorLoader("MUON"); | |
166 | if(whatToTest!=kOnlyTracker) MUONLoader->LoadRecPoints("READ"); | |
167 | if(whatToTest!=kOnlyTrigger) MUONLoader->LoadDigits("READ"); | |
6f6a3f3b | 168 | |
cc727848 | 169 | TTree *treeR = 0x0, *treeD = 0x0; |
170 | AliMUONVTriggerStore* trigStore = 0x0; | |
171 | AliMUONLocalTrigger* locTrg = 0x0; | |
6f6a3f3b | 172 | |
cc727848 | 173 | AliMUONVDigitStore* digitStoreTrack = 0x0; |
6f6a3f3b | 174 | AliMUONDigitStoreV2R digitStoreTrackCut; |
cc727848 | 175 | AliMUONVCluster* cluster = 0x0; |
6f6a3f3b | 176 | |
cc727848 | 177 | // Load segmentation |
178 | AliCDBManager::Instance()->SetDefaultStorage(cdbStorage.Data()); | |
179 | AliCDBManager::Instance()->SetRun(runNumber); | |
6f6a3f3b | 180 | |
cc727848 | 181 | AliMpCDB::LoadDDLStore(); |
6f6a3f3b | 182 | |
cc727848 | 183 | AliMUONGeometryTransformer* transformer = 0x0; |
6f6a3f3b | 184 | |
35be7ed7 | 185 | AliMUONRecoParam* recoParam = 0x0; |
6f6a3f3b | 186 | |
cc727848 | 187 | AliMUONClusterStoreV2 clusterStore; |
188 | AliMUONVClusterFinder* clusterFinder = 0x0; | |
189 | ||
190 | AliMUONSimpleClusterServer* clusterServer = 0x0; | |
6f6a3f3b | 191 | |
cc727848 | 192 | Int_t firstChamber(0); |
193 | Int_t lastChamber(9); | |
6f6a3f3b | 194 | |
cc727848 | 195 | if(whatToTest!=kOnlyTrigger){ |
196 | // Import TGeo geometry | |
197 | TString geoFilename = baseDir.Data(); | |
198 | geoFilename.Append("/geometry.root"); | |
199 | if ( ! AliGeomManager::GetGeometry() ) { | |
200 | AliGeomManager::LoadGeometry(geoFilename); | |
201 | if (! AliGeomManager::GetGeometry() ) { | |
6f6a3f3b | 202 | printf(">>> Error : getting geometry from file %s failed\n", geoFilename.Data()); |
203 | return; | |
cc727848 | 204 | } |
205 | } | |
206 | transformer = new AliMUONGeometryTransformer(); | |
207 | // Load geometry data | |
208 | transformer->LoadGeometryData(); | |
35be7ed7 | 209 | // Load reconstruction parameters |
210 | AliCDBPath path("MUON","Calib","RecoParam"); | |
211 | AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath()); | |
212 | if(entry) { | |
213 | recoParam = dynamic_cast<AliMUONRecoParam*>(entry->GetObject()); | |
214 | entry->SetOwner(0); | |
215 | AliCDBManager::Instance()->UnloadFromCache(path.GetPath()); | |
216 | } | |
217 | if (!recoParam) { | |
218 | printf("Couldn't find RecoParam object in OCDB: create default one"); | |
219 | recoParam = AliMUONRecoParam::GetLowFluxParam(); | |
220 | } | |
221 | recoParam->Print("FULL"); | |
cc727848 | 222 | clusterFinder = new AliMUONClusterFinderMLEM(kFALSE,new AliMUONPreClusterFinder); |
0f211038 | 223 | clusterServer = new AliMUONSimpleClusterServer(clusterFinder,*transformer); |
cc727848 | 224 | } |
6f6a3f3b | 225 | |
cc727848 | 226 | AliMUONDigitStoreV1 digitStore; |
227 | AliMUONVDigit* mDigit; | |
6f6a3f3b | 228 | |
cc727848 | 229 | Int_t clusterSize; |
6f6a3f3b | 230 | |
cc727848 | 231 | AliMUONDigitMaker digitMaker; |
232 | TList digitsList[kNplanes]; | |
6f6a3f3b | 233 | |
cc727848 | 234 | Int_t nevents = RunLoader->GetNumberOfEvents(); |
235 | Int_t analysisFrac = nevents/10+1; | |
6f6a3f3b | 236 | |
cc727848 | 237 | printf("\nNumber of events = %i\n\n",nevents); |
6f6a3f3b | 238 | |
cc727848 | 239 | for(Int_t ievent=0; ievent<nevents; ievent++){ |
240 | if(ievent%analysisFrac==0) printf("Analysing event = %i\n",ievent); | |
241 | RunLoader->GetEvent(ievent); | |
242 | if(whatToTest!=kOnlyTracker){ | |
243 | digitStore.Clear(); | |
244 | treeR = MUONLoader->TreeR(); | |
245 | trigStore = AliMUONVTriggerStore::Create(*treeR); | |
6f6a3f3b | 246 | |
cc727848 | 247 | if ( trigStore == 0x0 ) continue; |
248 | trigStore->Clear(); | |
249 | trigStore->Connect(*treeR); | |
250 | treeR->GetEvent(0); | |
6f6a3f3b | 251 | |
cc727848 | 252 | TIter nextLocal(trigStore->CreateLocalIterator()); |
253 | while ( (locTrg = static_cast<AliMUONLocalTrigger*>( nextLocal() )) != NULL ) | |
254 | { | |
6f6a3f3b | 255 | TArrayS xyPattern[2]; |
256 | locTrg->GetXPattern(xyPattern[0]); | |
257 | locTrg->GetYPattern(xyPattern[1]); | |
258 | ||
259 | Int_t nBoard = locTrg->LoCircuit(); | |
260 | digitMaker.TriggerDigits(nBoard, xyPattern, digitStore); | |
cc727848 | 261 | } |
6f6a3f3b | 262 | |
cc727848 | 263 | TIter next(digitStore.CreateIterator()); |
6f6a3f3b | 264 | |
cc727848 | 265 | while ( ( mDigit = static_cast<AliMUONVDigit*>(next()) ) ) |
266 | { | |
6f6a3f3b | 267 | Int_t detElemId = mDigit->DetElemId(); |
268 | Int_t ch = detElemId/100 - 11; | |
269 | Int_t cathode = mDigit->Cathode(); | |
270 | Int_t slat = detElemId%100; | |
271 | Int_t iplane = GetPlane(ch, cathode); | |
272 | trigStripMult[iplane]->Fill(slat); | |
273 | digitsList[iplane].Add(mDigit); | |
cc727848 | 274 | } // loop on digits |
6f6a3f3b | 275 | |
cc727848 | 276 | for(Int_t iplane=0; iplane<kNplanes; iplane++){ |
6f6a3f3b | 277 | while(digitsList[iplane].GetEntries()){ |
278 | clusterSize=1; | |
279 | mDigit = (AliMUONVDigit*)digitsList[iplane].At(0); | |
280 | digitsList[iplane].Remove(mDigit); | |
281 | ClusterSize(digitsList[iplane],mDigit,clusterSize,1); | |
282 | //if(clusterSize>1) printf("Cluster size = %i\n\n",clusterSize); | |
283 | trigClusterSize[iplane]->Fill(clusterSize); | |
284 | ||
285 | Int_t detElemId = mDigit->DetElemId(); | |
286 | Int_t slat = detElemId%100; | |
287 | trigClusterMult[iplane]->Fill(slat); | |
288 | } // loop o sorted digits | |
cc727848 | 289 | } // loop on planes |
290 | } // trigger part | |
6f6a3f3b | 291 | |
cc727848 | 292 | if(whatToTest!=kOnlyTrigger){ |
293 | clusterStore.Clear(); | |
294 | treeD = MUONLoader->TreeD(); | |
295 | digitStoreTrack = AliMUONVDigitStore::Create(*treeD); | |
6f6a3f3b | 296 | |
cc727848 | 297 | if ( digitStoreTrack == 0x0 ) continue; |
298 | digitStoreTrack->Clear(); | |
6f6a3f3b | 299 | digitStoreTrackCut.Clear(); |
cc727848 | 300 | digitStoreTrack->Connect(*treeD); |
301 | treeD->GetEvent(0); | |
6f6a3f3b | 302 | |
303 | // Cut low charge channels (pedestal subtraction) | |
304 | TIter nextDigitTrack(digitStoreTrack->CreateIterator()); | |
305 | while ( ( mDigit = static_cast<AliMUONVDigit*>(nextDigitTrack()) ) ) | |
306 | { | |
307 | //printf("Digit class = %s",mDigit->ClassName()); | |
308 | Float_t charge = mDigit->Charge(); | |
309 | if(charge<adcCut) continue; | |
310 | digitStoreTrackCut.Add(*mDigit, AliMUONVDigitStore::kDeny); | |
311 | } // loop on digits | |
312 | ||
0f211038 | 313 | TIter nextDigitTrackCut(digitStoreTrackCut.CreateIterator()); |
2e2d0c44 | 314 | clusterServer->UseDigits(nextDigitTrackCut,&digitStoreTrackCut); |
6f6a3f3b | 315 | |
cc727848 | 316 | for (Int_t ch = firstChamber; ch <= lastChamber; ++ch ){ |
35be7ed7 | 317 | clusterServer->Clusterize(ch, clusterStore, AliMpArea(),recoParam); |
cc727848 | 318 | } |
6f6a3f3b | 319 | |
cc727848 | 320 | TIter next(clusterStore.CreateIterator()); |
321 | while ( ( cluster = static_cast<AliMUONVCluster*>(next()) ) ) | |
322 | { | |
6f6a3f3b | 323 | Float_t charge = cluster->GetCharge(); |
324 | if(charge==0.) continue; | |
325 | Int_t ch = cluster->GetChamberId(); | |
326 | Int_t npads = cluster->GetNDigits(); | |
327 | Int_t detElemId = cluster->GetDetElemId()%100; | |
328 | //printf("charge = %f pads = %i detElemId = %i\n",charge,npads,detElemId); | |
329 | ||
330 | trackClusterSize[ch]->Fill(npads); | |
331 | trackClusterMult[ch]->Fill(charge); | |
332 | trackClusterYield[ch]->Fill(detElemId); | |
cc727848 | 333 | } // loop on clusters |
334 | } // tracker part | |
335 | } // loop on events | |
6f6a3f3b | 336 | |
cc727848 | 337 | MUONLoader->UnloadRecPoints(); |
338 | MUONLoader->UnloadDigits(); | |
339 | RunLoader->UnloadAll(); | |
340 | delete RunLoader; | |
341 | TString outFileName = outDir.Data(); | |
342 | outFileName.Append("/outTestRecPoints.root"); | |
343 | TFile* outFile = new TFile(outFileName.Data(), "RECREATE"); | |
344 | TDirectory* dir = 0x0; | |
345 | if(whatToTest!=kOnlyTracker){ | |
346 | outFile->cd(); | |
347 | dir = outFile->mkdir("Trigger"); | |
348 | dir->cd(); | |
349 | histoListTrig.Write(); | |
350 | } | |
351 | if(whatToTest!=kOnlyTrigger){ | |
352 | outFile->cd(); | |
353 | dir = outFile->mkdir("Tracker"); | |
354 | dir->cd(); | |
355 | histoListTrack.Write(); | |
356 | } | |
357 | outFile->Close(); | |
358 | printf("\nSee results in %s\n",outFileName.Data()); | |
359 | } | |
360 | ||
361 | void ClusterSize(TList& list, AliMUONVDigit* refDigit, Int_t& clusterSize, Int_t depthLevel) | |
362 | { | |
363 | AliMUONVDigit* mDigit = 0x0; | |
364 | for(Int_t idigit=0; idigit<list.GetEntries(); idigit++){ | |
365 | mDigit = (AliMUONVDigit*) list.At(idigit); | |
366 | if(mDigit->DetElemId() != refDigit->DetElemId()) continue; | |
367 | Int_t diffX = TMath::Abs(mDigit->PadX() - refDigit->PadX()); | |
368 | Int_t diffY = TMath::Abs(mDigit->PadY() - refDigit->PadY()); | |
369 | if(diffX>1) continue; | |
370 | if(diffY>1) continue; | |
371 | if(diffX*diffY != 0) continue; | |
372 | clusterSize++; | |
373 | list.Remove(mDigit); | |
374 | //printf("DetElemId = %i Level = %i Ref. (%2i,%2i) Found (%2i,%2i)\n",mDigit->DetElemId(),depthLevel,refDigit->PadX(),refDigit->PadY(),mDigit->PadX(),mDigit->PadY()); | |
375 | ClusterSize(list, mDigit, clusterSize, depthLevel+1); | |
376 | Int_t val = idigit + depthLevel - clusterSize; | |
377 | //printf("Level = %i Current digit = %i\t",depthLevel,idigit); | |
378 | idigit = TMath::Max(-1,val); | |
379 | //printf("New digit = %i\n",idigit); | |
380 | } | |
381 | } |