]> git.uio.no Git - u/mrichter/AliRoot.git/blob - STEER/AliAlignmentTracks.cxx
Bug fix
[u/mrichter/AliRoot.git] / STEER / AliAlignmentTracks.cxx
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>
26 #include <TFile.h>
27
28 #include "AliAlignmentTracks.h"
29 #include "AliTrackPointArray.h"
30 #include "AliAlignObjAngles.h"
31 #include "AliTrackFitterRieman.h"
32 #include "AliTrackResidualsChi2.h"
33 #include "AliESD.h"
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),
47   fAlignObjs(0),
48   fMisalignObjs(0),
49   fTrackFitter(0),
50   fMinimizer(0),
51   fDoUpdate(kTRUE)
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),
67   fAlignObjs(0),
68   fMisalignObjs(0),
69   fTrackFitter(0),
70   fMinimizer(0),
71   fDoUpdate(kTRUE)
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):
83   fESDChain(new TChain(esdtreename)),
84   fPointsFilename("AliTrackPoints.root"),
85   fPointsFile(0),
86   fPointsTree(0),
87   fLastIndex(0),
88   fArrayIndex(0),
89   fIsIndexBuilt(kFALSE),
90   fAlignObjs(0),
91   fMisalignObjs(0),
92   fTrackFitter(0),
93   fMinimizer(0),
94   fDoUpdate(kTRUE)
95 {
96   // Constructor in the case
97   // the user provides a single ESD file
98   // or a directory containing ESD files
99   fESDChain->Add(esdfilename);
100
101   InitIndex();
102   InitAlignObjs();
103 }
104
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
145 //______________________________________________________________________________
146 void AliAlignmentTracks::ProcessESD()
147 {
148   // Analyzes and filters ESD tracks
149   // Stores the selected track space points
150   // into the output file
151
152   if (!fESDChain) return;
153
154   AliESD *esd = 0;
155   fESDChain->SetBranchAddress("ESD",&esd);
156   AliESDfriend *esdf = 0; 
157   fESDChain->SetBranchStatus("ESDfriend*",1);
158   fESDChain->SetBranchAddress("ESDfriend.",&esdf);
159
160   // Open the output file
161   if (fPointsFilename.Data() == "") {
162     AliWarning("Incorrect output filename!");
163     return;
164   }
165
166   TFile *pointsFile = TFile::Open(fPointsFilename,"RECREATE");
167   if (!pointsFile || !pointsFile->IsOpen()) {
168     AliWarning(Form("Can't open %s !",fPointsFilename.Data()));
169     return;
170   }
171
172   TTree *pointsTree = new TTree("spTree", "Tree with track space point arrays");
173   const AliTrackPointArray *array = 0;
174   pointsTree->Branch("SP","AliTrackPointArray", &array);
175   Int_t ievent = 0;
176   while (fESDChain->GetEntry(ievent++)) {
177     if (!esd) break;
178
179     esd->SetESDfriend(esdf); //Attach the friend to the ESD
180
181     Int_t ntracks = esd->GetNumberOfTracks();
182     for (Int_t itrack=0; itrack < ntracks; itrack++) {
183       AliESDtrack * track = esd->GetTrack(itrack);
184       if (!track) continue;
185  
186      //  UInt_t status = AliESDtrack::kITSpid; 
187 //       status|=AliESDtrack::kTPCpid; 
188 //       status|=AliESDtrack::kTRDpid; 
189 //       if ((track->GetStatus() & status) != status) continue;
190
191       if (track->GetP() < 0.5) continue;
192
193       array = track->GetTrackPointArray();
194       pointsTree->Fill();
195     }
196   }
197
198   if (!pointsTree->Write()) {
199     AliWarning("Can't write the tree with track point arrays!");
200     return;
201   }
202
203   pointsFile->Close();
204 }
205
206 //______________________________________________________________________________
207 void AliAlignmentTracks::ProcessESD(TSelector *selector)
208 {
209   AliWarning(Form("ESD processing based on selector is not yet implemented (%p) !",selector));
210 }
211
212 //______________________________________________________________________________
213 void AliAlignmentTracks::BuildIndex()
214 {
215   // Build index of points tree entries
216   // Used for access based on the volume IDs
217   if (fIsIndexBuilt) return;
218
219   fIsIndexBuilt = kTRUE;
220
221   // Dummy object is created in order
222   // to initialize the volume paths
223   AliAlignObjAngles alobj;
224
225   fPointsFile = TFile::Open(fPointsFilename);
226   if (!fPointsFile || !fPointsFile->IsOpen()) {
227     AliWarning(Form("Can't open %s !",fPointsFilename.Data()));
228     return;
229   }
230   
231   //  AliTrackPointArray* array = new AliTrackPointArray;
232   AliTrackPointArray* array = 0;
233   fPointsTree = (TTree*) fPointsFile->Get("spTree");
234   if (!fPointsTree) {
235     AliWarning("No pointsTree found!");
236     return;
237   }
238   fPointsTree->SetBranchAddress("SP", &array);
239
240   Int_t nArrays = (Int_t)fPointsTree->GetEntries();
241   for (Int_t iArray = 0; iArray < nArrays; iArray++)
242     {
243       fPointsTree->GetEvent(iArray);
244       if (!array) continue;
245       for (Int_t ipoint = 0; ipoint < array->GetNPoints(); ipoint++) {
246         UShort_t volId = array->GetVolumeID()[ipoint];
247         // check if the volId is valid
248         if (!AliAlignObj::SymName(volId)) {
249           AliError(Form("The volume id %d has no default volume name !",
250                         volId));
251           continue;
252         }
253         Int_t modId;
254         Int_t layerId = AliAlignObj::VolUIDToLayer(volId,modId)
255                       - AliAlignObj::kFirstLayer;
256         if (!fArrayIndex[layerId][modId]) {
257           //first entry for this volume
258           fArrayIndex[layerId][modId] = new TArrayI(1000);
259         }
260         else {
261           Int_t size = fArrayIndex[layerId][modId]->GetSize();
262           // If needed allocate new size
263           if (fLastIndex[layerId][modId] >= size)
264             fArrayIndex[layerId][modId]->Set(size + 1000);
265         }
266
267         // Check if the index is already filled
268         Bool_t fillIndex = kTRUE;
269         if (fLastIndex[layerId][modId] != 0) {
270           if ((*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]-1] == iArray)
271             fillIndex = kFALSE;
272         }
273         // Fill the index array and store last filled index
274         if (fillIndex) {
275           (*fArrayIndex[layerId][modId])[fLastIndex[layerId][modId]] = iArray;
276           fLastIndex[layerId][modId]++;
277         }
278       }
279     }
280
281 }
282
283 //______________________________________________________________________________
284 void AliAlignmentTracks::InitIndex()
285 {
286   // Initialize the index arrays
287   Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
288   fLastIndex = new Int_t*[nLayers];
289   fArrayIndex = new TArrayI**[nLayers];
290   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
291     fLastIndex[iLayer] = new Int_t[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
292     fArrayIndex[iLayer] = new TArrayI*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
293     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
294       fLastIndex[iLayer][iModule] = 0;
295       fArrayIndex[iLayer][iModule] = 0;
296     }
297   }
298 }
299
300 //______________________________________________________________________________
301 void AliAlignmentTracks::ResetIndex()
302 {
303   // Reset the value of the last filled index
304   // Do not realocate memory
305
306   fIsIndexBuilt = kFALSE;
307   
308   for (Int_t iLayer = 0; iLayer < AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer; iLayer++) {
309     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
310       fLastIndex[iLayer][iModule] = 0;
311     }
312   }
313 }
314
315 //______________________________________________________________________________
316 void AliAlignmentTracks::DeleteIndex()
317 {
318   // Delete the index arrays
319   // Called by the destructor
320   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
321     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
322       if (fArrayIndex[iLayer][iModule]) {
323         delete fArrayIndex[iLayer][iModule];
324         fArrayIndex[iLayer][iModule] = 0;
325       }
326     }
327     delete [] fLastIndex[iLayer];
328     delete [] fArrayIndex[iLayer];
329   }
330   delete [] fLastIndex;
331   delete [] fArrayIndex;
332 }
333
334 //______________________________________________________________________________
335 Bool_t AliAlignmentTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName)
336 {
337   // Read alignment object from a file
338   // To be replaced by a call to CDB
339   AliWarning(Form("Method not yet implemented (%s in %s) !",arrayName,alignObjFileName));
340
341   return kFALSE;
342 }
343
344 //______________________________________________________________________________
345 void AliAlignmentTracks::InitAlignObjs()
346 {
347   // Initialize the alignment objects array
348   Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
349   fAlignObjs = new AliAlignObj**[nLayers];
350   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
351     fAlignObjs[iLayer] = new AliAlignObj*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
352     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++) {
353       UShort_t volid = AliAlignObj::LayerToVolUID(iLayer+ AliAlignObj::kFirstLayer,iModule);
354       fAlignObjs[iLayer][iModule] = new AliAlignObjAngles(AliAlignObj::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
355     }
356   }
357 }
358
359 //______________________________________________________________________________
360 void AliAlignmentTracks::ResetAlignObjs()
361 {
362   // Reset the alignment objects array
363   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
364     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
365       fAlignObjs[iLayer][iModule]->SetPars(0,0,0,0,0,0);
366   }
367 }
368
369 //______________________________________________________________________________
370 void AliAlignmentTracks::DeleteAlignObjs()
371 {
372   // Delete the alignment objects array
373   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
374     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
375       if (fAlignObjs[iLayer][iModule])
376         delete fAlignObjs[iLayer][iModule];
377     delete [] fAlignObjs[iLayer];
378   }
379   delete [] fAlignObjs;
380   fAlignObjs = 0;
381 }
382
383 void AliAlignmentTracks::AlignDetector(AliAlignObj::ELayerID firstLayer,
384                                        AliAlignObj::ELayerID lastLayer,
385                                        AliAlignObj::ELayerID layerRangeMin,
386                                        AliAlignObj::ELayerID layerRangeMax,
387                                        Int_t iterations)
388 {
389   // Align detector volumes within
390   // a given layer range
391   // (could be whole detector).
392   // Tracks are fitted only within
393   // the range defined by the user.
394   Int_t nModules = 0;
395   for (Int_t iLayer = firstLayer; iLayer < lastLayer; iLayer++)
396     nModules += AliAlignObj::LayerSize(iLayer);
397   TArrayI volIds(nModules);
398
399   Int_t modnum = 0;
400   for (Int_t iLayer = firstLayer; iLayer < lastLayer; iLayer++) {
401     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer); iModule++) {
402       UShort_t volId = AliAlignObj::LayerToVolUID(iLayer,iModule);
403       volIds.AddAt(volId,modnum);
404       modnum++;
405     }
406   }
407
408   while (iterations > 0) {
409     AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
410     iterations--;
411   }
412 }
413
414 //______________________________________________________________________________
415 void AliAlignmentTracks::AlignLayer(AliAlignObj::ELayerID layer,
416                                     AliAlignObj::ELayerID layerRangeMin,
417                                     AliAlignObj::ELayerID layerRangeMax,
418                                     Int_t iterations)
419 {
420   // Align detector volumes within
421   // a given layer.
422   // Tracks are fitted only within
423   // the range defined by the user.
424   Int_t nModules = AliAlignObj::LayerSize(layer);
425   TArrayI volIds(nModules);
426   for (Int_t iModule = 0; iModule < nModules; iModule++) {
427     UShort_t volId = AliAlignObj::LayerToVolUID(layer,iModule);
428     volIds.AddAt(volId,iModule);
429   }
430
431   while (iterations > 0) {
432     AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
433     iterations--;
434   }
435 }
436
437 //______________________________________________________________________________
438 void AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit,
439                                      Int_t iterations)
440 {
441   // Align single detector volume to
442   // another volume.
443   // Tracks are fitted only within
444   // the second volume.
445   TArrayI volIds(1);
446   volIds.AddAt(volId,0);
447   TArrayI volIdsFit(1);
448   volIdsFit.AddAt(volIdFit,0);
449
450   while (iterations > 0) {
451     AlignVolumes(&volIds,&volIdsFit);
452     iterations--;
453   }
454 }
455
456 //______________________________________________________________________________
457 void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *volidsfit,
458                                      AliAlignObj::ELayerID layerRangeMin,
459                                      AliAlignObj::ELayerID layerRangeMax,
460                                      Int_t iterations)
461 {
462   // Align a set of detector volumes.
463   // Tracks are fitted only within
464   // the range defined by the user
465   // (by layerRangeMin and layerRangeMax)
466   // or within the set of volidsfit
467   // Repeat the procedure 'iterations' times
468
469   Int_t nVolIds = volids->GetSize();
470   if (nVolIds == 0) {
471     AliError("Volume IDs array is empty!");
472     return;
473   }
474
475   // Load only the tracks with at least one
476   // space point in the set of volume (volids)
477   BuildIndex();
478   AliTrackPointArray **points;
479   // Start the iterations
480   while (iterations > 0) {
481     Int_t nArrays = LoadPoints(volids, points);
482     if (nArrays == 0) return;
483
484     AliTrackResiduals *minimizer = CreateMinimizer();
485     minimizer->SetNTracks(nArrays);
486     minimizer->InitAlignObj();
487     AliTrackFitter *fitter = CreateFitter();
488     for (Int_t iArray = 0; iArray < nArrays; iArray++) {
489       if (!points[iArray]) continue;
490       fitter->SetTrackPointArray(points[iArray], kFALSE);
491       if (fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE) continue;
492       AliTrackPointArray *pVolId,*pTrack;
493       fitter->GetTrackResiduals(pVolId,pTrack);
494       minimizer->AddTrackPointArrays(pVolId,pTrack);
495     }
496     minimizer->Minimize();
497
498     // Update the alignment object(s)
499     if (fDoUpdate) for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
500       UShort_t volid = (*volids)[iVolId];
501       Int_t iModule;
502       AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
503       AliAlignObj *alignObj = fAlignObjs[iLayer-AliAlignObj::kFirstLayer][iModule];      
504       *alignObj *= *minimizer->GetAlignObj();
505       alignObj->Print("");
506     }
507
508     UnloadPoints(nArrays, points);
509     
510     iterations--;
511   }
512 }
513   
514 //______________________________________________________________________________
515 Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray** &points)
516 {
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
521   // fast access.
522
523   if (!fPointsTree) {
524     AliError("Tree with the space point arrays not initialized!");
525     points = 0;
526     return 0;
527   }
528
529   Int_t nVolIds = volids->GetSize();
530   if (nVolIds == 0) {
531     AliError("Volume IDs array is empty!");
532     points = 0;
533     return 0;
534   }
535
536   Int_t nArrays = 0;
537   for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
538     UShort_t volid = (*volids)[iVolId];
539     Int_t iModule;
540     AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
541
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));
545       continue;
546     }
547     nArrays += fLastIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
548   }
549
550   if (nArrays == 0) {
551     AliError("There are no space-points belonging to all of the volumes which are to be aligned!");
552     points = 0;
553     return 0;
554   }
555
556   AliTrackPointArray* array = 0;
557   fPointsTree->SetBranchAddress("SP", &array);
558
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;
562
563   // Init the array used to flag already loaded tree entries
564   Bool_t *indexUsed = new Bool_t[(UInt_t)fPointsTree->GetEntries()];
565   for (Int_t i = 0; i < fPointsTree->GetEntries(); i++)
566     indexUsed[i] = kFALSE;
567
568   // Start the loop over the volume ids
569   Int_t iArray = 0;
570   for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
571     UShort_t volid = (*volids)[iVolId];
572     Int_t iModule;
573     AliAlignObj::ELayerID iLayer = AliAlignObj::VolUIDToLayer(volid,iModule);
574
575     Int_t nArraysId = fLastIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
576     TArrayI *index = fArrayIndex[iLayer-AliAlignObj::kFirstLayer][iModule];
577     AliTrackPoint p;
578
579     for (Int_t iArrayId = 0; iArrayId < nArraysId; iArrayId++) {
580
581       // Get tree entry
582       Int_t entry = (*index)[iArrayId];
583       if (indexUsed[entry] == kTRUE) continue;
584       fPointsTree->GetEvent(entry);
585       if (!array) {
586         AliWarning("Wrong space point array index!");
587         continue;
588       }
589       indexUsed[entry] = kTRUE;
590
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);
596         Int_t modnum;
597         AliAlignObj::ELayerID layer = AliAlignObj::VolUIDToLayer(p.GetVolumeID(),modnum);
598         // check if the layer id is valid
599         if ((layer < AliAlignObj::kFirstLayer) ||
600             (layer >= AliAlignObj::kLastLayer)) {
601           AliError(Form("Layer index is invalid: %d (%d -> %d) !",
602                         layer,AliAlignObj::kFirstLayer,AliAlignObj::kLastLayer-1));
603           continue;
604         }
605         if ((modnum >= AliAlignObj::LayerSize(layer)) ||
606             (modnum < 0)) {
607           AliError(Form("Module number inside layer %d is invalid: %d (0 -> %d)",
608                         layer,modnum,AliAlignObj::LayerSize(layer)));
609           continue;
610         }
611
612         // Misalignment is introduced here
613         // Switch it off in case of real
614         // alignment job!
615         if (fMisalignObjs) {
616           AliAlignObj *misalignObj = fMisalignObjs[layer-AliAlignObj::kFirstLayer][modnum];
617           if (misalignObj)
618             misalignObj->Transform(p);
619         }
620         // End of misalignment
621
622         AliAlignObj *alignObj = fAlignObjs[layer-AliAlignObj::kFirstLayer][modnum];
623         alignObj->Transform(p);
624         points[iArray]->AddPoint(iPoint,&p);
625       }
626       iArray++;
627     }
628   }
629
630
631   delete [] indexUsed;
632
633   return nArrays;
634 }
635
636 //______________________________________________________________________________
637 void AliAlignmentTracks::UnloadPoints(Int_t n, AliTrackPointArray **points)
638 {
639   // Unload track point arrays for a given
640   // detector volume
641   for (Int_t iArray = 0; iArray < n; iArray++)
642     delete points[iArray];
643   delete [] points;
644 }
645
646 //______________________________________________________________________________
647 AliTrackFitter *AliAlignmentTracks::CreateFitter()
648 {
649   // Check if the user has already supplied
650   // a track fitter object.
651   // If not, create a default one.
652   if (!fTrackFitter)
653     fTrackFitter = new AliTrackFitterRieman;
654
655   return fTrackFitter;
656 }
657
658 //______________________________________________________________________________
659 AliTrackResiduals *AliAlignmentTracks::CreateMinimizer()
660 {
661   // Check if the user has already supplied
662   // a track residuals minimizer object.
663   // If not, create a default one.
664   if (!fMinimizer)
665     fMinimizer = new AliTrackResidualsChi2;
666
667   return fMinimizer;
668 }
669
670 //______________________________________________________________________________
671 Bool_t AliAlignmentTracks::Misalign(const char *misalignObjFileName, const char* arrayName)
672 {
673   // The method reads from a file a set of AliAlignObj which are
674   // then used to apply misalignments directly on the track
675   // space-points. The method is supposed to be used only for
676   // fast development and debugging of the alignment algorithms.
677   // Be careful not to use it in the case of 'real' alignment
678   // scenario since it will bias the results.
679
680   // Initialize the misalignment objects array
681   Int_t nLayers = AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer;
682   fMisalignObjs = new AliAlignObj**[nLayers];
683   for (Int_t iLayer = 0; iLayer < (AliAlignObj::kLastLayer - AliAlignObj::kFirstLayer); iLayer++) {
684     fMisalignObjs[iLayer] = new AliAlignObj*[AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer)];
685     for (Int_t iModule = 0; iModule < AliAlignObj::LayerSize(iLayer + AliAlignObj::kFirstLayer); iModule++)
686       fMisalignObjs[iLayer][iModule] = 0x0;
687   }
688
689   // Open the misliagnment file and load the array with
690   // misalignment objects
691   TFile* inFile = TFile::Open(misalignObjFileName,"READ");
692   if (!inFile || !inFile->IsOpen()) {
693     AliError(Form("Could not open misalignment file %s !",misalignObjFileName));
694     return kFALSE;
695   }
696
697   TClonesArray* array = ((TClonesArray*) inFile->Get(arrayName));
698   if (!array) {
699     AliError(Form("Could not find misalignment array %s in the file %s !",arrayName,misalignObjFileName));
700     inFile->Close();
701     return kFALSE;
702   }
703   inFile->Close();
704
705   // Store the misalignment objects for further usage  
706   Int_t nObjs = array->GetEntriesFast();
707   AliAlignObj::ELayerID layerId; // volume layer
708   Int_t modId; // volume ID inside the layer
709   for(Int_t i=0; i<nObjs; i++)
710     {
711       AliAlignObj* alObj = (AliAlignObj*)array->UncheckedAt(i);
712       alObj->GetVolUID(layerId,modId);
713       fMisalignObjs[layerId-AliAlignObj::kFirstLayer][modId] = alObj;
714     }
715   return kTRUE;
716 }