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 **************************************************************************/
19 #include <TStopwatch.h>
22 #include <TClonesArray.h>
24 #include <TGeoManager.h>
26 #include <TGeoMatrix.h>
30 #include "AliITSRealignTracks.h"
31 #include "AliAlignmentTracks.h"
32 #include "AliAlignObjParams.h"
33 #include "AliAlignObj.h"
34 #include "AliGeomManager.h"
35 #include "AliTrackFitter.h"
36 #include "AliTrackFitterKalman.h"
37 #include "AliTrackFitterRieman.h"
38 #include "AliTrackResidualsFast.h"
39 #include "AliTrackResidualsChi2.h"
40 #include "AliTrackResidualsLinear.h"
46 ClassImp(AliITSRealignTracks)
48 const Int_t referSect=2;
51 AliITSRealignTracks::AliITSRealignTracks(TString minimizer,Int_t fit,Bool_t covUsed,TString fileintro,TString geometryfile,TString misalignmentFile,TString startingfile):
80 // minimizer="fast"->AliTrackResidualFast minimizer
81 // "minuit"->AliTrackResidualChi2 minimizer
82 // "minuitnorot"->AliTrackResidualChi2 minimizer without rotations degrees of freedom
83 // "linear"->AliTrackResidualLinear minimizer
84 //fit=0-> Riemann Fitter, fit=1->Kalman
85 //fileintro=file into which the Tree with the space points is stored
86 //geometryfile=file containing the geometry
89 SetPointsFilename(fileintro.Data());
90 SetGeomFilename(geometryfile);
92 if(!InitSurveyObjs(kFALSE))AliWarning("Unable to set Survey AlignObjs!");
94 if(startingfile=="")printf("Starting from default geometry \n");
95 else ReadAlignObjs(startingfile.Data(),"ITSAlignObjs");
97 if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT INTRODUCED \n");
99 Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
100 if(!misal)AliWarning("Incorrect fake misalignment filename!");;
103 if(!gGeoManager) AliGeomManager::LoadGeometry(fgeomfilename.Data());
104 if(covUsed)SetCovIsUsed(kTRUE);
105 if(!SelectFitter(fit))AliWarning("Incorrect fitter assignment!");
106 if(!SelectMinimizer(minimizer))AliWarning("Incorrect minimizer assignment!");
113 AliITSRealignTracks::AliITSRealignTracks(const AliITSRealignTracks &realignTracks):
114 AliAlignmentTracks(),
115 fSurveyObjs(new AliAlignObj**(*realignTracks.fSurveyObjs)),
116 fgeomfilename(realignTracks.fgeomfilename),
117 fmintracks(realignTracks.fmintracks),
118 fCovIsUsed(realignTracks.fCovIsUsed),
119 fUpdateCov(realignTracks.fUpdateCov),
120 fVarySigmaY(realignTracks.fVarySigmaY),
121 fCorrModules(new Double_t *(*realignTracks.fCorrModules)),
122 fLimitCorr(realignTracks.fLimitCorr),
123 fsigmaY(realignTracks.fsigmaY),
125 fAlignDrawObjs(realignTracks.fAlignDrawObjs),
126 fCanvPar(realignTracks.fCanvPar),
127 fCanvGr(realignTracks.fCanvGr),
128 fgrIterMeanX(realignTracks.fgrIterMeanX),
129 fgrIterRMSX(realignTracks.fgrIterRMSX),
130 fgrIterMeanY(realignTracks.fgrIterMeanY),
131 fgrIterRMSY(realignTracks.fgrIterRMSY),
132 fgrIterMeanZ(realignTracks.fgrIterMeanZ),
133 fgrIterRMSZ(realignTracks.fgrIterRMSZ),
134 fgrIterMeanPsi(realignTracks.fgrIterMeanPsi),
135 fgrIterRMSPsi(realignTracks.fgrIterRMSPsi),
136 fgrIterMeanTheta(realignTracks.fgrIterMeanTheta),
137 fgrIterRMSTheta(realignTracks.fgrIterRMSTheta),
138 fgrIterMeanPhi(realignTracks.fgrIterMeanPhi),
139 fgrIterRMSPhi(realignTracks.fgrIterRMSPhi)
142 AliWarning("Can't copy AliAlignmentTracks Data member!");
145 AliITSRealignTracks& AliITSRealignTracks::operator=(const AliITSRealignTracks &obj){
146 ////////////////////////
147 // Assignment operator
148 ////////////////////////
149 this->~AliITSRealignTracks();
150 new(this) AliITSRealignTracks(obj);
154 AliITSRealignTracks::~AliITSRealignTracks(){
157 if(fSurveyObjs) DeleteSurveyObjs();
158 if(fAlignDrawObjs) DeleteDrawHists();
159 //delete [] fSurveyObjs;
164 //_____________________________
165 Bool_t AliITSRealignTracks::SelectFitter(Int_t fit,Int_t minTrackPoint){
168 AliTrackFitterKalman *fitter= new AliTrackFitterKalman();
169 fitter->SetMinNPoints(minTrackPoint);
170 SetTrackFitter(fitter);
174 AliTrackFitterRieman *fitter=new AliTrackFitterRieman();
175 fitter=new AliTrackFitterRieman();
176 fitter->SetMinNPoints(minTrackPoint);
177 SetTrackFitter(fitter);
185 Bool_t AliITSRealignTracks::SelectMinimizer(TString minimizer,Int_t minpoints,const Bool_t *coord){
186 AliTrackResiduals *res;
187 if(minimizer=="minuit"){
188 res = new AliTrackResidualsChi2();
190 for(Int_t j=0;j<6;j++){
191 if(coord[j])res->FixParameter(j);
195 else if(minimizer=="minuitnorot"){
196 res = new AliTrackResidualsChi2();
197 res->FixParameter(3);
198 res->FixParameter(4);
199 res->FixParameter(5);
202 else if(minimizer=="fast"){
203 res = new AliTrackResidualsFast();
205 for(Int_t j=0;j<6;j++){
206 if(coord[j])res->FixParameter(j);
210 else if(minimizer=="linear"){
211 res = new AliTrackResidualsLinear();
215 printf("Trying to set a non existing minimizer! \n");
219 res->SetMinNPoints(minpoints);
225 //____________________________________
226 void AliITSRealignTracks::SetVarySigmaY(Bool_t varysigmay,Double_t sigmaYfixed){
228 fVarySigmaY=varysigmay;
230 if(sigmaYfixed>0.)fsigmaY=sigmaYfixed;
232 printf("Negative assignment to sigmaY! set it to default value (1 cm) \n");
238 //_______________________________________
239 void AliITSRealignTracks::RealignITSVolIndependent(Int_t iter1,Int_t iterations,Int_t minNtracks,Int_t layer,Int_t minTrackPoint){
242 //iter1=#iterations inside AliAlignmentTracks::AliAlignVolumesITS method
243 //iterations=#iterations on the all procedure
244 //layer=0->all ITS, otherways the usual notation is considered (1=SPD1,2=SPD2,3=SDD1,4=SDD2,5=SSD1,6=SSD2)
245 //minNtracks=minimun number of tracks passing through a module in order to try to realign the module itsself
246 // if minNtracks<0, minimun number of tracks is |minNtracks|*minNumPoint[layer]/fact (see the code below): this allows a different
247 // choice of the number of tracks required on different layers and to vary these numbers once tuned the relative proportions.
248 //minTrackPoint=minimun number of "good" points required to a track (THE POINT ON THE MODULE THAT IS GOING TO BE REALIGNED
249 //IS NEVER CONSIDERED->max number can be required is 11 for cosmics tracks) for the track being considered in the minimization
252 fTrackFitter->SetMinNPoints(minTrackPoint);
253 TArrayI volIDs2(2200);
259 Int_t iLayer,iLayerToAlign;
261 Int_t minNumPoint[6]={100,100,100,100,50,50};
266 Int_t layerNum,modNum,lastVolid=0;
267 TNtuple *ntVolumeAlign=new TNtuple("ntVolumeAlign","NTuple with volume tried to be realigned","layerNum:modNum:volumeIDnum");
269 TStopwatch *timer=new TStopwatch();
275 for(Int_t iter=0;iter<iterations;iter++){
277 //Starting Independent Modules Realignment
278 for(iLayerToAlign=(Int_t)AliGeomManager::kSPD1;iLayerToAlign<=(Int_t)AliGeomManager::kSSD2;iLayerToAlign++){
279 if(layer!=0&&iLayerToAlign!=layer)continue;
282 for(Int_t k=(Int_t)AliGeomManager::kSPD1;k<=(Int_t)AliGeomManager::kSSD2;k++){
283 size+=AliGeomManager::LayerSize(k);
284 printf("size: %d \n",size);
287 for (Int_t iModule=0;iModule<AliGeomManager::LayerSize(iLayerToAlign);iModule++){
290 if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNumPoint[iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer]*(-1*minNtracks/fact))continue; }
291 else if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNtracks)continue;
293 UShort_t volidAl = AliGeomManager::LayerToVolUID(iLayerToAlign,iModule);
296 volIDsFit.Set(size-1);
297 for (iLayer=AliGeomManager::kSPD1;iLayer<AliGeomManager::kTPC1;iLayer++){
298 for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
299 volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
300 if(AliGeomManager::LayerToVolUID(iLayer,iModule2)==volidAl)continue;
301 volIDsFit.AddAt(volid,j);
305 volIDs.AddAt((Int_t)volidAl,0);
306 if(iter==iterations-1){
307 volIDs2.AddAt(volidAl,lastVolid);
310 volIDs2.AddAt(volidAl,lastVolid);
311 AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iter1);
316 if((iter+1)%5==0||iter==0||iter==1||iter==2||iter==3||iter==iterations-1){
317 command="RealignObj";
319 command.Append(".root");
320 WriteRealignObjArray(command.Data(),AliGeomManager::kSPD1,AliGeomManager::kSSD2);
325 if(j==0){printf("j=0 \n");return;}
326 for(Int_t k=0;k<volIDs2.GetSize();k++){
327 if(volIDs2.At(k)==0)break;
328 layerNum=AliGeomManager::VolUIDToLayer(volIDs2.At(k),modNum);
329 ntVolumeAlign->Fill(layerNum,modNum,volIDs2.At(k));
331 printf("End of selecting modules cycle: %d modules selected \n",j);
332 TFile *f=new TFile("RealignVolNt.root","RECREATE");
334 ntVolumeAlign->Write();
343 void AliITSRealignTracks::RealignITStracks(TString minimizer,Int_t fit=0,Int_t iter1=1,Int_t iterations=5,Int_t minNtracks=-10,Int_t layer=0,Int_t minTrackPoint=6,Bool_t covUsed=kFALSE,TString misalignmentFile="",TString startingfile="",Int_t doGlobal=1)
346 // minimizer="fast"->AliTrackResidualFast minimizer
347 // "minuit"->AliTrackResidualChi2 minimizer
348 // "minuitnorot"->AliTrackResidualChi2 minimizer without rotations degrees of freedom
349 // "linear"->AliTrackResidualLinear minimizer
350 //fit=0-> Riemann Fitter, fit=1->Kalman
351 //iter1=#iterations inside AliAlignmentTracks::AliAlignVolumesITS method
352 //iterations=#iterations on the all procedure
353 //layer=0->all ITS, otherways the usual notation is considered (1=SPD1,2=SPD2,3=SDD1,4=SDD2,5=SSD1,6=SSD2)
354 //minNtracks=minimun number of tracks passing through a module in order to try to realign the module itsself
355 // if minNtracks<0, minimun number of tracks is |minNtracks|*minNumPoint[layer]/fact (see the code below): this allows a different
356 // choice of the number of tracks required on different layers and to vary these numbers once tuned the relative proportions.
357 //minTrackPoint=minimun number of "good" points required to a track (THE POINT ON THE MODULE THAT IS GOING TO BE REALIGNED
358 // IS NEVER CONSIDERED->max number that can be required is 11 for cosmics tracks) for the track being considered in the minimization
359 //doGlobal : do global realignment, 0=no, 1= yes, 2=only global
362 TArrayI volIDs2(2200);
368 Int_t iLayer,iLayerToAlign;
370 Int_t minNumPoint[6]={100,100,100,100,50,50};
375 Int_t layerNum,modNum,lastVolid=0;
376 TNtuple *ntVolumeAlign=new TNtuple("ntVolumeAlign","NTuple with volume tried to be realigned","layerNum:modNum:volumeIDnum");
379 if(!SelectFitter(fit))AliWarning("Incorrect fitter assignment!");
380 if(!SelectMinimizer(minimizer))AliWarning("Incorrect minimizer assignment!");
381 if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT INTRODUCED \n");
383 Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
388 TStopwatch *timer=new TStopwatch();
394 if(startingfile=="")printf("Starting from default geometry \n");
396 printf("Starting from AlignObjs file: %s",startingfile.Data());
397 ReadAlignObjs(startingfile.Data(),"ITSAlignObjs");
400 for(Int_t iter=0;iter<iterations;iter++){
401 if(covUsed)SetCovIsUsed(kTRUE);
404 //START HIERARCHY REALIGNMENT
406 if(layer==0&&(doGlobal==1||doGlobal==2)){
407 for(Int_t siter=0;siter<5;siter++){
408 fTrackFitter->SetMinNPoints(2);
409 SetCovUpdate(kFALSE);
410 AlignSPDHalfBarrelToSectorRef(referSect,3);
411 // AlignSPDBarrel(1);
412 // if(siter==0)SetCovUpdate(kFALSE);
413 // AlignSPDHalfBarrel(0,3);
414 // SetCovUpdate(kTRUE);
415 AlignSPDHalfBarrelToHalfBarrel(1,3);
416 // AlignSPDHalfBarrelToSectorRef(referSect,3);
417 for(Int_t sector=0;sector<10;sector++){
419 if(sector==referSect)continue;
420 AlignSPDSectorWithSectors(sector,1);
424 for(Int_t lay=1;lay<=6;lay++){
425 if(!AlignLayerToSector(lay,referSect,3))AlignLayerToSPDHalfBarrel(lay,0,3);
427 AlignSPDHalfBarrel(0,3);
429 Int_t layers[6]={2,2,1,0,0,0};
430 fTrackFitter->SetMinNPoints(4);
431 AlignLayersToLayers(layers,1);
433 fTrackFitter->SetMinNPoints(6);
435 layers[3]=1;//{2,2,2,1,0,0};
436 AlignLayersToLayers(layers,1);
438 fTrackFitter->SetMinNPoints(6);
440 layers[4]=1;//{2,2,2,2,1,0};
441 AlignLayersToLayers(layers,1);
443 fTrackFitter->SetMinNPoints(6);
445 layers[5]=1;//{2,2,2,2,2,1};
446 AlignLayersToLayers(layers,1);
449 for(Int_t sector=0;sector<10;sector++){
450 AlignSPDSectorToOuterLayers(sector,1);
452 WriteRealignObjArray("AfterGlobal.root",AliGeomManager::kSPD1,AliGeomManager::kSSD2);
456 if(doGlobal==2)return;
458 if(covUsed)SetCovUpdate(kTRUE);
462 // STARTS INDEPENDENT MOULES REALIGNMENT
464 fTrackFitter->SetMinNPoints(minTrackPoint);
465 for(iLayerToAlign=(Int_t)AliGeomManager::kSPD1;iLayerToAlign<=(Int_t)AliGeomManager::kSSD2;iLayerToAlign++){
466 if(layer!=0&&iLayerToAlign!=layer)continue;
469 for(Int_t k=(Int_t)AliGeomManager::kSPD1;k<=(Int_t)AliGeomManager::kSSD2;k++){
470 size+=AliGeomManager::LayerSize(k);
471 printf("size: %d \n",size);
474 for (Int_t iModule=0;iModule<AliGeomManager::LayerSize(iLayerToAlign);iModule++){
477 if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNumPoint[iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer]*(-1*minNtracks/fact))continue;
479 else if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNtracks)continue;
481 UShort_t volidAl = AliGeomManager::LayerToVolUID(iLayerToAlign,iModule);
484 volIDsFit.Set(size-1);
485 for (iLayer=(Int_t)AliGeomManager::kSPD1;iLayer<(Int_t)AliGeomManager::kTPC1;iLayer++){
486 for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
489 volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
491 if(AliGeomManager::LayerToVolUID(iLayer,iModule2)==volidAl)continue;
492 volIDsFit.AddAt(volid,count);
497 volIDs.AddAt((Int_t)volidAl,0);
498 if(iter==iterations-1){
499 volIDs2.AddAt(volidAl,lastVolid);
502 volIDs2.AddAt(volidAl,lastVolid);
503 AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iter1);
507 if((iter+1)%2==0||(iter+1)%5==0||iter==0||iter==1||iter==2||iter==3||iter==iterations-1){
508 command="RealignObj";
510 command.Append(".root");
511 WriteRealignObjArray(command.Data(),AliGeomManager::kSPD1,AliGeomManager::kSSD2);
516 if(count==0){printf("count=0 \n");return;}
517 for(Int_t k=0;k<volIDs2.GetSize();k++){
518 if(volIDs2.At(k)==0)break;
519 layerNum=AliGeomManager::VolUIDToLayer(volIDs2.At(k),modNum);
520 ntVolumeAlign->Fill(layerNum,modNum,volIDs2.At(k));
522 printf("End of selecting modules cycle: %d modules selected \n",count);
523 TFile *f=new TFile("RealignVolNt.root","RECREATE");
525 ntVolumeAlign->Write();
535 //______________________________________________________________________________
536 void AliITSRealignTracks::InitAlignObjs()
538 // Initialize the alignment objects array
541 for(Int_t i=0;i<21;i++)cov[i]=0.;
542 for(Int_t i=0;i<3;i++)cov[i*(i+1)/2+i]=0.05*0.05;//Set Default Error to 500 micron for Translations
543 for(Int_t i=3;i<6;i++)cov[i*(i+1)/2+i]=0.001*0.001*180*180/3.14/3.14;//and 1 mrad for rotations (global ref. sysytem->~40 micron for SPD1,~450 micron for SSD2)
545 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
546 fAlignObjs = new AliAlignObj**[nLayers];
547 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
548 fAlignObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
549 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
550 UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
551 fAlignObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
552 fAlignObjs[iLayer][iModule]->SetCorrMatrix(cov);
553 fAlignObjs[iLayer][iModule]->SetUniqueID(0);
558 //______________________________________________________________________________
559 void AliITSRealignTracks::ResetCorrModules(){
560 // Initialize and reset to 0 the array with the information on correlation
562 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
563 fCorrModules = new Double_t*[nLayers];
564 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
565 fCorrModules[iLayer] = new Double_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
566 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
567 fCorrModules[iLayer][iModule]=0.;
572 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
573 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
574 fCorrModules[iLayer][iModule]=0.;
580 //______________________________________________________________________________
581 Bool_t AliITSRealignTracks::InitSurveyObjs(Bool_t infinite,Double_t factor,TString filename,TString arrayName){
583 if(fSurveyObjs)DeleteSurveyObjs();
584 Bool_t fromfile=kFALSE;
586 TClonesArray *clnarray;
587 if(!filename.IsNull()){
588 //Initialize from file
589 if(gSystem->AccessPathName(filename.Data(),kFileExists)){
590 printf("Wrong Survey AlignObjs File Name \n");
593 if(arrayName.IsNull()){
594 printf("Null Survey Object Name! \n");
601 // Initialize the alignment objects array with default values
602 Double_t v=1.*factor;
603 if(infinite)v*=100000.;
606 for(Int_t i=0;i<21;i++)cov[i]=0.;
607 for(Int_t i=0;i<3;i++)cov[i*(i+1)/2+i]=0.1*0.1*v;//Set Default Error to 1 mm for Translation
608 for(Int_t i=3;i<6;i++)cov[i*(i+1)/2+i]=0.03*0.03*180.*180./3.14/3.14*v;//and 30 mrad (~1.7 degrees)for rotations (global ref. sysytem)
610 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
611 fSurveyObjs = new AliAlignObj**[nLayers];
612 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
613 fSurveyObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
614 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
615 fSurveyObjs[iLayer][iModule] = 0x0;
620 surveyObj=TFile::Open(filename.Data());
621 if (!surveyObj || !surveyObj->IsOpen()) {
622 AliError(Form("Could not open SurveyObjs file: %s !",filename.Data()));
625 printf("Getting TClonesArray \n");
626 clnarray=(TClonesArray*)surveyObj->Get(arrayName);
627 Int_t size=clnarray->GetSize();
629 for(Int_t ivol=0;ivol<size;ivol++){
630 AliAlignObjParams *a=(AliAlignObjParams*)clnarray->At(ivol);
631 volid=a->GetVolUID();
633 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
634 if(iLayer<=0)continue;
635 if(a->GetUniqueID()==0)continue;
636 printf("Updating survey for volume: %d ,layer: %d module: %d from file\n",volid,iLayer,iModule);
637 fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule] = new AliAlignObjParams(*a);
638 fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->SetUniqueID(a->GetUniqueID());
639 fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->Print("");
645 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
646 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
647 UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
648 if(!fSurveyObjs[iLayer][iModule]){
649 printf("Updating survey for volume: %d ,layer: %d module: %d with default values \n",volid,iLayer,iModule);
650 fSurveyObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
651 fSurveyObjs[iLayer][iModule]->SetCorrMatrix(cov);
652 fSurveyObjs[iLayer][iModule]->SetUniqueID(0);
663 //______________________________________________________________________________
664 Int_t AliITSRealignTracks::CheckWithSurvey(Double_t factor,TArrayI *volids){
666 // Check the parameters of the alignment objects in volids (or of all objects if volids is null)
667 // are into the boundaries set by the cov. matrix of the survey objs
668 // Returns the number of objects out of boudaries
669 AliAlignObj *alignObj;
672 Double_t surveycov[21],transl[3],rot[3],survtransl[3],survrot[3];
674 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
675 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++){
676 volid=AliGeomManager::LayerToVolUIDSafe(iLayer+AliGeomManager::kFirstLayer,iModule);
677 alignObj=GetAlignObj(volid);
678 alignObj->GetPars(transl,rot);
679 fSurveyObjs[iLayer][iModule]->GetCovMatrix(surveycov);
680 fSurveyObjs[iLayer][iModule]->GetPars(survtransl,survrot);
681 if(TMath::Sqrt(TMath::Abs(surveycov[0]))*factor<TMath::Abs(transl[0]-survtransl[0])||TMath::Sqrt(TMath::Abs(surveycov[2]))*factor<TMath::Abs(transl[1]-survtransl[1])||TMath::Sqrt(TMath::Abs(surveycov[5]))*factor<TMath::Abs(transl[2]-survtransl[2])||TMath::Sqrt(TMath::Abs(surveycov[9]))*factor<TMath::Abs(rot[0]-survrot[0])||TMath::Sqrt(TMath::Abs(surveycov[14]))*factor<TMath::Abs(rot[1]-survrot[1])||TMath::Sqrt(TMath::Abs(surveycov[20]))*factor<TMath::Abs(rot[2]-survrot[2])){
682 printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
683 // *alignObj = *alignObjSurv;
684 alignObj->SetPars(survtransl[0],survtransl[1],survtransl[2],survrot[0],survrot[1],survrot[2]);
685 alignObj->SetUniqueID(0);
686 if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
695 for(Int_t j=0;j<volids->GetSize();j++){
697 alignObj=GetAlignObj(volid);
698 alignObj->GetPars(transl,rot);
699 iLayer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volid,iModule)-(Int_t)AliGeomManager::kFirstLayer;
700 fSurveyObjs[iLayer][iModule]->GetCovMatrix(surveycov);
701 fSurveyObjs[iLayer][iModule]->GetPars(survtransl,survrot);
702 if(TMath::Sqrt(TMath::Abs(surveycov[0]))*factor<TMath::Abs(transl[0]-survtransl[0])||TMath::Sqrt(TMath::Abs(surveycov[2]))*factor<TMath::Abs(transl[1]-survtransl[1])||TMath::Sqrt(TMath::Abs(surveycov[5]))*factor<TMath::Abs(transl[2]-survtransl[2])||TMath::Sqrt(TMath::Abs(surveycov[9]))*factor<TMath::Abs(rot[0]-survrot[0])||TMath::Sqrt(TMath::Abs(surveycov[14]))*factor<TMath::Abs(rot[1]-survrot[1])||TMath::Sqrt(TMath::Abs(surveycov[20]))*factor<TMath::Abs(rot[2]-survrot[2])){
703 printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
704 // *alignObj = *alignObjSurv;
705 alignObj->SetPars(survtransl[0],survtransl[1],survtransl[2],survrot[0],survrot[1],survrot[2]);
706 alignObj->SetUniqueID(0);
707 if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
715 //___________________________________________________________________
717 void AliITSRealignTracks::ResetAlignObjs(Bool_t all,TArrayI *volids)
719 // Reset the alignment objects in volids or all if all=kTRUE
721 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
722 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++)
723 fAlignObjs[iLayer][iModule]->SetPars(0,0,0,0,0,0);
729 for(Int_t j=0;j<volids->GetSize();j++){
730 layer=(Int_t)AliGeomManager::VolUIDToLayer(volids->At(j),mod)-(Int_t)AliGeomManager::kFirstLayer;
731 fAlignObjs[layer][mod]->SetPars(0,0,0,0,0,0);
736 //______________________________________________-
737 void AliITSRealignTracks::DeleteSurveyObjs()
739 if(!fSurveyObjs)return;
740 // Delete the alignment objects array
741 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
742 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++){
743 if (fSurveyObjs[iLayer][iModule]) delete fSurveyObjs[iLayer][iModule];
746 if(fSurveyObjs[iLayer])delete [] fSurveyObjs[iLayer];
749 delete [] fSurveyObjs;
754 //______________________________________________________________________________
755 Bool_t AliITSRealignTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName){
757 // Read alignment object from a file: update the alignobj already present with the one in the file
758 // To be replaced by a call to CDB
760 if(gSystem->AccessPathName(alignObjFileName,kFileExists)){
761 printf("Wrong AlignObjs File Name \n");
765 TFile *fRealign=TFile::Open(alignObjFileName);
766 if (!fRealign || !fRealign->IsOpen()) {
767 AliError(Form("Could not open Align Obj File file %s !",alignObjFileName));
770 printf("Getting TClonesArray \n");
771 TClonesArray *clnarray=(TClonesArray*)fRealign->Get(arrayName);
772 Int_t size=clnarray->GetSize();
775 for(Int_t ivol=0;ivol<size;ivol++){
776 AliAlignObjParams *a=(AliAlignObjParams*)clnarray->At(ivol);
777 volid=a->GetVolUID();
779 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
780 if(iLayer<AliGeomManager::kFirstLayer||iLayer>AliGeomManager::kSSD2)continue;
781 printf("Updating volume: %d ,layer: %d module: %d \n",volid,iLayer,iModule);
782 *fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule] *= *a;
790 //_________________________________________
791 Bool_t AliITSRealignTracks::FirstAlignmentLayers(Bool_t *layers,Int_t minNtracks,Int_t iterations,Bool_t fitall,TArrayI *volidsSet){
793 //Align all modules in the set of layers independently according to a sequence based on the number of tracks passing through a given module
796 TString name="DrawFirstAlignment_Layers";
801 Int_t maxntr=0,nMod=0,modAligned=0,size=0;
802 for(Int_t i=0;i<6;i++){
804 size+=AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);
809 TArrayI *volFit=new TArrayI(size);
810 TArrayI *volFit2=new TArrayI(size-1);
811 TArrayI *sequence=new TArrayI(size);
812 TArrayI *volIn=new TArrayI(1);
814 // Initialize the index arrays
815 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
816 lastIndex = new Int_t*[nLayers];
817 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
818 lastIndex[iLayer] = new Int_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
819 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
820 lastIndex[iLayer][iModule] = fLastIndex[iLayer][iModule];
821 if(iLayer<=(AliGeomManager::kSSD2-AliGeomManager::kFirstLayer)&&layers[iLayer]==1){
822 volFit->AddAt(AliGeomManager::LayerToVolUID(iLayer+AliGeomManager::kFirstLayer,iModule),maxntr);
829 while (maxntr>minNtracks){
831 for (Int_t iLayer = 0; iLayer <= (AliGeomManager::kSSD2 - AliGeomManager::kFirstLayer); iLayer++) {
832 if(layers[iLayer]==0)continue;
833 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
834 if(lastIndex[iLayer][iModule]>maxntr){
835 maxntr=lastIndex[iLayer][iModule];
841 if(maxntr>minNtracks){
842 voluid=AliGeomManager::LayerToVolUID(laymax+AliGeomManager::kFirstLayer,modmax);
843 sequence->AddAt(voluid,nMod);
844 lastIndex[laymax][modmax]=0;
852 for(Int_t iter=0;iter<iterations;iter++){
853 if(iter>0&&fDraw)UpdateDraw(sequence,iter,iter);
855 for(Int_t k=0;k<nMod;k++){
857 voluid=sequence->At(k);
858 ilayer=AliGeomManager::VolUIDToLayer(voluid,imod);
859 volIn->AddAt(voluid,0);
862 for(Int_t j=0;j<nMod;j++){
867 else volFit2->AddAt(sequence->At(j),j-found);
869 volFit2->Set(nMod-1);
872 for(Int_t j=0;j<volFit->GetSize();j++){
873 if(volFit->At(j)!=volIn->At(0))volFit2->AddAt(volFit->At(j),j-found);
879 volFit3=IntersectVolArray(volidsSet,volFit2);
881 else volFit3=new TArrayI(*volFit2);
884 if(AlignVolumesITS(volIn,volFit3,AliGeomManager::kSPD1,AliGeomManager::kTPC1,2))modAligned++;
889 Int_t noutofsurv=CheckWithSurvey(2.,sequence);
890 printf("%d modules into the sequence \n %d modules re-aligned \n %d modules moved far away from survey (-> reset) \n",nMod,modAligned,noutofsurv);
891 name.Append("_iter");
893 name.Append(".root");
894 if(fDraw)WriteHists(name.Data());
903 //__________________________________________
904 Bool_t AliITSRealignTracks::FirstAlignmentSPD(Int_t minNtracks,Int_t iterations,Bool_t fitall,TArrayI *volidsSet){
912 Int_t maxntr=0,nMod=0,modAligned=0;
913 TArrayI *volFit=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2));
914 TArrayI *volFit2=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2)-1);
915 TArrayI *sequence=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2));
916 TArrayI *volIn=new TArrayI(1);
918 // Initialize the index arrays
919 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
920 lastIndex = new Int_t*[nLayers];
921 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
922 lastIndex[iLayer] = new Int_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
923 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
924 lastIndex[iLayer][iModule] = fLastIndex[iLayer][iModule];
925 if(iLayer+AliGeomManager::kFirstLayer<=AliGeomManager::kSPD2){
926 volFit->AddAt(AliGeomManager::LayerToVolUID(iLayer+AliGeomManager::kFirstLayer,iModule),maxntr);
933 while (maxntr>minNtracks){
935 for (Int_t iLayer = 0; iLayer <= (AliGeomManager::kSPD2 - AliGeomManager::kFirstLayer); iLayer++) {
936 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
937 if(lastIndex[iLayer][iModule]>maxntr){
940 maxntr=lastIndex[iLayer][iModule];
944 if(maxntr>minNtracks){
945 voluid=AliGeomManager::LayerToVolUID(laymax+AliGeomManager::kFirstLayer,modmax);
946 sequence->AddAt(voluid,nMod);
947 lastIndex[laymax][modmax]=0;
949 volIn->AddAt(voluid,0);
955 for(Int_t iter=0;iter<iterations;iter++){
957 for(Int_t k=0;k<nMod;k++){
959 voluid=sequence->At(k);
960 ilayer=AliGeomManager::VolUIDToLayer(voluid,imod);
961 volIn->AddAt(voluid,0);
964 for(Int_t j=0;j<nMod;j++){
969 else volFit2->AddAt(sequence->At(j),j-found);
971 volFit2->Set(nMod-1);
974 for(Int_t j=0;j<volFit->GetSize();j++){
975 if(volFit->At(j)!=volIn->At(0))volFit2->AddAt(volFit->At(j),j-found);
981 volFit3=IntersectVolArray(volidsSet,volFit2);
983 else volFit3=new TArrayI(*volFit2);
986 if(AlignVolumesITS(volIn,volFit3,AliGeomManager::kSPD1,AliGeomManager::kSDD1,2))modAligned++;
988 // if(volidsSet)delete volFit3;
991 Int_t noutofsurv=CheckWithSurvey(2.,sequence);
992 printf("%d modules into the sequence \n %d modules re-aligned \n %d modules moved far away from survey (-> reset) \n",nMod,modAligned,noutofsurv);
1001 //__________________________________
1002 Bool_t AliITSRealignTracks::SPDmodulesAlignToSSD(Int_t minNtracks,Int_t iterations){
1004 Int_t volSSD[6]={0,0,0,0,1,1};
1005 TArrayI *volOuter=GetLayersVolUID(volSSD);
1006 TArrayI *voluid=new TArrayI(1);
1007 for(Int_t iter=0;iter<iterations;iter++){
1009 for(Int_t imod=0;imod<AliGeomManager::LayerSize(AliGeomManager::kSPD1);imod++){ if(GetLastIndex(AliGeomManager::kSPD1-AliGeomManager::kFirstLayer,imod)<minNtracks){
1010 printf("Not enough tracks for module: lay %d mod %d \n",1,imod );
1013 voluid->AddAt(AliGeomManager::LayerToVolUID(AliGeomManager::kSPD1,imod),0);
1014 AlignVolumesITS(voluid,volOuter,AliGeomManager::kSSD1,AliGeomManager::kSSD2,2);
1017 for(Int_t imod=0;imod<AliGeomManager::LayerSize(AliGeomManager::kSPD2);imod++){
1018 if(GetLastIndex(AliGeomManager::kSPD2-AliGeomManager::kFirstLayer,imod)<minNtracks){
1019 printf("Not enough tracks for module: lay %d mod %d \n",2,imod );
1022 voluid->AddAt(AliGeomManager::LayerToVolUID(AliGeomManager::kSPD2,imod),0);
1023 AlignVolumesITS(voluid,volOuter,AliGeomManager::kSSD1,AliGeomManager::kSSD2,2);
1029 //______________________________________________________________________________
1030 Bool_t AliITSRealignTracks::AlignVolumesITS(const TArrayI *volids, const TArrayI *volidsfit,
1031 AliGeomManager::ELayerID layerRangeMin,
1032 AliGeomManager::ELayerID layerRangeMax,
1035 // Align a set of detector volumes.
1036 // Tracks are fitted only within
1037 // the range defined by the user
1038 // (by layerRangeMin and layerRangeMax)
1039 // or within the set of volidsfit
1040 // Repeat the procedure 'iterations' times
1042 Int_t nVolIds = volids->GetSize();
1044 AliError("Volume IDs array is empty!");
1047 Bool_t correlated=kFALSE;
1048 Double_t surveycov[21],transl[3],rot[3],survtransl[3],survrot[3];
1052 Double_t phiglob,smearing,rotorig[9],normplanevect[3]={0.,0.,0.},normplanevect2[3]={0.,0.,0.};
1054 TMatrixDSym covmatrx(6);
1055 AliAlignObj *modAlign;
1057 // Load only the tracks with at least one
1058 // space point in the set of volume (volids)
1060 AliTrackPointArray **points;
1061 Bool_t failed=kFALSE;
1062 Int_t pointsdim=0,skipped=0;
1063 // Start the iterations
1064 while (iterations > 0){
1065 normplanevect2[0]=0.;
1066 normplanevect2[1]=0.;
1067 normplanevect2[2]=0.;
1072 Int_t nArrays = LoadPoints(volids, points,pointsdim);
1074 if (nArrays < fmintracks||nArrays<=0){
1076 printf("Not enough tracks to try minimization: volUID %d and following in volids \n", volids->At(0));
1077 UnloadPoints(pointsdim, points);
1080 frac=1./(Double_t)nArrays;
1081 AliTrackResiduals *minimizer = CreateMinimizer();
1082 minimizer->SetNTracks(nArrays);
1083 minimizer->InitAlignObj();
1084 AliTrackFitter *fitter = CreateFitter();
1086 //Here prepare to set the plane for GetPCArot
1087 // if(volids->GetSize()==1){//TEMPORARY: to be improved
1088 AliGeomManager::GetOrigRotation(volids->At(0),rotorig);
1089 if((Int_t)AliGeomManager::VolUIDToLayer(volids->At(0))==1){//TEMPORARY: to be improved
1090 normplanevect[0]=-rotorig[1];
1091 normplanevect[1]=-rotorig[4];
1092 normplanevect[2]=0.;
1095 normplanevect[0]=rotorig[1];
1096 normplanevect[1]=rotorig[4];
1097 normplanevect[2]=0.;
1100 phiglob=TMath::ATan2(normplanevect[1],normplanevect[0]);
1102 modAlign=GetAlignObj(volids->At(0));
1103 modAlign->GetMatrix(hM);
1104 deltarot=hM.GetRotationMatrix();
1105 for(Int_t j=0;j<3;j++){
1106 for(Int_t i=0;i<3;i++){
1107 normplanevect2[j]+=deltarot[j*3+i]*normplanevect[i];
1109 // printf("Here the difference: norm1[%d]=%f norm2[%d]=%f \n",j,normplanevect[j],j,normplanevect2[j]);
1113 if(modAlign->GetUniqueID()==0)smearing=fsigmaY;
1115 modAlign->GetCovMatrix(covmatrx);
1116 smearing=5.*5.*(covmatrx(0,0)+covmatrx(1,1)+covmatrx(2,2)+10.*10.*covmatrx(3,3)+10.*10.*covmatrx(4,4)+10.*10.*covmatrx(5,5))/6.;
1117 //This is a sort of average: the trace with the variances of the angles
1118 //weighted with 10 cm divided per 6 and the result multiplied per 25
1119 // (the sqrt would be 5 times a sort of "mean sigma" )
1123 else smearing=fsigmaY;
1124 printf("This is the sigmaY value: %f \n",smearing);
1125 // the plane will be set into the loop on tracks
1128 for (Int_t iArray = 0; iArray < nArrays; iArray++) {
1129 if (!points[iArray]) continue;
1130 points[iArray]->Sort(kTRUE);
1131 fitter->SetTrackPointArray(points[iArray], kFALSE);
1132 // printf("Here normplane vect: %f \n",normplanevect2[1]); //TO BE REPLACED BY fitter->SetNormPlaneVect(normplanevect2);
1133 if (fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE) continue;
1139 TArrayI *volparray=new TArrayI(points[iArray]->GetNPoints());
1140 for(Int_t point=0;point<points[iArray]->GetNPoints();point++){
1141 points[iArray]->GetPoint(p,point);
1142 volparray->AddAt(p.GetVolumeID(),point);
1144 TArrayI *volpArray=ExcludeVolidsFromVolidsArray(volids,volparray);
1145 for(Int_t point=0;point<volpArray->GetSize();point++){
1146 layer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volpArray->At(point),module);
1147 if(fCorrModules[layer-AliGeomManager::kFirstLayer][module]>fLimitCorr){
1149 // printf("volid %d, iarray = %d : skipping %d for Volume: %d \n",volids->At(0),iArray,skipped,volpArray->At(point));
1155 for(Int_t point=0;point<volpArray->GetSize();point++){
1156 layer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volpArray->At(point),module);
1157 //printf("Number of common tracks: %d \n",fCorrModules[layer-AliGeomManager::kFirstLayer][module]);
1158 fCorrModules[layer-AliGeomManager::kFirstLayer][module]+=frac;
1170 AliTrackPointArray *pVolId,*pTrack;
1171 fitter->GetTrackResiduals(pVolId,pTrack);
1172 minimizer->AddTrackPointArrays(pVolId,pTrack);
1175 printf("Number of tracks considered: %d \n",nArrays);
1176 frac=(Double_t)skipped/(Double_t)nArrays;
1177 printf("Number of tracks skipped cause of correlation: %d (fraction: %f )\n",skipped,frac);
1179 Int_t ntracks=minimizer->GetNFilledTracks();
1180 frac=(Double_t)ntracks/(Double_t)nArrays;
1181 printf("Number of tracks into the minimizer: %d (fraction: %f )\n",ntracks,frac);
1182 if(ntracks<=fmintracks){
1183 printf("Not enough good tracks found: could not find parameter for volume %d (and following in volids)\n",volids->At(0));
1184 UnloadPoints(pointsdim, points);
1189 failed=(!minimizer->Minimize());
1191 // Update the alignment object(s)
1192 if (fDoUpdate) for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
1193 UShort_t volid = (*volids)[iVolId];
1196 AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
1197 AliAlignObj *alignObj = fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule];
1199 //Check the last minimization is not too large
1200 minimizer->GetAlignObj()->GetPars(transl,rot);
1201 fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->GetCovMatrix(surveycov);
1202 fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->GetPars(survtransl,survrot);
1203 if(TMath::Sqrt(TMath::Abs(surveycov[0]))*2<TMath::Abs(transl[0]-survtransl[0])||TMath::Sqrt(TMath::Abs(surveycov[2]))*2<TMath::Abs(transl[1]-survtransl[1])||TMath::Sqrt(TMath::Abs(surveycov[5]))*2<TMath::Abs(transl[2]-survtransl[2])||TMath::Sqrt(TMath::Abs(surveycov[9]))*2<TMath::Abs(rot[0]-survrot[0])||TMath::Sqrt(TMath::Abs(surveycov[14]))*2<TMath::Abs(rot[1]-survrot[1])||TMath::Sqrt(TMath::Abs(surveycov[20]))*2<TMath::Abs(rot[2]-survrot[2])){
1204 printf("Results for module %d too large: can't update them \n",volid);
1205 alignObj->SetUniqueID(2);
1212 *alignObj *= *minimizer->GetAlignObj();
1213 alignObj->SetUniqueID(1);
1216 alignObj->GetCovMatrix(covmatrx);
1217 *alignObj *= *minimizer->GetAlignObj();
1218 alignObj->SetCorrMatrix(covmatrx);
1219 alignObj->SetUniqueID(1);
1222 /*alignObj->GetPars(transl,rot);
1224 if(TMath::Sqrt(TMath::Abs(surveycov[0]))*20<TMath::Abs(transl[0])||TMath::Sqrt(TMath::Abs(surveycov[2]))*20<TMath::Abs(transl[1])||TMath::Sqrt(TMath::Abs(surveycov[5]))*20<TMath::Abs(transl[2])||TMath::Sqrt(TMath::Abs(surveycov[9]))*20<TMath::Abs(rot[0])||TMath::Sqrt(TMath::Abs(surveycov[14]))*20<TMath::Abs(rot[1])||TMath::Sqrt(TMath::Abs(surveycov[20]))*20<TMath::Abs(rot[2])){
1225 printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
1226 // *alignObj = *alignObjSurv;
1227 alignObj->SetPars(0.,0.,0.,0.,0.,0.);
1228 alignObj->SetUniqueID(0);
1229 if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
1235 if(iterations==1)alignObj->Print("");
1238 printf("Minimization failed: cannot update AlignObj for volume: %d \n",volid);
1241 UnloadPoints(pointsdim,points);
1243 minimizer->InitAlignObj();
1254 //______________________________________________
1255 Bool_t AliITSRealignTracks::AlignSPDBarrel(Int_t iterations){
1257 Int_t size=0,size2=0;
1258 Int_t layers[6]={1,1,0,0,0,0};
1259 for(Int_t k=1;k<=2;k++){
1260 size+=AliGeomManager::LayerSize(k);
1262 for(Int_t k=3;k<=6;k++){
1263 size2+=AliGeomManager::LayerSize(k);
1264 printf("size: %d \n",size2);
1267 printf("Aligning SPDBarrel: nmodules: %d \n",size);
1268 printf("Fitting modules: %d \n",size2);
1270 TArrayI *volIDs=GetLayersVolUID(layers);
1277 TArrayI *volIDsFit=GetLayersVolUID(layers);
1279 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSDD1,AliGeomManager::kTPC1,iterations);
1284 //______________________
1285 Bool_t AliITSRealignTracks::AlignSPDHalfBarrel(Int_t method,Int_t iterations){
1287 Int_t size=0,size2=0;
1288 Int_t layers[6]={0,0,1,1,1,1};
1289 Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1290 Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1293 if(method==0)updownstr="UpNoDown";
1294 else if (method==1)updownstr="DownNoUp";
1295 else if (method==10)updownstr="UpWithDown";
1296 else if (method==11)updownstr="DownWithUp";
1298 AliWarning("Wrong AlignSPDHalfBarrel method selected ");
1302 for(Int_t i=1;i<=2;i++){
1303 size+=AliGeomManager::LayerSize(i);
1306 for(Int_t i=3;i<=6;i++){
1307 size2+=AliGeomManager::LayerSize(i);
1311 if(method==10||method==11)size2+=size;
1313 printf("Aligning SPDHalfBarrel %s: nmodules: %d \n",updownstr.Data(),size);
1314 printf("Fitting modules: %d \n",size2);
1315 TArrayI *volIDsFit2;
1316 TArrayI *volids = NULL;
1317 TArrayI *volIDsFit=GetLayersVolUID(layers);
1318 if(method==0||method==10)volids=GetSPDSectorsVolids(sectorsUp);
1319 if(method==1||method==11)volids=GetSPDSectorsVolids(sectorsDown);
1321 if(method==10)volIDsFit2=JoinVolArrays(GetSPDSectorsVolids(sectorsDown),volIDsFit);
1322 else if(method==11)volIDsFit2=JoinVolArrays(GetSPDSectorsVolids(sectorsUp),volIDsFit);
1323 else volIDsFit2=volIDsFit;
1325 AlignVolumesITS(volids,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1331 //______________________________________________________
1332 Bool_t AliITSRealignTracks::AlignLayer(Int_t layer,Int_t iterations){
1334 Int_t size=0,size2=0;
1335 Int_t layers[6]={0,0,0,0,0,0};
1337 TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1338 for(Int_t k=1;k<=6;k++){
1339 if(k!=layer)size2+=AliGeomManager::LayerSize(k);
1341 size=AliGeomManager::LayerSize(layer);
1343 printf("Aligning layer %s, nmodules %d ,fitted modules %d \n",layerstr[layer-1].Data(),size,size2);
1346 TArrayI *volIDs=GetLayersVolUID(layers);
1354 TArrayI *volIDsFit=GetLayersVolUID(layers);
1356 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSDD1,AliGeomManager::kSSD2,iterations);
1361 //___________________________________________
1363 Bool_t AliITSRealignTracks::AlignLayersToLayers(Int_t *layer,Int_t iterations){
1366 Int_t size=0,size2=0,j=0,k=0;
1368 TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1369 TString command="",str;
1370 for(Int_t i=1;i<=6;i++){
1371 if(layer[i-1]==1||layer[i-1]==3){
1372 size+=AliGeomManager::LayerSize(i);
1373 command.Append(" ");
1374 command.Append(layerstr[i-1]);
1376 if(layer[i-1]==2||layer[i-1]==3){
1377 size2+=AliGeomManager::LayerSize(i);
1379 str.Append(layerstr[i-1]);
1383 printf("Aligning layers %s To layers %s, nmodules %d ,fitted modules %d \n",command.Data(),str.Data(),size,size2);
1386 TArrayI volIDs(size);
1387 TArrayI volIDsFit(size2);
1389 for (iLayer=(Int_t)AliGeomManager::kSPD1;iLayer<(Int_t)AliGeomManager::kTPC1;iLayer++){
1390 if(layer[iLayer-AliGeomManager::kFirstLayer]==0)continue;
1391 if(layer[iLayer-AliGeomManager::kFirstLayer]==1){
1392 for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1393 volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1394 volIDs.AddAt(volid,j);
1398 else if(layer[iLayer-AliGeomManager::kFirstLayer]==2){
1399 for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1400 volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1401 volIDsFit.AddAt(volid,k);
1405 else if(layer[iLayer-AliGeomManager::kFirstLayer]==3){
1406 for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1407 volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1408 volIDs.AddAt(volid,j);
1410 volIDsFit.AddAt(volid,k);
1416 AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1421 //______________________________________________
1423 Bool_t AliITSRealignTracks::AlignSPDSectorToOuterLayers(Int_t sector,Int_t iterations){
1426 Int_t layers[6]={0,0,1,1,1,1};
1428 Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1429 Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1435 sectorsIN[sector]=1;
1436 sectorsFit[sector]=0;
1437 TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1440 volIDsFit=JoinVolArrays(GetSPDSectorsVolids(sectorsFit),GetLayersVolUID(layers));
1442 else volIDsFit=GetLayersVolUID(layers);
1444 printf("Aligning SPD sector %d: nmodules: %d \n",sector,volIDs->GetSize());
1445 printf("Fitting modules: %d \n",volIDsFit->GetSize());
1447 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1452 //______________________________________________
1453 Bool_t AliITSRealignTracks::AlignSPDSectorWithSectors(Int_t sector,Int_t iterations){
1456 Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1457 Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1459 sectorsIN[sector]=1;
1460 sectorsFit[sector]=0;
1461 TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1462 TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);;
1464 printf("Aligning SPD sector %d: nmodules: %d \n",sector,volIDs->GetSize());
1465 printf("Fitting modules: %d \n",volIDsFit->GetSize());
1468 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSDD1,iterations);
1475 //___________________________________________________
1476 Bool_t AliITSRealignTracks::AlignSPDSectorsWithSectors(Int_t *sectorsIN,Int_t *sectorsFit,Int_t iterations){
1478 TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1479 TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);;
1481 printf("Aligning SPD sectors: modules: %d \n",volIDs->GetSize());
1482 printf("Fitting modules: %d \n",volIDsFit->GetSize());
1486 return AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSDD1,iterations);;
1489 //___________________________________________________
1490 Bool_t AliITSRealignTracks::AlignSPDStaves(Int_t *staves,Int_t *sectorsIN,Int_t *sectorsFit,Int_t iterations){
1492 TArrayI *volIDs=GetSPDStavesVolids(sectorsIN,staves);
1493 TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1495 if(volIDs->GetSize()==0){
1496 printf("EMPTY ARRAY !! \n");
1499 printf("Aligning SPD staves: modules: %d \n",volIDs->GetSize());
1500 printf("Fitting modules: %d \n",volIDsFit->GetSize());
1502 TArrayI *volIDsFit2=ExcludeVolidsFromVolidsArray(volIDs,volIDsFit);
1503 return AlignVolumesITS(volIDs,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSSD1,iterations);
1507 //___________________________________________
1509 Bool_t AliITSRealignTracks::AlignLayerToSPDHalfBarrel(Int_t layer,Int_t updown,Int_t iterations){
1512 Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1513 Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1514 TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1516 Int_t layers[6]={0,0,0,0,0,0};
1518 Int_t size=AliGeomManager::LayerSize(layer);
1519 TArrayI *volIDs=GetLayersVolUID(layers);
1522 volIDsFit=GetSPDSectorsVolids(sectorsUp);
1523 printf("Aligning layer %s, nmodules %d ,to half barrel Up \n",layerstr[layer-1].Data(),size);
1526 volIDsFit=GetSPDSectorsVolids(sectorsDown);
1527 printf("Aligning layer %s, nmodules %d ,to half barrel Down \n",layerstr[layer-1].Data(),size);
1530 printf("Wrong Half Barrel selection! \n");
1534 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1539 //___________________________________________
1541 Bool_t AliITSRealignTracks::AlignLayerToSector(Int_t layer,Int_t sector,Int_t iterations){
1544 printf("Wrong Sector selection! \n");
1547 Int_t sectors[10]={0,0,0,0,0,0,0,0,0,0};
1549 TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1551 Int_t layers[6]={0,0,0,0,0,0};
1553 TArrayI *volIDs=GetLayersVolUID(layers);
1554 Int_t size=AliGeomManager::LayerSize(layer);
1557 volIDsFit=GetSPDSectorsVolids(sectors);
1558 printf("Aligning layer %s, nmodules %d ,to half barrel Up \n",layerstr[layer-1].Data(),size);
1560 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1565 //_______________________________________________
1567 Bool_t AliITSRealignTracks::AlignSPDHalfBarrelToHalfBarrel(Int_t updown,Int_t iterations){
1570 Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1571 Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1573 TArrayI *volIDsUp=GetSPDSectorsVolids(sectorsUp);
1574 TArrayI *volIDsDown=GetSPDSectorsVolids(sectorsDown);
1577 printf("Aligning SPD HalfBarrel up to half Barrel down : nmodules: %d \n",volIDsUp->GetSize());
1578 printf("Fitting modules: %d \n",volIDsDown->GetSize());
1579 AlignVolumesITS(volIDsUp,volIDsDown,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1582 printf("Aligning SPD HalfBarrel down to half Barrel Up : nmodules: %d \n",volIDsDown->GetSize());
1583 printf("Fitting modules: %d \n",volIDsUp->GetSize());
1584 AlignVolumesITS(volIDsDown,volIDsUp,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1587 printf("Wrong Half Barrel selection! \n");
1595 //_______________________
1596 Bool_t AliITSRealignTracks::AlignSPDHalfBarrelToSectorRef(Int_t sector,Int_t iterations){
1598 Int_t sectorsIN[10]={0,0,0,0,0,1,1,1,1,1};
1599 Int_t sectorsFit[10]={0,0,0,0,0,0,0,0,0,0};
1601 sectorsFit[sector]=1;
1603 TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1604 TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1606 printf("Aligning SPD HalfBarrel to sector 0 %d: nmodules: %d \n",sector,volIDs->GetSize());
1607 printf("Fitting modules: %d \n",volIDsFit->GetSize());
1610 AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1615 //_________________________________________
1616 Bool_t AliITSRealignTracks::AlignSPD1SectorRef(Int_t sector,Int_t iterations){
1618 Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1619 Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1620 sectorsIN[sector]=1;
1621 sectorsFit[sector]=0;
1622 TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1623 TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1624 Int_t size=volIDs->GetSize();
1625 Int_t size2=volIDsFit->GetSize();
1629 TArrayI *volIDsSPD1=new TArrayI(size-8);
1630 TArrayI *volIDsFit2=new TArrayI(size2+8);
1632 for(Int_t j=0;j<size;j++){
1633 volID=volIDs->At(j);
1634 if(AliGeomManager::VolUIDToLayer(volID)==AliGeomManager::kSPD1){
1635 volIDsSPD1->AddAt(volID,size2+k);
1638 else volIDsFit2->AddAt(volID,j-k);
1642 for(Int_t j=0;j<size2;j++){
1643 volID=volIDsFit->At(j);
1644 volIDsFit2->AddAt(volID,size-k+j);
1647 printf("Aligning SPD Sector %d: nmodules: %d \n",sector,volIDsSPD1->GetSize());
1648 printf("Fitting modules: %d \n",volIDsFit2->GetSize());
1650 AlignVolumesITS(volIDsSPD1,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1655 //_____________________________________________
1657 AliAlignObjParams* AliITSRealignTracks::MediateAlignObj(TArrayI *volIDs,Int_t lastVolid){
1661 Double_t *rot,*transl;
1662 Double_t rotSum[9],translSum[3]={0.,0.,0.};
1663 for(Int_t k=0;k<8;k++)rotSum[k]=0.;
1666 for(Int_t ivol=0;ivol<lastVolid;ivol++){
1667 volid=volIDs->At(ivol);
1669 GetAlignObj(volIDs->At(ivol))->GetMatrix(hm);
1671 rot=hm.GetRotationMatrix();
1672 transl=hm.GetTranslation();
1674 for(Int_t j=0;j<9;j++)rotSum[j]+=rot[j];
1675 for(Int_t jt=0;jt<3;jt++)translSum[jt]+=transl[jt];
1678 for(Int_t j=0;j<9;j++)rotSum[j]=rotSum[j]/lastVolid;
1679 for(Int_t jt=0;jt<3;jt++)translSum[jt]=translSum[jt]/lastVolid;
1681 else printf("Try to mediate results for zero modules \n");
1683 hm.SetRotation(rotSum);
1684 hm.SetTranslation(translSum);
1688 AliAlignObjParams *alignObj=new AliAlignObjParams("average", 0,hm, kTRUE);
1694 //________________________________________________
1695 TArrayI* AliITSRealignTracks::GetSPDStavesVolids(Int_t *sectors,Int_t* staves){
1698 // This method gets the volID Array for the chosen staves into the
1699 // chosen sectors. You have to pass an array (10 dim) with a 1 for each
1700 // selected sector and an array (6 dim) with a 1 for each chosen stave.
1701 // The staves are numbered in this way: 0,1 for SPD1 and 2,3,4,5 for SPD2
1702 // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
1703 // staves[6]={0,1,1,0,0,1} -> Staves 1 on SPD1 and 0 and 3 on SPD2 selected
1705 Int_t nSect=0,nStaves=0;
1709 for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1710 if(sectors[co]==1) nSect++;
1713 for(Int_t co=0;co<6;co++){ //counts the number of sectors chosen
1714 if(staves[co]==1) nStaves++;
1717 if(nSect<1||nStaves<1){ //if no sector chosen -> exit
1718 Printf("Error! No Sector/s or staves Selected!");
1722 TArrayI *volIDs = new TArrayI(nSect*nStaves*4);
1723 TString stave="/Stave",str,symn,laystr;
1725 TArrayI *sectvol=GetSPDSectorsVolids(sectors);
1728 for(Int_t k=0;k<2;k++){
1732 for(Int_t i=0;i<sectvol->GetSize();i++){
1733 symn=AliGeomManager::SymName(sectvol->At(i));
1734 if(symn.Contains(str)&&symn.Contains(laystr)){
1735 // printf("Adding: %s \n",symn.Data());
1736 volIDs->AddAt(sectvol->At(i),last);
1744 for(Int_t k=2;k<6;k++){
1748 for(Int_t i=0;i<sectvol->GetSize();i++){
1749 symn=AliGeomManager::SymName(sectvol->At(i));
1750 if(symn.Contains(str)&&symn.Contains(laystr)){
1751 volIDs->AddAt(sectvol->At(i),last);
1752 printf("Adding: %s \n",symn.Data());
1763 //________________________________________________
1764 TArrayI* AliITSRealignTracks::GetSPDSectorsVolids(Int_t *sectors)
1767 // This method gets the volID Array for the chosen sectors.
1768 // You have to pass an array with a 1 for each selected sector.
1769 // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
1776 for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1777 if(sectors[co]==1) nSect++;
1780 if(nSect<1){ //if no sector chosen -> exit
1781 Printf("Error! No Sector/s Selected!");
1785 TArrayI *volIDs = new TArrayI(nSect*24);
1787 if(sectors[0]==1){ //--->cSect = 0 <---
1788 volIDs->AddAt(2048,iModule); iModule++;
1789 volIDs->AddAt(2049,iModule); iModule++;
1790 volIDs->AddAt(2050,iModule); iModule++;
1791 volIDs->AddAt(2051,iModule); iModule++;
1792 volIDs->AddAt(2052,iModule); iModule++;
1793 volIDs->AddAt(2053,iModule); iModule++;
1794 volIDs->AddAt(2054,iModule); iModule++;
1795 volIDs->AddAt(2055,iModule); iModule++;
1796 volIDs->AddAt(4096,iModule); iModule++;
1797 volIDs->AddAt(4097,iModule); iModule++;
1798 volIDs->AddAt(4098,iModule); iModule++;
1799 volIDs->AddAt(4099,iModule); iModule++;
1800 volIDs->AddAt(4100,iModule); iModule++;
1801 volIDs->AddAt(4101,iModule); iModule++;
1802 volIDs->AddAt(4102,iModule); iModule++;
1803 volIDs->AddAt(4103,iModule); iModule++;
1804 volIDs->AddAt(4104,iModule); iModule++;
1805 volIDs->AddAt(4105,iModule); iModule++;
1806 volIDs->AddAt(4106,iModule); iModule++;
1807 volIDs->AddAt(4107,iModule); iModule++;
1808 volIDs->AddAt(4108,iModule); iModule++;
1809 volIDs->AddAt(4109,iModule); iModule++;
1810 volIDs->AddAt(4110,iModule); iModule++;
1811 volIDs->AddAt(4111,iModule); iModule++;
1813 if(sectors[1]==1){ //--->cSect = 1 <//---
1814 volIDs->AddAt(2056,iModule); iModule++;
1815 volIDs->AddAt(2057,iModule); iModule++;
1816 volIDs->AddAt(2058,iModule); iModule++;
1817 volIDs->AddAt(2059,iModule); iModule++;
1818 volIDs->AddAt(2060,iModule); iModule++;
1819 volIDs->AddAt(2061,iModule); iModule++;
1820 volIDs->AddAt(2062,iModule); iModule++;
1821 volIDs->AddAt(2063,iModule); iModule++;
1822 volIDs->AddAt(4112,iModule); iModule++;
1823 volIDs->AddAt(4113,iModule); iModule++;
1824 volIDs->AddAt(4114,iModule); iModule++;
1825 volIDs->AddAt(4115,iModule); iModule++;
1826 volIDs->AddAt(4116,iModule); iModule++;
1827 volIDs->AddAt(4117,iModule); iModule++;
1828 volIDs->AddAt(4118,iModule); iModule++;
1829 volIDs->AddAt(4119,iModule); iModule++;
1830 volIDs->AddAt(4120,iModule); iModule++;
1831 volIDs->AddAt(4121,iModule); iModule++;
1832 volIDs->AddAt(4122,iModule); iModule++;
1833 volIDs->AddAt(4123,iModule); iModule++;
1834 volIDs->AddAt(4124,iModule); iModule++;
1835 volIDs->AddAt(4125,iModule); iModule++;
1836 volIDs->AddAt(4126,iModule); iModule++;
1837 volIDs->AddAt(4127,iModule); iModule++;
1839 if(sectors[2]==1){//--->cSect = 2 <//---
1840 volIDs->AddAt(2064,iModule); iModule++;
1841 volIDs->AddAt(2065,iModule); iModule++;
1842 volIDs->AddAt(2066,iModule); iModule++;
1843 volIDs->AddAt(2067,iModule); iModule++;
1844 volIDs->AddAt(2068,iModule); iModule++;
1845 volIDs->AddAt(2069,iModule); iModule++;
1846 volIDs->AddAt(2070,iModule); iModule++;
1847 volIDs->AddAt(2071,iModule); iModule++;
1848 volIDs->AddAt(4128,iModule); iModule++;
1849 volIDs->AddAt(4129,iModule); iModule++;
1850 volIDs->AddAt(4130,iModule); iModule++;
1851 volIDs->AddAt(4131,iModule); iModule++;
1852 volIDs->AddAt(4132,iModule); iModule++;
1853 volIDs->AddAt(4133,iModule); iModule++;
1854 volIDs->AddAt(4134,iModule); iModule++;
1855 volIDs->AddAt(4135,iModule); iModule++;
1856 volIDs->AddAt(4136,iModule); iModule++;
1857 volIDs->AddAt(4137,iModule); iModule++;
1858 volIDs->AddAt(4138,iModule); iModule++;
1859 volIDs->AddAt(4139,iModule); iModule++;
1860 volIDs->AddAt(4140,iModule); iModule++;
1861 volIDs->AddAt(4141,iModule); iModule++;
1862 volIDs->AddAt(4142,iModule); iModule++;
1863 volIDs->AddAt(4143,iModule); iModule++;
1865 if(sectors[3]==1){//--->cSect = 3 <//---
1866 volIDs->AddAt(2072,iModule); iModule++;
1867 volIDs->AddAt(2073,iModule); iModule++;
1868 volIDs->AddAt(2074,iModule); iModule++;
1869 volIDs->AddAt(2075,iModule); iModule++;
1870 volIDs->AddAt(2076,iModule); iModule++;
1871 volIDs->AddAt(2077,iModule); iModule++;
1872 volIDs->AddAt(2078,iModule); iModule++;
1873 volIDs->AddAt(2079,iModule); iModule++;
1874 volIDs->AddAt(4144,iModule); iModule++;
1875 volIDs->AddAt(4145,iModule); iModule++;
1876 volIDs->AddAt(4146,iModule); iModule++;
1877 volIDs->AddAt(4147,iModule); iModule++;
1878 volIDs->AddAt(4148,iModule); iModule++;
1879 volIDs->AddAt(4149,iModule); iModule++;
1880 volIDs->AddAt(4150,iModule); iModule++;
1881 volIDs->AddAt(4151,iModule); iModule++;
1882 volIDs->AddAt(4152,iModule); iModule++;
1883 volIDs->AddAt(4153,iModule); iModule++;
1884 volIDs->AddAt(4154,iModule); iModule++;
1885 volIDs->AddAt(4155,iModule); iModule++;
1886 volIDs->AddAt(4156,iModule); iModule++;
1887 volIDs->AddAt(4157,iModule); iModule++;
1888 volIDs->AddAt(4158,iModule); iModule++;
1889 volIDs->AddAt(4159,iModule); iModule++;
1891 if(sectors[4]==1){//--->cSect = 4 <//---
1892 volIDs->AddAt(2080,iModule); iModule++;
1893 volIDs->AddAt(2081,iModule); iModule++;
1894 volIDs->AddAt(2082,iModule); iModule++;
1895 volIDs->AddAt(2083,iModule); iModule++;
1896 volIDs->AddAt(2084,iModule); iModule++;
1897 volIDs->AddAt(2085,iModule); iModule++;
1898 volIDs->AddAt(2086,iModule); iModule++;
1899 volIDs->AddAt(2087,iModule); iModule++;
1900 volIDs->AddAt(4160,iModule); iModule++;
1901 volIDs->AddAt(4161,iModule); iModule++;
1902 volIDs->AddAt(4162,iModule); iModule++;
1903 volIDs->AddAt(4163,iModule); iModule++;
1904 volIDs->AddAt(4164,iModule); iModule++;
1905 volIDs->AddAt(4165,iModule); iModule++;
1906 volIDs->AddAt(4166,iModule); iModule++;
1907 volIDs->AddAt(4167,iModule); iModule++;
1908 volIDs->AddAt(4168,iModule); iModule++;
1909 volIDs->AddAt(4169,iModule); iModule++;
1910 volIDs->AddAt(4170,iModule); iModule++;
1911 volIDs->AddAt(4171,iModule); iModule++;
1912 volIDs->AddAt(4172,iModule); iModule++;
1913 volIDs->AddAt(4173,iModule); iModule++;
1914 volIDs->AddAt(4174,iModule); iModule++;
1915 volIDs->AddAt(4175,iModule); iModule++;
1917 if(sectors[5]==1){//--->cSect = 5 <//---
1918 volIDs->AddAt(2088,iModule); iModule++;
1919 volIDs->AddAt(2089,iModule); iModule++;
1920 volIDs->AddAt(2090,iModule); iModule++;
1921 volIDs->AddAt(2091,iModule); iModule++;
1922 volIDs->AddAt(2092,iModule); iModule++;
1923 volIDs->AddAt(2093,iModule); iModule++;
1924 volIDs->AddAt(2094,iModule); iModule++;
1925 volIDs->AddAt(2095,iModule); iModule++;
1926 volIDs->AddAt(4176,iModule); iModule++;
1927 volIDs->AddAt(4177,iModule); iModule++;
1928 volIDs->AddAt(4178,iModule); iModule++;
1929 volIDs->AddAt(4179,iModule); iModule++;
1930 volIDs->AddAt(4180,iModule); iModule++;
1931 volIDs->AddAt(4181,iModule); iModule++;
1932 volIDs->AddAt(4182,iModule); iModule++;
1933 volIDs->AddAt(4183,iModule); iModule++;
1934 volIDs->AddAt(4184,iModule); iModule++;
1935 volIDs->AddAt(4185,iModule); iModule++;
1936 volIDs->AddAt(4186,iModule); iModule++;
1937 volIDs->AddAt(4187,iModule); iModule++;
1938 volIDs->AddAt(4188,iModule); iModule++;
1939 volIDs->AddAt(4189,iModule); iModule++;
1940 volIDs->AddAt(4190,iModule); iModule++;
1941 volIDs->AddAt(4191,iModule); iModule++;
1943 if(sectors[6]==1){//--->cSect = 6 <//---
1944 volIDs->AddAt(2096,iModule); iModule++;
1945 volIDs->AddAt(2097,iModule); iModule++;
1946 volIDs->AddAt(2098,iModule); iModule++;
1947 volIDs->AddAt(2099,iModule); iModule++;
1948 volIDs->AddAt(2100,iModule); iModule++;
1949 volIDs->AddAt(2101,iModule); iModule++;
1950 volIDs->AddAt(2102,iModule); iModule++;
1951 volIDs->AddAt(2103,iModule); iModule++;
1952 volIDs->AddAt(4192,iModule); iModule++;
1953 volIDs->AddAt(4193,iModule); iModule++;
1954 volIDs->AddAt(4194,iModule); iModule++;
1955 volIDs->AddAt(4195,iModule); iModule++;
1956 volIDs->AddAt(4196,iModule); iModule++;
1957 volIDs->AddAt(4197,iModule); iModule++;
1958 volIDs->AddAt(4198,iModule); iModule++;
1959 volIDs->AddAt(4199,iModule); iModule++;
1960 volIDs->AddAt(4200,iModule); iModule++;
1961 volIDs->AddAt(4201,iModule); iModule++;
1962 volIDs->AddAt(4202,iModule); iModule++;
1963 volIDs->AddAt(4203,iModule); iModule++;
1964 volIDs->AddAt(4204,iModule); iModule++;
1965 volIDs->AddAt(4205,iModule); iModule++;
1966 volIDs->AddAt(4206,iModule); iModule++;
1967 volIDs->AddAt(4207,iModule); iModule++;
1969 if(sectors[7]==1){ //--->cSect = 7 <//---
1970 volIDs->AddAt(2104,iModule); iModule++;
1971 volIDs->AddAt(2105,iModule); iModule++;
1972 volIDs->AddAt(2106,iModule); iModule++;
1973 volIDs->AddAt(2107,iModule); iModule++;
1974 volIDs->AddAt(2108,iModule); iModule++;
1975 volIDs->AddAt(2109,iModule); iModule++;
1976 volIDs->AddAt(2110,iModule); iModule++;
1977 volIDs->AddAt(2111,iModule); iModule++;
1978 volIDs->AddAt(4208,iModule); iModule++;
1979 volIDs->AddAt(4209,iModule); iModule++;
1980 volIDs->AddAt(4210,iModule); iModule++;
1981 volIDs->AddAt(4211,iModule); iModule++;
1982 volIDs->AddAt(4212,iModule); iModule++;
1983 volIDs->AddAt(4213,iModule); iModule++;
1984 volIDs->AddAt(4214,iModule); iModule++;
1985 volIDs->AddAt(4215,iModule); iModule++;
1986 volIDs->AddAt(4216,iModule); iModule++;
1987 volIDs->AddAt(4217,iModule); iModule++;
1988 volIDs->AddAt(4218,iModule); iModule++;
1989 volIDs->AddAt(4219,iModule); iModule++;
1990 volIDs->AddAt(4220,iModule); iModule++;
1991 volIDs->AddAt(4221,iModule); iModule++;
1992 volIDs->AddAt(4222,iModule); iModule++;
1993 volIDs->AddAt(4223,iModule); iModule++;
1995 if(sectors[8]==1){//--->cSect = 8 <//---
1996 volIDs->AddAt(2112,iModule); iModule++;
1997 volIDs->AddAt(2113,iModule); iModule++;
1998 volIDs->AddAt(2114,iModule); iModule++;
1999 volIDs->AddAt(2115,iModule); iModule++;
2000 volIDs->AddAt(2116,iModule); iModule++;
2001 volIDs->AddAt(2117,iModule); iModule++;
2002 volIDs->AddAt(2118,iModule); iModule++;
2003 volIDs->AddAt(2119,iModule); iModule++;
2004 volIDs->AddAt(4224,iModule); iModule++;
2005 volIDs->AddAt(4225,iModule); iModule++;
2006 volIDs->AddAt(4226,iModule); iModule++;
2007 volIDs->AddAt(4227,iModule); iModule++;
2008 volIDs->AddAt(4228,iModule); iModule++;
2009 volIDs->AddAt(4229,iModule); iModule++;
2010 volIDs->AddAt(4230,iModule); iModule++;
2011 volIDs->AddAt(4231,iModule); iModule++;
2012 volIDs->AddAt(4232,iModule); iModule++;
2013 volIDs->AddAt(4233,iModule); iModule++;
2014 volIDs->AddAt(4234,iModule); iModule++;
2015 volIDs->AddAt(4235,iModule); iModule++;
2016 volIDs->AddAt(4236,iModule); iModule++;
2017 volIDs->AddAt(4237,iModule); iModule++;
2018 volIDs->AddAt(4238,iModule); iModule++;
2019 volIDs->AddAt(4239,iModule); iModule++;
2021 if(sectors[9]==1){//--->cSect = 9 <//---
2022 volIDs->AddAt(2120,iModule); iModule++;
2023 volIDs->AddAt(2121,iModule); iModule++;
2024 volIDs->AddAt(2122,iModule); iModule++;
2025 volIDs->AddAt(2123,iModule); iModule++;
2026 volIDs->AddAt(2124,iModule); iModule++;
2027 volIDs->AddAt(2125,iModule); iModule++;
2028 volIDs->AddAt(2126,iModule); iModule++;
2029 volIDs->AddAt(2127,iModule); iModule++;
2030 volIDs->AddAt(4240,iModule); iModule++;
2031 volIDs->AddAt(4241,iModule); iModule++;
2032 volIDs->AddAt(4242,iModule); iModule++;
2033 volIDs->AddAt(4243,iModule); iModule++;
2034 volIDs->AddAt(4244,iModule); iModule++;
2035 volIDs->AddAt(4245,iModule); iModule++;
2036 volIDs->AddAt(4246,iModule); iModule++;
2037 volIDs->AddAt(4247,iModule); iModule++;
2038 volIDs->AddAt(4248,iModule); iModule++;
2039 volIDs->AddAt(4249,iModule); iModule++;
2040 volIDs->AddAt(4250,iModule); iModule++;
2041 volIDs->AddAt(4251,iModule); iModule++;
2042 volIDs->AddAt(4252,iModule); iModule++;
2043 volIDs->AddAt(4253,iModule); iModule++;
2044 volIDs->AddAt(4254,iModule); iModule++;
2045 volIDs->AddAt(4255,iModule); iModule++;
2051 //___________________________________
2052 TArrayI* AliITSRealignTracks::GetLayersVolUID(Int_t *layer){
2054 TArrayI *out=new TArrayI(2198);
2057 for(Int_t i=0;i<6;i++){
2059 for(Int_t mod=0;mod<AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);mod++){
2060 voluid=AliGeomManager::LayerToVolUID(i+AliGeomManager::kFirstLayer,mod);
2061 out->AddAt(voluid,last);
2062 // printf("voluid %d at position %d \n",out->At(last),last);
2072 TArrayI* AliITSRealignTracks::SelectLayerInVolids(const TArrayI *volidsIN,AliGeomManager::ELayerID layer){
2074 Int_t size=volidsIN->GetSize();
2076 for(Int_t j=0;j<size;j++){
2077 if(AliGeomManager::VolUIDToLayer(volidsIN->At(j))==layer)count++;
2079 TArrayI *volidsOUT=new TArrayI(count);
2081 for(Int_t j=0;j<size;j++){
2082 if(AliGeomManager::VolUIDToLayer(volidsIN->At(j))==layer){
2083 volidsOUT->AddAt(volidsIN->At(j),count);
2090 //______________________________________________
2092 TArrayI* AliITSRealignTracks::IntersectVolArray(const TArrayI *vol1,const TArrayI *vol2){
2094 Int_t size1=vol1->GetSize();
2095 Int_t size2=vol2->GetSize();
2098 TArrayI *volidOut=new TArrayI(size1+size2);
2100 for(Int_t k=0;k<size1;k++){
2103 for(Int_t j=0;j<size2;j++){
2104 if(vol2->At(j)==volid)found=kTRUE;
2107 volidOut->AddAt(volid,last);
2111 volidOut->Set(last);
2114 //_________________________________________
2116 TArrayI* AliITSRealignTracks::JoinVolArrays(const TArrayI *vol1,const TArrayI *vol2){
2117 //!BE CAREFUL: If an index is repeated in vol1 or vol2 will be repeated also in the final array
2119 Int_t size1=vol1->GetSize();
2120 Int_t size2=vol2->GetSize();
2124 TArrayI *volidOut=new TArrayI(size1+size2);
2126 for(Int_t k=0;k<size1;k++){
2128 volidOut->AddAt(volid,k);
2131 for(Int_t k=0;k<size1;k++){
2134 for(Int_t j=0;j<size2;j++){
2135 if(vol2->At(j)==volid)found=kTRUE;
2138 volidOut->AddAt(volid,size1+count);
2142 volidOut->Set(size1+count);
2146 //______________________________________
2148 TArrayI* AliITSRealignTracks::ExcludeVolidsFromVolidsArray(const TArrayI *volidsToExclude,const TArrayI *volStart){
2150 Int_t size1=volidsToExclude->GetSize();
2151 Int_t size2=volStart->GetSize();
2155 TArrayI *volidOut=new TArrayI(size2);
2157 for(Int_t k=0;k<size2;k++){
2159 volid=volStart->At(k);
2160 for(Int_t j=0;j<size1;j++){
2161 if(volidsToExclude->At(j)==volid){
2167 volidOut->AddAt(volid,last);
2171 volidOut->Set(last);
2176 //________________________________________
2178 TArrayI* AliITSRealignTracks::GetLayerVolumes(Int_t *layer){
2180 TArrayI *out=new TArrayI(2198);
2183 for(Int_t i=0;i<6;i++){
2185 for(Int_t mod=0;mod<AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);mod++){
2186 voluid=AliGeomManager::LayerToVolUID(i+AliGeomManager::kFirstLayer,mod);
2187 out->AddAt(voluid,last);
2188 // printf("voluid %d at position %d \n",out->At(last),last);
2199 //______________________________
2200 TArrayI* AliITSRealignTracks::GetAlignedVolumes(char *filename){
2201 if(gSystem->AccessPathName(filename)){
2202 printf("Wrong Realignment file name \n");
2205 TFile *f=TFile::Open(filename,"READ");
2206 TClonesArray *array=(TClonesArray*)f->Get("ITSAlignObjs");
2207 AliAlignObjParams *a;
2209 TArrayI *volidOut=new TArrayI(2200);
2210 for(Int_t j=0;j<array->GetSize();j++){
2211 a=(AliAlignObjParams*)array->At(j);
2212 if(a->GetUniqueID()==0)continue;
2215 volidOut->AddAt(a->GetVolUID(),last);
2219 volidOut->Set(last);
2225 //________________________________________
2226 void AliITSRealignTracks::SetDraw(Bool_t draw,Bool_t refresh){
2230 if(fAlignDrawObjs)DeleteDrawHists();
2237 void AliITSRealignTracks::DeleteDrawHists(){
2238 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
2239 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
2240 delete fAlignDrawObjs[iLayer][iModule];
2242 if(fAlignDrawObjs[iLayer])delete [] fAlignDrawObjs[iLayer];
2245 delete [] fAlignDrawObjs;
2251 delete fgrIterMeanX;
2253 delete fgrIterMeanY;
2255 delete fgrIterMeanZ;
2257 delete fgrIterMeanPsi;
2258 delete fgrIterRMSPsi;
2259 delete fgrIterMeanTheta;
2260 delete fgrIterRMSTheta;
2261 delete fgrIterMeanPhi;
2262 delete fgrIterRMSPhi;
2266 void AliITSRealignTracks::InitDrawHists(){
2267 Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
2268 fAlignDrawObjs = new AliAlignObj**[nLayers];
2269 for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
2270 fAlignDrawObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
2271 for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
2272 UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
2273 fAlignDrawObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
2274 fAlignDrawObjs[iLayer][iModule]->SetUniqueID(1);
2279 TH1F *hX=new TH1F("hX","hX",1000,-10000.,10000.);
2280 TH1F *hY=new TH1F("hY","hY",1000,-10000.,10000.);
2281 TH1F *hZ=new TH1F("hZ","hZ",1000,-10000.,10000.);
2282 TH1F *hPsi=new TH1F("hPsi","hPsi",1000,-5000.,5000.);
2283 TH1F *hTheta=new TH1F("hTheta","hTheta",1000,-5000.,5000.);
2284 TH1F *hPhi=new TH1F("hPhi","hPhi",1000,-5000.,5000.);
2286 fCanvPar=new TCanvas("fCanvPar","Parameters trend during iterations: Convergence \n");
2287 fCanvPar->Divide(3,2);
2290 hX->SetXTitle("#mum");
2293 hY->SetXTitle("#mum");
2295 hZ->SetXTitle("#mum");
2298 hPsi->SetXTitle("mdeg");
2301 hTheta->SetXTitle("mdeg");
2304 hPhi->SetXTitle("mdeg");
2309 fCanvGr=new TCanvas("fCanvGr","Parameters trend during iterations: Convergence \n");
2310 fCanvGr->Divide(3,2);
2313 fgrIterMeanX=new TGraph(1);
2314 fgrIterRMSX=new TGraph(1);
2315 fgrIterRMSX->GetYaxis()->SetRangeUser(-1000.,1000.);
2316 fgrIterRMSX->SetName("fgrIterRMSX");
2317 fgrIterRMSX->SetLineColor(2);
2318 fgrIterMeanX->SetName("fgrIterMeanX");
2319 fgrIterMeanX->SetTitle("Convergence of #deltaX \n");
2320 fgrIterMeanX->GetXaxis()->SetTitle("#mum");
2321 fgrIterRMSX->Draw("acp");
2322 fgrIterMeanX->Draw("cp");
2325 fgrIterMeanY=new TGraph(1);
2326 fgrIterRMSY=new TGraph(1);
2327 fgrIterRMSY->GetYaxis()->SetRangeUser(-1000.,1000.);
2328 fgrIterRMSY->SetName("fgrIterRMSY");
2329 fgrIterRMSY->SetLineColor(2);
2330 fgrIterMeanY->SetName("fgrIterMeanY");
2331 fgrIterMeanY->SetTitle("Convergence of #deltaY \n");
2332 fgrIterMeanY->GetXaxis()->SetTitle("#mum");
2333 fgrIterRMSY->Draw("acp");
2334 fgrIterMeanY->Draw("cp");
2337 fgrIterMeanZ=new TGraph(1);
2338 fgrIterRMSZ=new TGraph(1);
2339 fgrIterRMSZ->GetYaxis()->SetRangeUser(-1000.,1000.);
2340 fgrIterRMSZ->SetName("fgrIterRMSZ");
2341 fgrIterRMSZ->SetLineColor(2);
2342 fgrIterMeanZ->SetName("fgrIterMeanZ");
2343 fgrIterMeanZ->SetTitle("Convergence of #deltaZ \n");
2344 fgrIterMeanZ->GetXaxis()->SetTitle("#mum");
2345 fgrIterRMSZ->Draw("acp");
2346 fgrIterMeanZ->Draw("cp");
2349 fgrIterMeanPsi=new TGraph(1);
2350 fgrIterRMSPsi=new TGraph(1);
2351 fgrIterRMSPsi->GetYaxis()->SetRangeUser(-1000.,1000.);
2352 fgrIterRMSPsi->SetName("fgrIterRMSPsi");
2353 fgrIterRMSPsi->SetLineColor(2);
2354 fgrIterMeanPsi->SetName("fgrIterMeanPsi");
2355 fgrIterMeanPsi->SetTitle("Convergence of #deltaPsi \n");
2356 fgrIterMeanPsi->GetXaxis()->SetTitle("mdeg");
2357 fgrIterRMSPsi->Draw("acp");
2358 fgrIterMeanPsi->Draw("cp");
2361 fgrIterMeanTheta=new TGraph(1);
2362 fgrIterRMSTheta=new TGraph(1);
2363 fgrIterRMSTheta->GetYaxis()->SetRangeUser(-1000.,1000.);
2364 fgrIterRMSTheta->SetName("fgrIterRMSTheta");
2365 fgrIterRMSTheta->SetLineColor(2);
2366 fgrIterMeanTheta->SetName("fgrIterMeanTheta");
2367 fgrIterMeanTheta->SetTitle("Convergence of #deltaTheta \n");
2368 fgrIterMeanTheta->GetXaxis()->SetTitle("mdeg");
2369 fgrIterRMSTheta->Draw("acp");
2370 fgrIterMeanTheta->Draw("cp");
2373 fgrIterMeanPhi=new TGraph(1);
2374 fgrIterRMSPhi=new TGraph(1);
2375 fgrIterRMSPhi->GetYaxis()->SetRangeUser(-1000.,1000.);
2376 fgrIterRMSPhi->SetName("fgrIterRMSPhi");
2377 fgrIterRMSPhi->SetLineColor(2);
2378 fgrIterMeanPhi->SetName("fgrIterMeanPhi");
2379 fgrIterMeanPhi->SetTitle("Convergence of #deltaPhi \n");
2380 fgrIterMeanPhi->GetXaxis()->SetTitle("mdeg");
2381 fgrIterRMSPhi->Draw("acp");
2382 fgrIterMeanPhi->Draw("cp");
2388 void AliITSRealignTracks::UpdateDraw(TArrayI *volids,Int_t iter,Int_t color){
2393 name.Append("iter");
2394 TH1F *hX=new TH1F("hX",name.Data(),1000,-10000.,10000.);
2398 name.Append("iter");
2399 TH1F *hY=new TH1F("hY",name.Data(),1000,-10000.,10000.);
2403 name.Append("iter");
2404 TH1F *hZ=new TH1F("hZ",name.Data(),1000,-10000.,10000.);
2408 name.Append("iter");
2409 TH1F *hPsi=new TH1F("hPsi",name.Data(),1000,-5000.,5000.);
2413 name.Append("iter");
2414 TH1F *hTheta=new TH1F("hTheta",name.Data(),1000,-5000.,5000.);
2418 name.Append("iter");
2419 TH1F *hPhi=new TH1F("hPhi",name.Data(),1000,-5000.,5000.);
2422 Double_t transl[3],rot[3],transldr[3],rotdr[3];
2424 for(Int_t i=0;i<volids->GetSize();i++){
2425 layer=AliGeomManager::VolUIDToLayer(volids->At(i),mod);
2426 fAlignObjs[layer-AliGeomManager::kFirstLayer][mod]->GetPars(transl,rot);
2427 fAlignDrawObjs[layer-AliGeomManager::kFirstLayer][mod]->GetPars(transldr,rotdr);
2429 hX->Fill(10000.*(transl[0]-transldr[0]));
2430 hY->Fill(10000.*(transl[1]-transldr[1]));
2431 hZ->Fill(10000.*(transl[2]-transldr[2]));
2432 hPsi->Fill(1000.*(rot[0]-rotdr[0]));
2433 hTheta->Fill(1000.*(rot[1]-rotdr[1]));
2434 hPhi->Fill(1000.*(rot[1]-rotdr[2]));
2435 //Update the pars of the draw object
2436 fAlignDrawObjs[layer-AliGeomManager::kFirstLayer][mod]->SetPars(transl[0],transl[1],transl[2],rot[0],rot[1],rot[2]);
2439 hX->SetLineColor(color);
2440 hY->SetLineColor(color);
2441 hZ->SetLineColor(color);
2442 hPsi->SetLineColor(color);
2443 hTheta->SetLineColor(color);
2444 hPhi->SetLineColor(color);
2456 hTheta->Draw("Same");
2461 fCanvPar->Modified();
2463 fgrIterMeanX->SetPoint(fgrIterMeanX->GetN()+1,iter,hX->GetMean());
2464 fgrIterRMSX->SetPoint(fgrIterRMSX->GetN()+1,iter,hX->GetRMS());
2465 fgrIterMeanY->SetPoint(fgrIterMeanY->GetN()+1,iter,hY->GetMean());
2466 fgrIterRMSY->SetPoint(fgrIterRMSY->GetN()+1,iter,hY->GetRMS());
2467 fgrIterMeanZ->SetPoint(fgrIterMeanZ->GetN()+1,iter,hZ->GetMean());
2468 fgrIterRMSZ->SetPoint(fgrIterRMSZ->GetN()+1,iter,hZ->GetRMS());
2469 fgrIterMeanPsi->SetPoint(fgrIterMeanPsi->GetN()+1,iter,hPsi->GetMean());
2470 fgrIterRMSPsi->SetPoint(fgrIterRMSPsi->GetN()+1,iter,hPsi->GetRMS());
2471 fgrIterMeanTheta->SetPoint(fgrIterMeanTheta->GetN()+1,iter,hTheta->GetMean());
2472 fgrIterRMSTheta->SetPoint(fgrIterRMSTheta->GetN()+1,iter,hTheta->GetRMS());
2473 fgrIterMeanPhi->SetPoint(fgrIterMeanPhi->GetN()+1,iter,hPhi->GetMean());
2474 fgrIterRMSPhi->SetPoint(fgrIterRMSPhi->GetN()+1,iter,hPhi->GetRMS());
2481 void AliITSRealignTracks::WriteHists(const char *outfile){
2483 TFile *f=new TFile(outfile,"RECREATE");
2487 fgrIterMeanX->Write();
2488 fgrIterRMSX->Write();
2489 fgrIterMeanY->Write();
2490 fgrIterRMSY->Write();
2491 fgrIterMeanZ->Write();
2492 fgrIterRMSZ->Write();
2493 fgrIterMeanPsi->Write();
2494 fgrIterRMSPsi->Write();
2495 fgrIterMeanTheta->Write();
2496 fgrIterRMSTheta->Write();
2497 fgrIterMeanPhi->Write();
2498 fgrIterRMSPhi->Write();