1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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 //-----------------------------------------------------------------
28 #include "AliAlignmentTracks.h"
29 #include "AliTrackPointArray.h"
30 #include "AliAlignObjAngles.h"
31 #include "AliTrackFitterRieman.h"
32 #include "AliTrackResidualsChi2.h"
36 ClassImp(AliAlignmentTracks)
38 //______________________________________________________________________________
39 AliAlignmentTracks::AliAlignmentTracks():
41 fPointsFilename("AliTrackPoints.root"),
46 fIsIndexBuilt(kFALSE),
51 // Default constructor
56 //______________________________________________________________________________
57 AliAlignmentTracks::AliAlignmentTracks(TChain *esdchain):
59 fPointsFilename("AliTrackPoints.root"),
64 fIsIndexBuilt(kFALSE),
69 // Constructor in the case
70 // the user provides an already
71 // built TChain with ESD trees
77 //______________________________________________________________________________
78 AliAlignmentTracks::AliAlignmentTracks(const char *esdfilename, const char *esdtreename):
79 fPointsFilename("AliTrackPoints.root"),
84 fIsIndexBuilt(kFALSE),
89 // Constructor in the case
90 // the user provides a single ESD file
91 // or a directory containing ESD files
92 fESDChain = new TChain(esdtreename);
93 fESDChain->Add(esdfilename);
99 //______________________________________________________________________________
100 AliAlignmentTracks::AliAlignmentTracks(const AliAlignmentTracks &alignment):
105 AliWarning("Copy constructor not implemented!");
108 //______________________________________________________________________________
109 AliAlignmentTracks& AliAlignmentTracks::operator= (const AliAlignmentTracks& alignment)
111 // Asignment operator
113 if(this==&alignment) return *this;
115 AliWarning("Asignment operator not implemented!");
117 ((TObject *)this)->operator=(alignment);
122 //______________________________________________________________________________
123 AliAlignmentTracks::~AliAlignmentTracks()
126 if (fESDChain) delete fESDChain;
134 if (fPointsFile) fPointsFile->Close();
137 //______________________________________________________________________________
138 void AliAlignmentTracks::AddESD(TChain *esdchain)
140 // Add a chain with ESD files
142 fESDChain->Add(esdchain);
144 fESDChain = esdchain;
147 //______________________________________________________________________________
148 void AliAlignmentTracks::AddESD(const char *esdfilename, const char *esdtreename)
150 // Add a single file or
151 // a directory to the chain
152 // with the ESD files
154 fESDChain->AddFile(esdfilename,TChain::kBigNumber,esdtreename);
156 fESDChain = new TChain(esdtreename);
157 fESDChain->Add(esdfilename);
161 //______________________________________________________________________________
162 void AliAlignmentTracks::ProcessESD()
164 // Analyzes and filters ESD tracks
165 // Stores the selected track space points
166 // into the output file
168 if (!fESDChain) return;
171 fESDChain->SetBranchAddress("ESD",&esd);
173 // Open the output file
174 if (fPointsFilename.Data() == "") {
175 AliWarning("Incorrect output filename!");
179 TFile *pointsFile = TFile::Open(fPointsFilename,"RECREATE");
180 if (!pointsFile || !pointsFile->IsOpen()) {
181 AliWarning(Form("Can't open %s !",fPointsFilename.Data()));
185 TTree *pointsTree = new TTree("spTree", "Tree with track space point arrays");
186 AliTrackPointArray *array = 0;
187 pointsTree->Branch("SP","AliTrackPointArray", &array);
190 while (fESDChain->GetEntry(ievent++)) {
192 Int_t ntracks = esd->GetNumberOfTracks();
193 for (Int_t itrack=0; itrack < ntracks; itrack++) {
194 AliESDtrack * track = esd->GetTrack(itrack);
195 if (!track) continue;
197 // UInt_t status = AliESDtrack::kITSpid;
198 // status|=AliESDtrack::kTPCpid;
199 // status|=AliESDtrack::kTRDpid;
200 // if ((track->GetStatus() & status) != status) continue;
202 if (track->GetP() < 0.5) continue;
204 array = track->GetTrackPointArray();
209 if (!pointsTree->Write()) {
210 AliWarning("Can't write the tree with track point arrays!");
217 //______________________________________________________________________________
218 void AliAlignmentTracks::ProcessESD(TSelector *selector)
220 AliWarning(Form("ESD processing based on selector is not yet implemented (%p) !",selector));
223 //______________________________________________________________________________
224 void AliAlignmentTracks::BuildIndex()
226 // Build index of points tree entries
227 // Used for access based on the volume IDs
228 if (fIsIndexBuilt) return;
230 fIsIndexBuilt = kTRUE;
232 TFile *fPointsFile = TFile::Open(fPointsFilename);
233 if (!fPointsFile || !fPointsFile->IsOpen()) {
234 AliWarning(Form("Can't open %s !",fPointsFilename.Data()));
238 // AliTrackPointArray* array = new AliTrackPointArray;
239 AliTrackPointArray* array = 0;
240 fPointsTree = (TTree*) fPointsFile->Get("spTree");
242 AliWarning("No pointsTree found!");
245 fPointsTree->SetBranchAddress("SP", &array);
247 Int_t nArrays = fPointsTree->GetEntries();
248 for (Int_t iArray = 0; iArray < nArrays; iArray++)
250 fPointsTree->GetEvent(iArray);
251 if (!array) continue;
252 for (Int_t ipoint = 0; ipoint < array->GetNPoints(); ipoint++) {
253 UShort_t volId = array->GetVolumeID()[ipoint];
255 Int_t layerId = AliAlignObj::VolUIDToLayer(volId,modId)
256 - AliAlignObj::kFirstLayer;
257 if (!fArrayIndex[layerId][modId]) {
258 //first entry for this volume
259 fArrayIndex[layerId][modId] = new TArrayI(1000);
262 Int_t size = fArrayIndex[layerId][modId]->GetSize();
263 // If needed allocate new size
264 if (fLastIndex[layerId][modId] >= size)
265 fArrayIndex[layerId][modId]->Set(size + 1000);
268 // Check if the index is already filled
269 Bool_t fillIndex = kTRUE;
270 if (fLastIndex[layerId][modId] != 0) {
271 if ((*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]-1] == iArray)
274 // Fill the index array and store last filled index
276 (*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]] = iArray;
277 fLastIndex[layerId][modId]++;
284 //______________________________________________________________________________
285 void AliAlignmentTracks::InitIndex()
287 // Initialize the index arrays
288 Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
289 fLastIndex = new Int_t*[nLayers];
290 fArrayIndex = new TArrayI**[nLayers];
291 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
292 fLastIndex[iLayer] = new Int_t[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
293 fArrayIndex[iLayer] = new TArrayI*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
294 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
295 fLastIndex[iLayer][iModule] = 0;
296 fArrayIndex[iLayer][iModule] = 0;
301 //______________________________________________________________________________
302 void AliAlignmentTracks::ResetIndex()
304 // Reset the value of the last filled index
305 // Do not realocate memory
307 fIsIndexBuilt = kFALSE;
309 for (Int_t iLayer = 0; iLayer < AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer; iLayer++) {
310 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
311 fLastIndex[iLayer][iModule] = 0;
316 //______________________________________________________________________________
317 void AliAlignmentTracks::DeleteIndex()
319 // Delete the index arrays
320 // Called by the destructor
321 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
322 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
323 if (fArrayIndex[iLayer][iModule]) {
324 delete fArrayIndex[iLayer][iModule];
325 fArrayIndex[iLayer][iModule] = 0;
328 delete [] fLastIndex[iLayer];
329 delete [] fArrayIndex[iLayer];
331 delete [] fLastIndex;
332 delete [] fArrayIndex;
335 //______________________________________________________________________________
336 Bool_t AliAlignmentTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName)
338 // Read alignment object from a file
339 // To be replaced by a call to CDB
340 AliWarning(Form("Method not yet implemented (%s in %s) !",arrayName,alignObjFileName));
345 //______________________________________________________________________________
346 void AliAlignmentTracks::InitAlignObjs()
348 // Initialize the alignment objects array
349 Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
350 fAlignObjs = new AliAlignObj**[nLayers];
351 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
352 fAlignObjs[iLayer] = new AliAlignObj*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
353 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
354 UShort_t volid = AliAlignObj::LayerToVolUID(iLayer+ AliAlignObj::kFirstLayer,iModule);
355 fAlignObjs[iLayer][iModule] = new AliAlignObjAngles("",volid,0,0,0,0,0,0);
360 //______________________________________________________________________________
361 void AliAlignmentTracks::ResetAlignObjs()
363 // Reset the alignment objects array
364 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
365 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
366 fAlignObjs[iLayer][iModule]->SetPars(0,0,0,0,0,0);
370 //______________________________________________________________________________
371 void AliAlignmentTracks::DeleteAlignObjs()
373 // Delete the alignment objects array
374 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
375 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
376 if (fAlignObjs[iLayer][iModule])
377 delete fAlignObjs[iLayer][iModule];
378 delete [] fAlignObjs[iLayer];
380 delete [] fAlignObjs;
384 void AliAlignmentTracks::AlignDetector(AliAlignObj::ELayerID firstLayer,
385 AliAlignObj::ELayerID lastLayer,
386 AliAlignObj::ELayerID layerRangeMin,
387 AliAlignObj::ELayerID layerRangeMax,
390 // Align detector volumes within
391 // a given layer range
392 // (could be whole detector).
393 // Tracks are fitted only within
394 // the range defined by the user.
396 for (Int_t iLayer = firstLayer; iLayer < lastLayer; iLayer++)
397 nModules += AliAlignObj::LayerSize(iLayer);
398 TArrayI volIds(nModules);
401 for (Int_t iLayer = firstLayer; iLayer < lastLayer; iLayer++) {
402 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer); iModule++) {
403 UShort_t volId = AliAlignObj::LayerToVolUID(iLayer,iModule);
404 volIds.AddAt(volId,modnum);
409 while (iterations > 0) {
410 AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
415 //______________________________________________________________________________
416 void AliAlignmentTracks::AlignLayer(AliAlignObj::ELayerID layer,
417 AliAlignObj::ELayerID layerRangeMin,
418 AliAlignObj::ELayerID layerRangeMax,
421 // Align detector volumes within
423 // Tracks are fitted only within
424 // the range defined by the user.
425 Int_t nModules = AliAlignObj::LayerSize(layer);
426 TArrayI volIds(nModules);
427 for (Int_t iModule = 0; iModule < nModules; iModule++) {
428 UShort_t volId = AliAlignObj::LayerToVolUID(layer,iModule);
429 volIds.AddAt(volId,iModule);
432 while (iterations > 0) {
433 AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
438 //______________________________________________________________________________
439 void AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit,
442 // Align single detector volume to
444 // Tracks are fitted only within
445 // the second volume.
447 volIds.AddAt(volId,0);
448 TArrayI volIdsFit(1);
449 volIdsFit.AddAt(volIdFit,0);
451 while (iterations > 0) {
452 AlignVolumes(&volIds,&volIdsFit);
457 //______________________________________________________________________________
458 void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *volidsfit,
459 AliAlignObj::ELayerID layerRangeMin,
460 AliAlignObj::ELayerID layerRangeMax,
463 // Align a set of detector volumes.
464 // Tracks are fitted only within
465 // the range defined by the user
466 // (by layerRangeMin and layerRangeMax)
467 // or within the set of volidsfit
468 // Repeat the procedure 'iterations' times
470 Int_t nVolIds = volids->GetSize();
472 AliError("Volume IDs array is empty!");
476 // Load only the tracks with at least one
477 // space point in the set of volume (volids)
479 AliTrackPointArray **points;
480 // Start the iterations
481 while (iterations > 0) {
482 Int_t nArrays = LoadPoints(volids, points);
483 if (nArrays == 0) return;
485 AliTrackResiduals *minimizer = CreateMinimizer();
486 minimizer->SetNTracks(nArrays);
487 minimizer->InitAlignObj();
488 AliTrackFitter *fitter = CreateFitter();
489 for (Int_t iArray = 0; iArray < nArrays; iArray++) {
490 fitter->SetTrackPointArray(points[iArray], kFALSE);
491 fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax);
492 AliTrackPointArray *pVolId,*pTrack;
493 fitter->GetTrackResiduals(pVolId,pTrack);
494 minimizer->AddTrackPointArrays(pVolId,pTrack);
496 minimizer->Minimize();
498 // Update the alignment object(s)
499 for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
500 UShort_t volid = (*volids)[iVolId];
502 AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
503 AliAlignObj *alignObj = fAlignObjs[iLayer-AliAlignObj::kFirstLayer][iModule];
504 *alignObj *= *minimizer->GetAlignObj();
508 UnloadPoints(nArrays, points);
514 //______________________________________________________________________________
515 Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray** &points)
517 // Load track point arrays with at least
518 // one space point in a given set of detector
519 // volumes (array volids).
520 // Use the already created tree index for
524 AliError("Tree with the space point arrays not initialized!");
529 Int_t nVolIds = volids->GetSize();
531 AliError("Volume IDs array is empty!");
537 for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
538 UShort_t volid = (*volids)[iVolId];
540 AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
542 // In case of empty index
543 if (fLastIndex[iLayer-AliAlignObj::kFirstLayer][iModule] == 0) {
544 AliWarning(Form("There are no space-points belonging to the volume which is to be aligned (Volume ID =%d)!",volid));
547 nArrays += fLastIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
551 AliError("There are no space-points belonging to all of the volumes which are to be aligned!");
556 AliTrackPointArray* array = 0;
557 fPointsTree->SetBranchAddress("SP", &array);
559 // Allocate the pointer to the space-point arrays
560 points = new AliTrackPointArray*[nArrays];
561 for (Int_t i = 0; i < nArrays; i++) points[i] = 0x0;
563 // Init the array used to flag already loaded tree entries
564 Bool_t *indexUsed = new Bool_t[fPointsTree->GetEntries()];
565 for (Int_t i = 0; i < fPointsTree->GetEntries(); i++)
566 indexUsed[i] = kFALSE;
568 // Start the loop over the volume ids
570 for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
571 UShort_t volid = (*volids)[iVolId];
573 AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
575 Int_t nArraysId = fLastIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
576 TArrayI *index = fArrayIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
579 for (Int_t iArrayId = 0; iArrayId < nArraysId; iArrayId++) {
582 Int_t entry = (*index)[iArrayId];
583 if (indexUsed[entry] == kTRUE) continue;
584 fPointsTree->GetEvent(entry);
586 AliWarning("Wrong space point array index!");
589 indexUsed[entry] = kTRUE;
591 // Get the space-point array
592 Int_t nPoints = array->GetNPoints();
593 points[iArray] = new AliTrackPointArray(nPoints);
594 for (Int_t iPoint = 0; iPoint < nPoints; iPoint++) {
595 array->GetPoint(p,iPoint);
597 AliAlignObj::ELayerID layer = AliAlignObj::VolUIDToLayer(p.GetVolumeID(),modnum);
599 // Misalignment is introduced here
600 // Switch it off in case of real
603 AliAlignObj *misalignObj = fMisalignObjs[layer-AliAlignObj::kFirstLayer][modnum];
605 misalignObj->Transform(p);
607 // End of misalignment
609 AliAlignObj *alignObj = fAlignObjs[layer-AliAlignObj::kFirstLayer][modnum];
610 alignObj->Transform(p);
611 points[iArray]->AddPoint(iPoint,&p);
623 //______________________________________________________________________________
624 void AliAlignmentTracks::UnloadPoints(Int_t n, AliTrackPointArray **points)
626 // Unload track point arrays for a given
628 for (Int_t iArray = 0; iArray < n; iArray++)
629 delete points[iArray];
633 //______________________________________________________________________________
634 AliTrackFitter *AliAlignmentTracks::CreateFitter()
636 // Check if the user has already supplied
637 // a track fitter object.
638 // If not, create a default one.
640 fTrackFitter = new AliTrackFitterRieman;
645 //______________________________________________________________________________
646 AliTrackResiduals *AliAlignmentTracks::CreateMinimizer()
648 // Check if the user has already supplied
649 // a track residuals minimizer object.
650 // If not, create a default one.
652 fMinimizer = new AliTrackResidualsChi2;
657 //______________________________________________________________________________
658 Bool_t AliAlignmentTracks::Misalign(const char *misalignObjFileName, const char* arrayName)
660 // The method reads from a file a set of AliAlignObj which are
661 // then used to apply misalignments directly on the track
662 // space-points. The method is supposed to be used only for
663 // fast development and debugging of the alignment algorithms.
664 // Be careful not to use it in the case of 'real' alignment
665 // scenario since it will bias the results.
667 // Initialize the misalignment objects array
668 Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
669 fMisalignObjs = new AliAlignObj**[nLayers];
670 for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
671 fMisalignObjs[iLayer] = new AliAlignObj*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
672 for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
673 fMisalignObjs[iLayer][iModule] = 0x0;
676 // Open the misliagnment file and load the array with
677 // misalignment objects
678 TFile* inFile = TFile::Open(misalignObjFileName,"READ");
679 if (!inFile || !inFile->IsOpen()) {
680 AliError(Form("Could not open misalignment file %s !",misalignObjFileName));
684 TClonesArray* array = ((TClonesArray*) inFile->Get(arrayName));
686 AliError(Form("Could not find misalignment array %s in the file %s !",arrayName,misalignObjFileName));
692 // Store the misalignment objects for further usage
693 Int_t nObjs = array->GetEntriesFast();
694 AliAlignObj::ELayerID layerId; // volume layer
695 Int_t modId; // volume ID inside the layer
696 for(Int_t i=0; i<nObjs; i++)
698 AliAlignObj* alObj = (AliAlignObj*)array->UncheckedAt(i);
699 alObj->GetVolUID(layerId,modId);
700 fMisalignObjs[layerId-AliAlignObj::kFirstLayer][modId] = alObj;