]>
Commit | Line | Data |
---|---|---|
98937d93 | 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 | //----------------------------------------------------------------- | |
17 | // Implementation of the alignment steering class | |
18 | // It provides an access to the track space points | |
19 | // written along the esd tracks. The class enables | |
20 | // the user to plug any track fitter (deriving from | |
21 | // AliTrackFitter class) and minimization fo the | |
22 | // track residual sums (deriving from the AliTrackResiduals). | |
23 | //----------------------------------------------------------------- | |
24 | ||
25 | #include <TChain.h> | |
bfbd5665 | 26 | #include <TFile.h> |
98937d93 | 27 | |
28 | #include "AliAlignmentTracks.h" | |
29 | #include "AliTrackPointArray.h" | |
90dbf5fb | 30 | #include "AliAlignObjParams.h" |
98937d93 | 31 | #include "AliTrackFitterRieman.h" |
32 | #include "AliTrackResidualsChi2.h" | |
9a1304f0 | 33 | #include "AliESDEvent.h" |
98937d93 | 34 | #include "AliLog.h" |
35 | ||
36 | ClassImp(AliAlignmentTracks) | |
37 | ||
38 | //______________________________________________________________________________ | |
39 | AliAlignmentTracks::AliAlignmentTracks(): | |
40 | fESDChain(0), | |
41 | fPointsFilename("AliTrackPoints.root"), | |
42 | fPointsFile(0), | |
43 | fPointsTree(0), | |
44 | fLastIndex(0), | |
45 | fArrayIndex(0), | |
46 | fIsIndexBuilt(kFALSE), | |
75e3794b | 47 | fAlignObjs(0), |
46ae650f | 48 | fMisalignObjs(0), |
98937d93 | 49 | fTrackFitter(0), |
91264393 | 50 | fMinimizer(0), |
51 | fDoUpdate(kTRUE) | |
98937d93 | 52 | { |
53 | // Default constructor | |
54 | InitIndex(); | |
55 | InitAlignObjs(); | |
56 | } | |
57 | ||
58 | //______________________________________________________________________________ | |
59 | AliAlignmentTracks::AliAlignmentTracks(TChain *esdchain): | |
60 | fESDChain(esdchain), | |
61 | fPointsFilename("AliTrackPoints.root"), | |
62 | fPointsFile(0), | |
63 | fPointsTree(0), | |
64 | fLastIndex(0), | |
65 | fArrayIndex(0), | |
66 | fIsIndexBuilt(kFALSE), | |
75e3794b | 67 | fAlignObjs(0), |
46ae650f | 68 | fMisalignObjs(0), |
98937d93 | 69 | fTrackFitter(0), |
91264393 | 70 | fMinimizer(0), |
71 | fDoUpdate(kTRUE) | |
98937d93 | 72 | { |
73 | // Constructor in the case | |
74 | // the user provides an already | |
75 | // built TChain with ESD trees | |
76 | InitIndex(); | |
77 | InitAlignObjs(); | |
78 | } | |
79 | ||
80 | ||
81 | //______________________________________________________________________________ | |
82 | AliAlignmentTracks::AliAlignmentTracks(const char *esdfilename, const char *esdtreename): | |
75e3794b | 83 | fESDChain(new TChain(esdtreename)), |
98937d93 | 84 | fPointsFilename("AliTrackPoints.root"), |
85 | fPointsFile(0), | |
86 | fPointsTree(0), | |
87 | fLastIndex(0), | |
88 | fArrayIndex(0), | |
89 | fIsIndexBuilt(kFALSE), | |
75e3794b | 90 | fAlignObjs(0), |
46ae650f | 91 | fMisalignObjs(0), |
98937d93 | 92 | fTrackFitter(0), |
91264393 | 93 | fMinimizer(0), |
94 | fDoUpdate(kTRUE) | |
98937d93 | 95 | { |
96 | // Constructor in the case | |
97 | // the user provides a single ESD file | |
98 | // or a directory containing ESD files | |
98937d93 | 99 | fESDChain->Add(esdfilename); |
100 | ||
101 | InitIndex(); | |
102 | InitAlignObjs(); | |
103 | } | |
104 | ||
98937d93 | 105 | |
106 | //______________________________________________________________________________ | |
107 | AliAlignmentTracks::~AliAlignmentTracks() | |
108 | { | |
109 | // Destructor | |
110 | if (fESDChain) delete fESDChain; | |
111 | ||
112 | DeleteIndex(); | |
113 | DeleteAlignObjs(); | |
114 | ||
115 | delete fTrackFitter; | |
116 | delete fMinimizer; | |
117 | ||
118 | if (fPointsFile) fPointsFile->Close(); | |
119 | } | |
120 | ||
121 | //______________________________________________________________________________ | |
122 | void AliAlignmentTracks::AddESD(TChain *esdchain) | |
123 | { | |
124 | // Add a chain with ESD files | |
125 | if (fESDChain) | |
126 | fESDChain->Add(esdchain); | |
127 | else | |
128 | fESDChain = esdchain; | |
129 | } | |
130 | ||
131 | //______________________________________________________________________________ | |
132 | void AliAlignmentTracks::AddESD(const char *esdfilename, const char *esdtreename) | |
133 | { | |
134 | // Add a single file or | |
135 | // a directory to the chain | |
136 | // with the ESD files | |
137 | if (fESDChain) | |
138 | fESDChain->AddFile(esdfilename,TChain::kBigNumber,esdtreename); | |
139 | else { | |
140 | fESDChain = new TChain(esdtreename); | |
141 | fESDChain->Add(esdfilename); | |
142 | } | |
143 | } | |
144 | ||
babbb915 | 145 | //________________________________________________________________________ |
146 | ||
147 | void AliAlignmentTracks::ProcessESD(Bool_t onlyITS, | |
148 | Int_t minITSpts, | |
149 | Bool_t cuts, | |
150 | Float_t minMom,Float_t maxMom, | |
151 | Float_t minAbsSinPhi,Float_t maxAbsSinPhi, | |
152 | Float_t minSinTheta,Float_t maxSinTheta) | |
98937d93 | 153 | { |
154 | // Analyzes and filters ESD tracks | |
155 | // Stores the selected track space points | |
156 | // into the output file | |
157 | ||
158 | if (!fESDChain) return; | |
159 | ||
9a1304f0 | 160 | AliESDEvent *esd = new AliESDEvent(); |
161 | esd->ReadFromTree(fESDChain); | |
cf0f66c2 | 162 | AliESDfriend *esdf = 0; |
163 | fESDChain->SetBranchStatus("ESDfriend*",1); | |
164 | fESDChain->SetBranchAddress("ESDfriend.",&esdf); | |
98937d93 | 165 | |
166 | // Open the output file | |
167 | if (fPointsFilename.Data() == "") { | |
168 | AliWarning("Incorrect output filename!"); | |
169 | return; | |
170 | } | |
171 | ||
172 | TFile *pointsFile = TFile::Open(fPointsFilename,"RECREATE"); | |
173 | if (!pointsFile || !pointsFile->IsOpen()) { | |
174 | AliWarning(Form("Can't open %s !",fPointsFilename.Data())); | |
175 | return; | |
176 | } | |
177 | ||
178 | TTree *pointsTree = new TTree("spTree", "Tree with track space point arrays"); | |
15e85efa | 179 | const AliTrackPointArray *array = 0; |
babbb915 | 180 | AliTrackPointArray *array2 = 0; |
181 | if(onlyITS) { // only ITS AliTrackPoints | |
182 | pointsTree->Branch("SP","AliTrackPointArray", &array2); | |
183 | } else { | |
184 | pointsTree->Branch("SP","AliTrackPointArray", &array); | |
185 | } | |
186 | ||
98937d93 | 187 | Int_t ievent = 0; |
188 | while (fESDChain->GetEntry(ievent++)) { | |
189 | if (!esd) break; | |
cf0f66c2 | 190 | |
191 | esd->SetESDfriend(esdf); //Attach the friend to the ESD | |
192 | ||
98937d93 | 193 | Int_t ntracks = esd->GetNumberOfTracks(); |
194 | for (Int_t itrack=0; itrack < ntracks; itrack++) { | |
195 | AliESDtrack * track = esd->GetTrack(itrack); | |
196 | if (!track) continue; | |
98937d93 | 197 | |
babbb915 | 198 | if(track->GetNcls(0) < minITSpts) continue; |
199 | if(cuts) { | |
200 | if(track->GetP()<minMom || track->GetP()>maxMom) continue; | |
201 | Float_t abssinphi = TMath::Abs(TMath::Sin(track->GetAlpha()+TMath::ASin(track->GetSnp()))); | |
202 | if(abssinphi<minAbsSinPhi || abssinphi>maxAbsSinPhi) continue; | |
203 | Float_t sintheta = TMath::Sin(0.5*TMath::Pi()-TMath::ATan(track->GetTgl())); | |
204 | if(sintheta<minSinTheta || sintheta>maxSinTheta) continue; | |
205 | } | |
206 | // UInt_t status = AliESDtrack::kITSpid; | |
207 | // status|=AliESDtrack::kTPCpid; | |
208 | // status|=AliESDtrack::kTRDpid; | |
209 | // if ((track->GetStatus() & status) != status) continue; | |
98937d93 | 210 | |
babbb915 | 211 | AliTrackPoint point; |
98937d93 | 212 | array = track->GetTrackPointArray(); |
babbb915 | 213 | |
214 | if(onlyITS) { | |
215 | array2 = new AliTrackPointArray(track->GetNcls(0)); | |
216 | Int_t ipt,volId,modId,layerId; | |
217 | Int_t jpt=0; | |
218 | for(ipt=0; ipt<array->GetNPoints(); ipt++) { | |
219 | array->GetPoint(point,ipt); | |
220 | volId = point.GetVolumeID(); | |
221 | layerId = AliGeomManager::VolUIDToLayer(volId,modId); | |
222 | if(layerId<7) { | |
223 | array2->AddPoint(jpt,&point); | |
224 | jpt++; | |
225 | } | |
226 | } | |
227 | } // end if(onlyITS) | |
228 | ||
98937d93 | 229 | pointsTree->Fill(); |
230 | } | |
231 | } | |
232 | ||
233 | if (!pointsTree->Write()) { | |
234 | AliWarning("Can't write the tree with track point arrays!"); | |
235 | return; | |
236 | } | |
237 | ||
238 | pointsFile->Close(); | |
239 | } | |
240 | ||
babbb915 | 241 | //______________________________________________________________________________ |
242 | void AliAlignmentTracks::ProcessESDCosmics(Bool_t onlyITS, | |
9a1304f0 | 243 | Int_t minITSpts,Float_t maxMatchingAngle, |
babbb915 | 244 | Bool_t cuts, |
245 | Float_t minMom,Float_t maxMom, | |
246 | Float_t minAbsSinPhi,Float_t maxAbsSinPhi, | |
247 | Float_t minSinTheta,Float_t maxSinTheta) | |
248 | { | |
249 | // Analyzes and filters ESD tracks | |
250 | // Merges inward and outward tracks in one single track | |
251 | // Stores the selected track space points | |
252 | // into the output file | |
253 | ||
254 | if (!fESDChain) return; | |
255 | ||
9a1304f0 | 256 | AliESDEvent *esd = new AliESDEvent(); |
257 | esd->ReadFromTree(fESDChain); | |
babbb915 | 258 | AliESDfriend *esdf = 0; |
259 | fESDChain->SetBranchStatus("ESDfriend*",1); | |
260 | fESDChain->SetBranchAddress("ESDfriend.",&esdf); | |
261 | ||
262 | // Open the output file | |
263 | if (fPointsFilename.Data() == "") { | |
264 | AliWarning("Incorrect output filename!"); | |
265 | return; | |
266 | } | |
267 | ||
268 | TFile *pointsFile = TFile::Open(fPointsFilename,"RECREATE"); | |
269 | if (!pointsFile || !pointsFile->IsOpen()) { | |
270 | AliWarning(Form("Can't open %s !",fPointsFilename.Data())); | |
271 | return; | |
272 | } | |
273 | ||
274 | TTree *pointsTree = new TTree("spTree", "Tree with track space point arrays"); | |
275 | const AliTrackPointArray *array = 0; | |
276 | AliTrackPointArray *array2 = 0; | |
277 | pointsTree->Branch("SP","AliTrackPointArray", &array2); | |
278 | ||
279 | Int_t ievent = 0; | |
280 | while (fESDChain->GetEntry(ievent++)) { | |
281 | if (!esd) break; | |
282 | ||
283 | esd->SetESDfriend(esdf); //Attach the friend to the ESD | |
284 | ||
285 | Int_t ntracks = esd->GetNumberOfTracks(); | |
9a1304f0 | 286 | if(ntracks<2) continue; |
287 | Int_t *goodtracksArray = new Int_t[ntracks]; | |
288 | Float_t *phiArray = new Float_t[ntracks]; | |
289 | Float_t *thetaArray = new Float_t[ntracks]; | |
babbb915 | 290 | Int_t ngt=0; |
babbb915 | 291 | for (Int_t itrack=0; itrack < ntracks; itrack++) { |
292 | AliESDtrack * track = esd->GetTrack(itrack); | |
293 | if (!track) continue; | |
294 | ||
295 | if(track->GetNcls(0) < minITSpts) continue; | |
9a1304f0 | 296 | Float_t phi = track->GetAlpha()+TMath::ASin(track->GetSnp()); |
297 | Float_t theta = 0.5*TMath::Pi()-TMath::ATan(track->GetTgl()); | |
babbb915 | 298 | if(cuts) { |
299 | if(track->GetP()<minMom || track->GetP()>maxMom) continue; | |
9a1304f0 | 300 | Float_t abssinphi = TMath::Abs(TMath::Sin(phi)); |
babbb915 | 301 | if(abssinphi<minAbsSinPhi || abssinphi>maxAbsSinPhi) continue; |
9a1304f0 | 302 | Float_t sintheta = TMath::Sin(theta); |
babbb915 | 303 | if(sintheta<minSinTheta || sintheta>maxSinTheta) continue; |
304 | } | |
9a1304f0 | 305 | goodtracksArray[ngt]=itrack; |
306 | phiArray[ngt]=phi; | |
307 | thetaArray[ngt]=theta; | |
babbb915 | 308 | ngt++; |
309 | } | |
310 | ||
9a1304f0 | 311 | if(ngt<2) { |
312 | delete [] goodtracksArray; goodtracksArray=0; | |
313 | delete [] phiArray; phiArray=0; | |
314 | delete [] thetaArray; thetaArray=0; | |
315 | continue; | |
316 | } | |
317 | ||
318 | // check matching of the two tracks from the muon | |
319 | Float_t min = 10000000.; | |
320 | Int_t good1 = -1, good2 = -1; | |
321 | for(Int_t itr1=0; itr1<ngt-1; itr1++) { | |
322 | for(Int_t itr2=itr1+1; itr2<ngt; itr2++) { | |
323 | Float_t deltatheta = TMath::Abs(TMath::Pi()-thetaArray[itr1]-thetaArray[itr2]); | |
324 | if(deltatheta>maxMatchingAngle) continue; | |
325 | Float_t deltaphi = TMath::Abs(TMath::Abs(phiArray[itr1]-phiArray[itr2])-TMath::Pi()); | |
326 | if(deltaphi>maxMatchingAngle) continue; | |
327 | printf("%f %f %f %f\n",deltaphi,deltatheta,thetaArray[itr1],thetaArray[itr2]); | |
328 | if(deltatheta+deltaphi<min) { | |
329 | min=deltatheta+deltaphi; | |
330 | good1 = goodtracksArray[itr1]; | |
331 | good2 = goodtracksArray[itr2]; | |
332 | } | |
333 | } | |
334 | } | |
335 | ||
336 | delete [] goodtracksArray; goodtracksArray=0; | |
337 | delete [] phiArray; phiArray=0; | |
338 | delete [] thetaArray; thetaArray=0; | |
339 | ||
340 | if(good1<0) continue; | |
babbb915 | 341 | |
9a1304f0 | 342 | AliESDtrack * track1 = esd->GetTrack(good1); |
343 | AliESDtrack * track2 = esd->GetTrack(good2); | |
babbb915 | 344 | |
345 | Int_t ntotpts; | |
346 | if(onlyITS) { | |
347 | ntotpts = track1->GetNcls(0)+track2->GetNcls(0); | |
348 | } else { | |
349 | ntotpts = track1->GetTrackPointArray()->GetNPoints()+ | |
350 | track2->GetTrackPointArray()->GetNPoints(); | |
351 | } | |
352 | array2 = new AliTrackPointArray(ntotpts); | |
353 | AliTrackPoint point; | |
354 | Int_t ipt,volId,modId,layerId; | |
355 | Int_t jpt=0; | |
356 | array = track1->GetTrackPointArray(); | |
357 | for(ipt=0; ipt<array->GetNPoints(); ipt++) { | |
358 | array->GetPoint(point,ipt); | |
359 | if(onlyITS) { | |
360 | volId = point.GetVolumeID(); | |
361 | layerId = AliGeomManager::VolUIDToLayer(volId,modId); | |
362 | if(layerId>6) continue; | |
363 | } | |
364 | array2->AddPoint(jpt,&point); | |
365 | jpt++; | |
366 | } | |
367 | array = track2->GetTrackPointArray(); | |
368 | for(ipt=0; ipt<array->GetNPoints(); ipt++) { | |
369 | array->GetPoint(point,ipt); | |
370 | if(onlyITS) { | |
371 | volId = point.GetVolumeID(); | |
372 | layerId = AliGeomManager::VolUIDToLayer(volId,modId); | |
373 | if(layerId>6) continue; | |
374 | } | |
375 | array2->AddPoint(jpt,&point); | |
376 | jpt++; | |
377 | } | |
378 | ||
379 | pointsTree->Fill(); | |
380 | } | |
381 | ||
382 | if (!pointsTree->Write()) { | |
383 | AliWarning("Can't write the tree with track point arrays!"); | |
384 | return; | |
385 | } | |
386 | ||
387 | pointsFile->Close(); | |
388 | } | |
389 | ||
98937d93 | 390 | //______________________________________________________________________________ |
391 | void AliAlignmentTracks::ProcessESD(TSelector *selector) | |
392 | { | |
393 | AliWarning(Form("ESD processing based on selector is not yet implemented (%p) !",selector)); | |
394 | } | |
395 | ||
396 | //______________________________________________________________________________ | |
397 | void AliAlignmentTracks::BuildIndex() | |
398 | { | |
399 | // Build index of points tree entries | |
400 | // Used for access based on the volume IDs | |
46ae650f | 401 | if (fIsIndexBuilt) return; |
402 | ||
403 | fIsIndexBuilt = kTRUE; | |
98937d93 | 404 | |
33bebad6 | 405 | // Dummy object is created in order |
406 | // to initialize the volume paths | |
90dbf5fb | 407 | AliAlignObjParams alobj; |
33bebad6 | 408 | |
7fa594d2 | 409 | fPointsFile = TFile::Open(fPointsFilename); |
98937d93 | 410 | if (!fPointsFile || !fPointsFile->IsOpen()) { |
411 | AliWarning(Form("Can't open %s !",fPointsFilename.Data())); | |
412 | return; | |
413 | } | |
414 | ||
415 | // AliTrackPointArray* array = new AliTrackPointArray; | |
416 | AliTrackPointArray* array = 0; | |
46ae650f | 417 | fPointsTree = (TTree*) fPointsFile->Get("spTree"); |
418 | if (!fPointsTree) { | |
98937d93 | 419 | AliWarning("No pointsTree found!"); |
420 | return; | |
421 | } | |
46ae650f | 422 | fPointsTree->SetBranchAddress("SP", &array); |
98937d93 | 423 | |
7fa594d2 | 424 | Int_t nArrays = (Int_t)fPointsTree->GetEntries(); |
98937d93 | 425 | for (Int_t iArray = 0; iArray < nArrays; iArray++) |
426 | { | |
46ae650f | 427 | fPointsTree->GetEvent(iArray); |
98937d93 | 428 | if (!array) continue; |
429 | for (Int_t ipoint = 0; ipoint < array->GetNPoints(); ipoint++) { | |
430 | UShort_t volId = array->GetVolumeID()[ipoint]; | |
33bebad6 | 431 | // check if the volId is valid |
25be1e5c | 432 | if (!AliGeomManager::SymName(volId)) { |
b760c02e | 433 | AliError(Form("The volume id %d has no default volume name !", |
6baf24ce | 434 | volId)); |
435 | continue; | |
436 | } | |
98937d93 | 437 | Int_t modId; |
25be1e5c | 438 | Int_t layerId = AliGeomManager::VolUIDToLayer(volId,modId) |
439 | - AliGeomManager::kFirstLayer; | |
98937d93 | 440 | if (!fArrayIndex[layerId][modId]) { |
441 | //first entry for this volume | |
442 | fArrayIndex[layerId][modId] = new TArrayI(1000); | |
443 | } | |
444 | else { | |
445 | Int_t size = fArrayIndex[layerId][modId]->GetSize(); | |
446 | // If needed allocate new size | |
447 | if (fLastIndex[layerId][modId] >= size) | |
448 | fArrayIndex[layerId][modId]->Set(size + 1000); | |
449 | } | |
450 | ||
451 | // Check if the index is already filled | |
452 | Bool_t fillIndex = kTRUE; | |
453 | if (fLastIndex[layerId][modId] != 0) { | |
454 | if ((*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]-1] == iArray) | |
455 | fillIndex = kFALSE; | |
456 | } | |
457 | // Fill the index array and store last filled index | |
458 | if (fillIndex) { | |
459 | (*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]] = iArray; | |
460 | fLastIndex[layerId][modId]++; | |
461 | } | |
462 | } | |
463 | } | |
464 | ||
465 | } | |
466 | ||
98937d93 | 467 | //______________________________________________________________________________ |
468 | void AliAlignmentTracks::InitIndex() | |
469 | { | |
470 | // Initialize the index arrays | |
25be1e5c | 471 | Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer; |
98937d93 | 472 | fLastIndex = new Int_t*[nLayers]; |
473 | fArrayIndex = new TArrayI**[nLayers]; | |
25be1e5c | 474 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
475 | fLastIndex[iLayer] = new Int_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)]; | |
476 | fArrayIndex[iLayer] = new TArrayI*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)]; | |
477 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) { | |
98937d93 | 478 | fLastIndex[iLayer][iModule] = 0; |
479 | fArrayIndex[iLayer][iModule] = 0; | |
480 | } | |
481 | } | |
482 | } | |
483 | ||
484 | //______________________________________________________________________________ | |
485 | void AliAlignmentTracks::ResetIndex() | |
486 | { | |
487 | // Reset the value of the last filled index | |
488 | // Do not realocate memory | |
46ae650f | 489 | |
490 | fIsIndexBuilt = kFALSE; | |
491 | ||
25be1e5c | 492 | for (Int_t iLayer = 0; iLayer < AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer; iLayer++) { |
493 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) { | |
98937d93 | 494 | fLastIndex[iLayer][iModule] = 0; |
495 | } | |
496 | } | |
497 | } | |
498 | ||
499 | //______________________________________________________________________________ | |
500 | void AliAlignmentTracks::DeleteIndex() | |
501 | { | |
502 | // Delete the index arrays | |
503 | // Called by the destructor | |
25be1e5c | 504 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
505 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) { | |
98937d93 | 506 | if (fArrayIndex[iLayer][iModule]) { |
507 | delete fArrayIndex[iLayer][iModule]; | |
508 | fArrayIndex[iLayer][iModule] = 0; | |
509 | } | |
510 | } | |
511 | delete [] fLastIndex[iLayer]; | |
512 | delete [] fArrayIndex[iLayer]; | |
513 | } | |
514 | delete [] fLastIndex; | |
515 | delete [] fArrayIndex; | |
516 | } | |
517 | ||
518 | //______________________________________________________________________________ | |
46ae650f | 519 | Bool_t AliAlignmentTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName) |
98937d93 | 520 | { |
521 | // Read alignment object from a file | |
522 | // To be replaced by a call to CDB | |
46ae650f | 523 | AliWarning(Form("Method not yet implemented (%s in %s) !",arrayName,alignObjFileName)); |
98937d93 | 524 | |
525 | return kFALSE; | |
526 | } | |
527 | ||
528 | //______________________________________________________________________________ | |
529 | void AliAlignmentTracks::InitAlignObjs() | |
530 | { | |
531 | // Initialize the alignment objects array | |
25be1e5c | 532 | Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer; |
98937d93 | 533 | fAlignObjs = new AliAlignObj**[nLayers]; |
25be1e5c | 534 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
535 | fAlignObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)]; | |
536 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) { | |
537 | UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule); | |
90dbf5fb | 538 | fAlignObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE); |
46ae650f | 539 | } |
98937d93 | 540 | } |
541 | } | |
542 | ||
543 | //______________________________________________________________________________ | |
544 | void AliAlignmentTracks::ResetAlignObjs() | |
545 | { | |
546 | // Reset the alignment objects array | |
25be1e5c | 547 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
548 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) | |
98937d93 | 549 | fAlignObjs[iLayer][iModule]->SetPars(0,0,0,0,0,0); |
550 | } | |
551 | } | |
552 | ||
553 | //______________________________________________________________________________ | |
554 | void AliAlignmentTracks::DeleteAlignObjs() | |
555 | { | |
556 | // Delete the alignment objects array | |
25be1e5c | 557 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
558 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) | |
98937d93 | 559 | if (fAlignObjs[iLayer][iModule]) |
560 | delete fAlignObjs[iLayer][iModule]; | |
561 | delete [] fAlignObjs[iLayer]; | |
562 | } | |
563 | delete [] fAlignObjs; | |
46ae650f | 564 | fAlignObjs = 0; |
98937d93 | 565 | } |
566 | ||
25be1e5c | 567 | void AliAlignmentTracks::AlignDetector(AliGeomManager::ELayerID firstLayer, |
568 | AliGeomManager::ELayerID lastLayer, | |
569 | AliGeomManager::ELayerID layerRangeMin, | |
570 | AliGeomManager::ELayerID layerRangeMax, | |
cc345ce3 | 571 | Int_t iterations) |
98937d93 | 572 | { |
cc345ce3 | 573 | // Align detector volumes within |
574 | // a given layer range | |
575 | // (could be whole detector). | |
576 | // Tracks are fitted only within | |
577 | // the range defined by the user. | |
578 | Int_t nModules = 0; | |
98f0f87b | 579 | for (Int_t iLayer = firstLayer; iLayer <= lastLayer; iLayer++) |
25be1e5c | 580 | nModules += AliGeomManager::LayerSize(iLayer); |
cc345ce3 | 581 | TArrayI volIds(nModules); |
582 | ||
583 | Int_t modnum = 0; | |
98f0f87b | 584 | for (Int_t iLayer = firstLayer; iLayer <= lastLayer; iLayer++) { |
25be1e5c | 585 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer); iModule++) { |
586 | UShort_t volId = AliGeomManager::LayerToVolUID(iLayer,iModule); | |
cc345ce3 | 587 | volIds.AddAt(volId,modnum); |
588 | modnum++; | |
589 | } | |
590 | } | |
591 | ||
98937d93 | 592 | while (iterations > 0) { |
cc345ce3 | 593 | AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax); |
594 | iterations--; | |
98937d93 | 595 | } |
596 | } | |
597 | ||
598 | //______________________________________________________________________________ | |
25be1e5c | 599 | void AliAlignmentTracks::AlignLayer(AliGeomManager::ELayerID layer, |
600 | AliGeomManager::ELayerID layerRangeMin, | |
601 | AliGeomManager::ELayerID layerRangeMax, | |
98937d93 | 602 | Int_t iterations) |
603 | { | |
604 | // Align detector volumes within | |
605 | // a given layer. | |
606 | // Tracks are fitted only within | |
607 | // the range defined by the user. | |
25be1e5c | 608 | Int_t nModules = AliGeomManager::LayerSize(layer); |
cc345ce3 | 609 | TArrayI volIds(nModules); |
610 | for (Int_t iModule = 0; iModule < nModules; iModule++) { | |
25be1e5c | 611 | UShort_t volId = AliGeomManager::LayerToVolUID(layer,iModule); |
cc345ce3 | 612 | volIds.AddAt(volId,iModule); |
613 | } | |
614 | ||
98937d93 | 615 | while (iterations > 0) { |
cc345ce3 | 616 | AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax); |
617 | iterations--; | |
618 | } | |
619 | } | |
620 | ||
621 | //______________________________________________________________________________ | |
622 | void AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit, | |
623 | Int_t iterations) | |
624 | { | |
625 | // Align single detector volume to | |
626 | // another volume. | |
627 | // Tracks are fitted only within | |
628 | // the second volume. | |
629 | TArrayI volIds(1); | |
630 | volIds.AddAt(volId,0); | |
631 | TArrayI volIdsFit(1); | |
632 | volIdsFit.AddAt(volIdFit,0); | |
633 | ||
634 | while (iterations > 0) { | |
635 | AlignVolumes(&volIds,&volIdsFit); | |
98937d93 | 636 | iterations--; |
637 | } | |
638 | } | |
639 | ||
640 | //______________________________________________________________________________ | |
cc345ce3 | 641 | void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *volidsfit, |
25be1e5c | 642 | AliGeomManager::ELayerID layerRangeMin, |
643 | AliGeomManager::ELayerID layerRangeMax, | |
46ae650f | 644 | Int_t iterations) |
98937d93 | 645 | { |
cc345ce3 | 646 | // Align a set of detector volumes. |
98937d93 | 647 | // Tracks are fitted only within |
46ae650f | 648 | // the range defined by the user |
649 | // (by layerRangeMin and layerRangeMax) | |
cc345ce3 | 650 | // or within the set of volidsfit |
46ae650f | 651 | // Repeat the procedure 'iterations' times |
98937d93 | 652 | |
cc345ce3 | 653 | Int_t nVolIds = volids->GetSize(); |
654 | if (nVolIds == 0) { | |
655 | AliError("Volume IDs array is empty!"); | |
656 | return; | |
657 | } | |
98937d93 | 658 | |
cc345ce3 | 659 | // Load only the tracks with at least one |
660 | // space point in the set of volume (volids) | |
98937d93 | 661 | BuildIndex(); |
662 | AliTrackPointArray **points; | |
46ae650f | 663 | // Start the iterations |
664 | while (iterations > 0) { | |
cc345ce3 | 665 | Int_t nArrays = LoadPoints(volids, points); |
46ae650f | 666 | if (nArrays == 0) return; |
667 | ||
668 | AliTrackResiduals *minimizer = CreateMinimizer(); | |
669 | minimizer->SetNTracks(nArrays); | |
cc345ce3 | 670 | minimizer->InitAlignObj(); |
46ae650f | 671 | AliTrackFitter *fitter = CreateFitter(); |
672 | for (Int_t iArray = 0; iArray < nArrays; iArray++) { | |
33bebad6 | 673 | if (!points[iArray]) continue; |
46ae650f | 674 | fitter->SetTrackPointArray(points[iArray], kFALSE); |
33bebad6 | 675 | if (fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE) continue; |
46ae650f | 676 | AliTrackPointArray *pVolId,*pTrack; |
677 | fitter->GetTrackResiduals(pVolId,pTrack); | |
678 | minimizer->AddTrackPointArrays(pVolId,pTrack); | |
679 | } | |
680 | minimizer->Minimize(); | |
cc345ce3 | 681 | |
682 | // Update the alignment object(s) | |
91264393 | 683 | if (fDoUpdate) for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) { |
cc345ce3 | 684 | UShort_t volid = (*volids)[iVolId]; |
685 | Int_t iModule; | |
25be1e5c | 686 | AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule); |
687 | AliAlignObj *alignObj = fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule]; | |
cc345ce3 | 688 | *alignObj *= *minimizer->GetAlignObj(); |
689 | alignObj->Print(""); | |
690 | } | |
98937d93 | 691 | |
46ae650f | 692 | UnloadPoints(nArrays, points); |
693 | ||
694 | iterations--; | |
98937d93 | 695 | } |
98937d93 | 696 | } |
697 | ||
698 | //______________________________________________________________________________ | |
cc345ce3 | 699 | Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray** &points) |
98937d93 | 700 | { |
701 | // Load track point arrays with at least | |
cc345ce3 | 702 | // one space point in a given set of detector |
703 | // volumes (array volids). | |
98937d93 | 704 | // Use the already created tree index for |
705 | // fast access. | |
98937d93 | 706 | |
46ae650f | 707 | if (!fPointsTree) { |
cc345ce3 | 708 | AliError("Tree with the space point arrays not initialized!"); |
98937d93 | 709 | points = 0; |
710 | return 0; | |
711 | } | |
712 | ||
cc345ce3 | 713 | Int_t nVolIds = volids->GetSize(); |
714 | if (nVolIds == 0) { | |
715 | AliError("Volume IDs array is empty!"); | |
716 | points = 0; | |
717 | return 0; | |
718 | } | |
719 | ||
720 | Int_t nArrays = 0; | |
721 | for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) { | |
722 | UShort_t volid = (*volids)[iVolId]; | |
723 | Int_t iModule; | |
25be1e5c | 724 | AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule); |
cc345ce3 | 725 | |
726 | // In case of empty index | |
25be1e5c | 727 | if (fLastIndex[iLayer-AliGeomManager::kFirstLayer][iModule] == 0) { |
cc345ce3 | 728 | AliWarning(Form("There are no space-points belonging to the volume which is to be aligned (Volume ID =%d)!",volid)); |
729 | continue; | |
730 | } | |
25be1e5c | 731 | nArrays += fLastIndex[iLayer-AliGeomManager::kFirstLayer][iModule]; |
cc345ce3 | 732 | } |
46ae650f | 733 | |
46ae650f | 734 | if (nArrays == 0) { |
cc345ce3 | 735 | AliError("There are no space-points belonging to all of the volumes which are to be aligned!"); |
98937d93 | 736 | points = 0; |
737 | return 0; | |
738 | } | |
739 | ||
98937d93 | 740 | AliTrackPointArray* array = 0; |
741 | fPointsTree->SetBranchAddress("SP", &array); | |
742 | ||
cc345ce3 | 743 | // Allocate the pointer to the space-point arrays |
98937d93 | 744 | points = new AliTrackPointArray*[nArrays]; |
5caa5567 | 745 | for (Int_t i = 0; i < nArrays; i++) points[i] = 0x0; |
cc345ce3 | 746 | |
747 | // Init the array used to flag already loaded tree entries | |
7fa594d2 | 748 | Bool_t *indexUsed = new Bool_t[(UInt_t)fPointsTree->GetEntries()]; |
cc345ce3 | 749 | for (Int_t i = 0; i < fPointsTree->GetEntries(); i++) |
750 | indexUsed[i] = kFALSE; | |
751 | ||
752 | // Start the loop over the volume ids | |
753 | Int_t iArray = 0; | |
754 | for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) { | |
755 | UShort_t volid = (*volids)[iVolId]; | |
756 | Int_t iModule; | |
25be1e5c | 757 | AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule); |
cc345ce3 | 758 | |
25be1e5c | 759 | Int_t nArraysId = fLastIndex[iLayer-AliGeomManager::kFirstLayer][iModule]; |
760 | TArrayI *index = fArrayIndex[iLayer-AliGeomManager::kFirstLayer][iModule]; | |
cc345ce3 | 761 | AliTrackPoint p; |
762 | ||
763 | for (Int_t iArrayId = 0; iArrayId < nArraysId; iArrayId++) { | |
764 | ||
765 | // Get tree entry | |
766 | Int_t entry = (*index)[iArrayId]; | |
767 | if (indexUsed[entry] == kTRUE) continue; | |
768 | fPointsTree->GetEvent(entry); | |
769 | if (!array) { | |
770 | AliWarning("Wrong space point array index!"); | |
771 | continue; | |
46ae650f | 772 | } |
cc345ce3 | 773 | indexUsed[entry] = kTRUE; |
774 | ||
775 | // Get the space-point array | |
776 | Int_t nPoints = array->GetNPoints(); | |
777 | points[iArray] = new AliTrackPointArray(nPoints); | |
778 | for (Int_t iPoint = 0; iPoint < nPoints; iPoint++) { | |
779 | array->GetPoint(p,iPoint); | |
780 | Int_t modnum; | |
25be1e5c | 781 | AliGeomManager::ELayerID layer = AliGeomManager::VolUIDToLayer(p.GetVolumeID(),modnum); |
33bebad6 | 782 | // check if the layer id is valid |
25be1e5c | 783 | if ((layer < AliGeomManager::kFirstLayer) || |
784 | (layer >= AliGeomManager::kLastLayer)) { | |
6baf24ce | 785 | AliError(Form("Layer index is invalid: %d (%d -> %d) !", |
25be1e5c | 786 | layer,AliGeomManager::kFirstLayer,AliGeomManager::kLastLayer-1)); |
6baf24ce | 787 | continue; |
788 | } | |
25be1e5c | 789 | if ((modnum >= AliGeomManager::LayerSize(layer)) || |
6baf24ce | 790 | (modnum < 0)) { |
791 | AliError(Form("Module number inside layer %d is invalid: %d (0 -> %d)", | |
25be1e5c | 792 | layer,modnum,AliGeomManager::LayerSize(layer))); |
6baf24ce | 793 | continue; |
794 | } | |
cc345ce3 | 795 | |
796 | // Misalignment is introduced here | |
797 | // Switch it off in case of real | |
798 | // alignment job! | |
799 | if (fMisalignObjs) { | |
25be1e5c | 800 | AliAlignObj *misalignObj = fMisalignObjs[layer-AliGeomManager::kFirstLayer][modnum]; |
cc345ce3 | 801 | if (misalignObj) |
802 | misalignObj->Transform(p); | |
803 | } | |
804 | // End of misalignment | |
46ae650f | 805 | |
25be1e5c | 806 | AliAlignObj *alignObj = fAlignObjs[layer-AliGeomManager::kFirstLayer][modnum]; |
cc345ce3 | 807 | alignObj->Transform(p); |
808 | points[iArray]->AddPoint(iPoint,&p); | |
809 | } | |
810 | iArray++; | |
98937d93 | 811 | } |
812 | } | |
813 | ||
cc345ce3 | 814 | |
815 | delete [] indexUsed; | |
816 | ||
98937d93 | 817 | return nArrays; |
818 | } | |
819 | ||
820 | //______________________________________________________________________________ | |
821 | void AliAlignmentTracks::UnloadPoints(Int_t n, AliTrackPointArray **points) | |
822 | { | |
823 | // Unload track point arrays for a given | |
824 | // detector volume | |
825 | for (Int_t iArray = 0; iArray < n; iArray++) | |
826 | delete points[iArray]; | |
827 | delete [] points; | |
828 | } | |
829 | ||
830 | //______________________________________________________________________________ | |
831 | AliTrackFitter *AliAlignmentTracks::CreateFitter() | |
832 | { | |
833 | // Check if the user has already supplied | |
834 | // a track fitter object. | |
835 | // If not, create a default one. | |
836 | if (!fTrackFitter) | |
837 | fTrackFitter = new AliTrackFitterRieman; | |
838 | ||
839 | return fTrackFitter; | |
840 | } | |
841 | ||
842 | //______________________________________________________________________________ | |
843 | AliTrackResiduals *AliAlignmentTracks::CreateMinimizer() | |
844 | { | |
845 | // Check if the user has already supplied | |
846 | // a track residuals minimizer object. | |
847 | // If not, create a default one. | |
848 | if (!fMinimizer) | |
849 | fMinimizer = new AliTrackResidualsChi2; | |
850 | ||
851 | return fMinimizer; | |
852 | } | |
46ae650f | 853 | |
854 | //______________________________________________________________________________ | |
855 | Bool_t AliAlignmentTracks::Misalign(const char *misalignObjFileName, const char* arrayName) | |
856 | { | |
857 | // The method reads from a file a set of AliAlignObj which are | |
858 | // then used to apply misalignments directly on the track | |
859 | // space-points. The method is supposed to be used only for | |
860 | // fast development and debugging of the alignment algorithms. | |
861 | // Be careful not to use it in the case of 'real' alignment | |
862 | // scenario since it will bias the results. | |
863 | ||
864 | // Initialize the misalignment objects array | |
25be1e5c | 865 | Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer; |
46ae650f | 866 | fMisalignObjs = new AliAlignObj**[nLayers]; |
25be1e5c | 867 | for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) { |
868 | fMisalignObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)]; | |
869 | for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) | |
46ae650f | 870 | fMisalignObjs[iLayer][iModule] = 0x0; |
871 | } | |
872 | ||
873 | // Open the misliagnment file and load the array with | |
874 | // misalignment objects | |
875 | TFile* inFile = TFile::Open(misalignObjFileName,"READ"); | |
876 | if (!inFile || !inFile->IsOpen()) { | |
877 | AliError(Form("Could not open misalignment file %s !",misalignObjFileName)); | |
878 | return kFALSE; | |
879 | } | |
880 | ||
881 | TClonesArray* array = ((TClonesArray*) inFile->Get(arrayName)); | |
882 | if (!array) { | |
883 | AliError(Form("Could not find misalignment array %s in the file %s !",arrayName,misalignObjFileName)); | |
884 | inFile->Close(); | |
885 | return kFALSE; | |
886 | } | |
887 | inFile->Close(); | |
888 | ||
889 | // Store the misalignment objects for further usage | |
890 | Int_t nObjs = array->GetEntriesFast(); | |
25be1e5c | 891 | AliGeomManager::ELayerID layerId; // volume layer |
46ae650f | 892 | Int_t modId; // volume ID inside the layer |
893 | for(Int_t i=0; i<nObjs; i++) | |
894 | { | |
895 | AliAlignObj* alObj = (AliAlignObj*)array->UncheckedAt(i); | |
896 | alObj->GetVolUID(layerId,modId); | |
25be1e5c | 897 | fMisalignObjs[layerId-AliGeomManager::kFirstLayer][modId] = alObj; |
46ae650f | 898 | } |
899 | return kTRUE; | |
900 | } |