ITS/AliITSRealignTracks.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 //Class to perform the realignment if the Inner Tracking System
  18 //with an iterative approach based on track to cluster residuals
  19 // minimization. A chi2 function of the residuals is minimized with
  20 //respect to alignment parameters. The class allows both single module
  21 //realignment and set of modules realignment. Tracks are fitted with
  22 //AliTrackFitter* fitters. AliTrackFitterKalman is more suited for
  23 //straight line (e.g. cosmic in the absence of magnetic field) but can't
  24 //work with helixes. AliTrackFitterRieman is suited for helixes.
  25 //The minimization is performed by AliTrackResiduals* classes: default
  26 //one is AliTrackResidualsFast (analytic minimization by inversion).
  27 //For numerical minimization using MINUIT, use AliTrackResidualChi2.
  28 //Methods are present to defined both the set of modules where the tracks
  29 //are fittef and the set of modules which  are to be realigned
  30 //The main method is  AlignVolumesITS
  31 //
  32 //Class by: A. Rossi, andrea,rossi@ts.infn.it
  33
  34 #include <TFile.h>
  35 #include <TStopwatch.h>
  36 #include <TNtuple.h>
  37 #include <TClonesArray.h>
  38 #include <TMath.h>
  39 #include <TGraph.h>
  40 #include <TCanvas.h>
  41 #include <TH1F.h>
  42 #include "AliITSRealignTracks.h"
  43 #include "AliAlignObjParams.h"
  44 #include "AliAlignObj.h"
  45 #include "AliGeomManager.h"
  46 #include "AliTrackFitter.h"
  47 #include "AliTrackFitterKalman.h"
  48 #include "AliTrackFitterRieman.h"
  49 #include "AliTrackResidualsFast.h"
  50 #include "AliTrackResidualsChi2.h"
  51 #include "AliTrackResidualsLinear.h"
  52
  53 #include "AliLog.h"
  54 #include <TSystem.h>
  55 #include <TGeoManager.h>
  56
  57 class AliAlignmentTracks;
  58 class TGeoMatrix;
  59 class TArray;
  60
  61
  62 /* $Id$ */
  63
  64
  65 ClassImp(AliITSRealignTracks)
  66
  67 const Int_t kreferSect=2;
  68
  69
  70 AliITSRealignTracks::AliITSRealignTracks(TString minimizer,Int_t fit,Bool_t covUsed,TString fileintro,TString geometryfile,TString misalignmentFile,TString startingfile):
  71   AliAlignmentTracks(),
  72   fSurveyObjs(0),
  73   fgeomfilename(),
  74   fmintracks(),
  75   fUpdateCov(kFALSE),
  76   fVarySigmaY(kFALSE),
  77   fCorrModules(0),
  78   fLimitCorr(0),
  79   fsigmaY(),
  80   fDraw(kFALSE),
  81   fAlignDrawObjs(0),
  82   fCanvPar(0),
  83   fCanvGr(0),
  84   fgrIterMeanX(0),
  85   fgrIterRMSX(0),
  86   fgrIterMeanY(0),
  87   fgrIterRMSY(0),
  88   fgrIterMeanZ(0),
  89   fgrIterRMSZ(0),
  90   fgrIterMeanPsi(0),
  91   fgrIterRMSPsi(0),
  92   fgrIterMeanTheta(0),
  93   fgrIterRMSTheta(0),
  94   fgrIterMeanPhi(0),
  95   fgrIterRMSPhi(0)
  96 {
  97
  98   // minimizer="fast"->AliTrackResidualFast minimizer
  99   //           "minuit"->AliTrackResidualChi2 minimizer
 100   //           "minuitnorot"->AliTrackResidualChi2 minimizer without rotations degrees of freedom
 101   //           "linear"->AliTrackResidualLinear minimizer
 102   //fit=0-> Riemann Fitter, fit=1->Kalman
 103   //fileintro=file into which the Tree with the space points is stored
 104   //geometryfile=file containing the geometry
 105
 106
 107   SetPointsFilename(fileintro.Data());
 108   SetGeomFilename(geometryfile);
 109   InitAlignObjs();
 110   if(!InitSurveyObjs(kFALSE))AliWarning("Unable to set Survey AlignObjs!");
 111
 112   if(startingfile=="")printf("Starting from default geometry \n");
 113   else ReadAlignObjs(startingfile.Data(),"ITSAlignObjs");
 114
 115   if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT INTRODUCED \n");
 116   else {
 117     Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
 118     if(!misal)AliWarning("Incorrect fake misalignment filename!");;
 119   }
 120
 121   if(!gGeoManager) AliGeomManager::LoadGeometry(fgeomfilename.Data());
 122   if(covUsed)SetCovIsUsed(kTRUE);
 123   if(!SelectFitter(fit))AliWarning("Incorrect fitter assignment!");
 124   if(!SelectMinimizer(minimizer))AliWarning("Incorrect minimizer assignment!");
 125   fsigmaY=1.;
 126   fmintracks=1;
 127   BuildIndex();
 128
 129 }
 130
 131 AliITSRealignTracks::AliITSRealignTracks(const AliITSRealignTracks &realignTracks):
 132   AliAlignmentTracks(),
 133   fSurveyObjs(new AliAlignObj**(*realignTracks.fSurveyObjs)),
 134   fgeomfilename(realignTracks.fgeomfilename),
 135   fmintracks(realignTracks.fmintracks),
 136   fUpdateCov(realignTracks.fUpdateCov),
 137   fVarySigmaY(realignTracks.fVarySigmaY),
 138   fCorrModules(new Double_t *(*realignTracks.fCorrModules)),
 139   fLimitCorr(realignTracks.fLimitCorr),
 140   fsigmaY(realignTracks.fsigmaY),
 141   fDraw(kFALSE),
 142   fAlignDrawObjs(realignTracks.fAlignDrawObjs),
 143   fCanvPar(realignTracks.fCanvPar),
 144   fCanvGr(realignTracks.fCanvGr),
 145   fgrIterMeanX(realignTracks.fgrIterMeanX),
 146   fgrIterRMSX(realignTracks.fgrIterRMSX),
 147   fgrIterMeanY(realignTracks.fgrIterMeanY),
 148   fgrIterRMSY(realignTracks.fgrIterRMSY),
 149   fgrIterMeanZ(realignTracks.fgrIterMeanZ),
 150   fgrIterRMSZ(realignTracks.fgrIterRMSZ),
 151   fgrIterMeanPsi(realignTracks.fgrIterMeanPsi),
 152   fgrIterRMSPsi(realignTracks.fgrIterRMSPsi),
 153   fgrIterMeanTheta(realignTracks.fgrIterMeanTheta),
 154   fgrIterRMSTheta(realignTracks.fgrIterRMSTheta),
 155   fgrIterMeanPhi(realignTracks.fgrIterMeanPhi),
 156   fgrIterRMSPhi(realignTracks.fgrIterRMSPhi)
 157
 158 {//Copy Constructor
 159   AliWarning("Can't copy AliAlignmentTracks Data member!");
 160 }
 161
 162 AliITSRealignTracks& AliITSRealignTracks::operator=(const AliITSRealignTracks &obj){
 163   ////////////////////////
 164   // Assignment operator
 165   ////////////////////////
 166   this->~AliITSRealignTracks();
 167   new(this) AliITSRealignTracks(obj);
 168   return *this;
 169 }
 170
 171 AliITSRealignTracks::~AliITSRealignTracks(){
 172   //destructor
 173
 174   if(fSurveyObjs)   DeleteSurveyObjs();
 175   if(fAlignDrawObjs) DeleteDrawHists();
 176   //delete [] fSurveyObjs;
 177
 178 }
 179
 180
 181 //_____________________________
 182 Bool_t AliITSRealignTracks::SelectFitter(Int_t fit,Int_t minTrackPoint){
 183   //Method to select the fitter: 0 for AliTrackFitterRieman (use this for helixes)
 184   //                             1 for AliTrackFitterKalman
 185   //minTrackPoint defines the minimum number of points (not rejected by the fit itself)
 186   //a track should have to fit it
 187
 188   if(fit==1){
 189      AliTrackFitterKalman *fitter= new AliTrackFitterKalman();
 190      fitter->SetMinNPoints(minTrackPoint);
 191      SetTrackFitter(fitter);
 192   }
 193
 194   else if(fit==0){
 195     AliTrackFitterRieman *fitter=new AliTrackFitterRieman();
 196     fitter->SetMinNPoints(minTrackPoint);
 197     SetTrackFitter(fitter);
 198   }
 199   else return kFALSE;
 200
 201   return kTRUE;
 202 }
 203
 204
 205 Bool_t AliITSRealignTracks::SelectMinimizer(TString minimizer,Int_t minpoints,const Bool_t *coord){
 206   //Method to select the minimizer: "minuit" for AliTrackFitterChi2 (numerical minimization by MINUIT)
 207   //                                "fast" for AliTrackResidualsFast
 208   //                                "linear" for AliTrackResidualsLinear
 209   //                                "minuitnorot" for AliTrackFitterChi2 by
 210   // coord[6] allows to fix the degrees of freedom in the minimization (e.g. look only for tranlsations,
 211   //       or only for rotations). The coord are: dTx,dTy,dTz,dPsi,dTheta,dPhi where "d" stands for
 212   //       "differential" and "T" for translation. If coord[i] is set to kTRUE then the i coord is fixed
 213   //       When a coordinate is fixed the value returnd for it is 0
 214   // minnpoints fix the minimum number of residuals to perform the minimization: it's not safe
 215   //     to align a module with a small number of track passing through it since the results could be
 216   //     not reliable. For single modules the number of residuals and the number of tracks passing
 217   //     through it and accepted bu the fit procedure coincide. This is not the case for sets of modules,
 218   //     since a given track can pass through two or more modules in the set (e.g a cosmic track can give
 219   //     two cluster on a layer)
 220
 221   AliTrackResiduals *res;
 222   if(minimizer=="minuit"){
 223     res = new AliTrackResidualsChi2();
 224     if(coord){
 225       for(Int_t j=0;j<6;j++){
 226         if(coord[j])res->FixParameter(j);
 227       }
 228     }
 229   }
 230   else if(minimizer=="minuitnorot"){
 231     res = new AliTrackResidualsChi2();
 232     res->FixParameter(3);
 233     res->FixParameter(4);
 234     res->FixParameter(5);
 235   }
 236
 237   else if(minimizer=="fast"){
 238     res = new AliTrackResidualsFast();
 239     if(coord){
 240       for(Int_t j=0;j<6;j++){
 241         if(coord[j])res->FixParameter(j);
 242       }
 243     }
 244   }
 245   else if(minimizer=="linear"){
 246     res = new AliTrackResidualsLinear();
 247   }
 248
 249   else {
 250     printf("Trying to set a non existing minimizer! \n");
 251     return kFALSE;
 252   }
 253
 254   res->SetMinNPoints(minpoints);
 255   SetMinimizer(res);
 256
 257   return kTRUE;
 258 }
 259
 260 //____________________________________
 261 void AliITSRealignTracks::SetVarySigmaY(Bool_t varysigmay,Double_t sigmaYfixed){
 262   //SigmaY is the value of the error along the track direction assigned
 263   //to the AliTrackPoint constructed from the extrapolation of the fit of a track
 264   //to the module one is realigning. This error simulate the uncertainty on the
 265   //position of the cluster in space due to the misalingment of the module (i.e.
 266   //you don't the real position and orientation of the plane of the desired extrapolation
 267   //but you rely on the fit and consider the point to lie along the track)
 268
 269
 270   fVarySigmaY=varysigmay;
 271   if(!varysigmay){
 272     if(sigmaYfixed>0.)fsigmaY=sigmaYfixed;
 273     else {
 274       printf("Negative assignment to sigmaY! set it to default value (1 cm) \n");
 275       fsigmaY=1.;
 276     }
 277   }
 278 }
 279
 280 //_______________________________________
 281 void AliITSRealignTracks::RealignITSVolIndependent(Int_t iter1,Int_t iterations,Int_t minNtracks,Int_t layer,Int_t minTrackPoint){
 282
 283
 284   //iter1=#iterations inside AliAlignmentTracks::AliAlignVolumesITS method
 285   //iterations=#iterations on the all procedure
 286   //layer=0->all ITS, otherways the usual notation is considered (1=SPD1,2=SPD2,3=SDD1,4=SDD2,5=SSD1,6=SSD2)
 287   //minNtracks=minimun number of tracks passing through a module in order to try to realign the module itsself
 288   //           if minNtracks<0, minimun number of tracks is |minNtracks|*minNumPoint[layer]/fact (see the code below): this allows a different
 289   //           choice of the number of tracks required on different layers and to vary these numbers once tuned the relative proportions.
 290   //minTrackPoint=minimun number of "good" points required to a track (THE POINT ON THE MODULE THAT IS GOING TO BE REALIGNED
 291   //IS NEVER CONSIDERED->max number can be required is 11 for cosmics tracks) for the track being considered in the minimization
 292
 293
 294   fTrackFitter->SetMinNPoints(minTrackPoint);
 295   TArrayI volIDs2(2200);
 296   volIDs2.Reset(0);
 297   TArrayI volIDs(1);
 298   TString command;
 299   TArrayI volIDsFit;
 300
 301   Int_t iLayer,iLayerToAlign;
 302
 303   Int_t minNumPoint[6]={100,100,100,100,50,50};
 304   Double_t fact=10;
 305   Int_t j=0;
 306
 307   Int_t size=0;
 308   Int_t layerNum,modNum,lastVolid=0;
 309   TNtuple *ntVolumeAlign=new TNtuple("ntVolumeAlign","NTuple with volume tried to be realigned","layerNum:modNum:volumeIDnum");
 310
 311   TStopwatch *timer=new TStopwatch();
 312   timer->Start();
 313   BuildIndex();
 314   j=0;
 315   UShort_t volid;
 316
 317   for(Int_t iter=0;iter<iterations;iter++){
 318
 319     //Starting Independent Modules Realignment
 320     for(iLayerToAlign=(Int_t)AliGeomManager::kSPD1;iLayerToAlign<=(Int_t)AliGeomManager::kSSD2;iLayerToAlign++){
 321       if(layer!=0&&iLayerToAlign!=layer)continue;
 322       j=0;
 323       size=0;
 324       for(Int_t k=(Int_t)AliGeomManager::kSPD1;k<=(Int_t)AliGeomManager::kSSD2;k++){
 325         size+=AliGeomManager::LayerSize(k);
 326         printf("size: %d \n",size);
 327       }
 328
 329       for (Int_t iModule=0;iModule<AliGeomManager::LayerSize(iLayerToAlign);iModule++){
 330         j=0;
 331         if(minNtracks<0){
 332           if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNumPoint[iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer]*(-1*minNtracks/fact))continue;        }
 333         else if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNtracks)continue;
 334
 335         UShort_t volidAl = AliGeomManager::LayerToVolUID(iLayerToAlign,iModule);
 336
 337         volIDsFit.Reset(0);
 338         volIDsFit.Set(size-1);
 339         for (iLayer=AliGeomManager::kSPD1;iLayer<AliGeomManager::kTPC1;iLayer++){
 340           for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
 341             volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
 342             if(AliGeomManager::LayerToVolUID(iLayer,iModule2)==volidAl)continue;
 343             volIDsFit.AddAt(volid,j);
 344             j++;
 345           }
 346         }
 347         volIDs.AddAt((Int_t)volidAl,0);
 348         if(iter==iterations-1){
 349           volIDs2.AddAt(volidAl,lastVolid);
 350           lastVolid++;
 351         }
 352         volIDs2.AddAt(volidAl,lastVolid);
 353         AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iter1);
 354       }
 355     }
 356
 357
 358     if((iter+1)%5==0||iter==0||iter==1||iter==2||iter==3||iter==iterations-1){
 359       command="RealignObj";
 360       command+=(iter+1);
 361       command.Append(".root");
 362       WriteRealignObjArray(command.Data(),AliGeomManager::kSPD1,AliGeomManager::kSSD2);
 363     }
 364   }
 365
 366
 367   if(j==0){printf("j=0 \n");return;}
 368   for(Int_t k=0;k<volIDs2.GetSize();k++){
 369     if(volIDs2.At(k)==0)break;
 370     layerNum=AliGeomManager::VolUIDToLayer(volIDs2.At(k),modNum);
 371     ntVolumeAlign->Fill(layerNum,modNum,volIDs2.At(k));
 372   }
 373   printf("End of selecting modules cycle: %d modules selected \n",j);
 374   TFile *f=new TFile("RealignVolNt.root","RECREATE");
 375   f->cd();
 376   ntVolumeAlign->Write();
 377   f->Close();
 378
 379   timer->Stop();
 380   timer->Print();
 381   return;
 382 }
 383
 384
 385 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)
 386 {
 387
 388   // minimizer="fast"->AliTrackResidualFast minimizer
 389   //           "minuit"->AliTrackResidualChi2 minimizer
 390   //           "minuitnorot"->AliTrackResidualChi2 minimizer without rotations degrees of freedom
 391   //           "linear"->AliTrackResidualLinear minimizer
 392   //fit=0-> Riemann Fitter, fit=1->Kalman
 393   //iter1=#iterations inside AliAlignmentTracks::AliAlignVolumesITS method
 394   //iterations=#iterations on the all procedure
 395   //layer=0->all ITS, otherways the usual notation is considered (1=SPD1,2=SPD2,3=SDD1,4=SDD2,5=SSD1,6=SSD2)
 396   //minNtracks=minimun number of tracks passing through a module in order to try to realign the module itsself
 397   //           if minNtracks<0, minimun number of tracks is |minNtracks|*minNumPoint[layer]/fact (see the code below): this allows a different
 398   //           choice of the number of tracks required on different layers and to vary these numbers once tuned the relative proportions.
 399   //minTrackPoint=minimun number of "good" points required to a track (THE POINT ON THE MODULE THAT IS GOING TO BE REALIGNED
 400   //              IS NEVER CONSIDERED->max number that can be required is 11 for cosmics tracks) for the track being considered in the minimization
 401   //doGlobal : do global realignment, 0=no, 1= yes, 2=only global
 402
 403
 404   TArrayI volIDs2(2200);
 405   volIDs2.Reset(0);
 406   TArrayI volIDs(1);
 407   TString command;
 408   TArrayI volIDsFit;
 409
 410   Int_t iLayer,iLayerToAlign;
 411
 412   Int_t minNumPoint[6]={100,100,100,100,50,50};
 413   Double_t fact=10;
 414   Int_t count=0;
 415
 416   Int_t size=0;
 417   Int_t layerNum,modNum,lastVolid=0;
 418   TNtuple *ntVolumeAlign=new TNtuple("ntVolumeAlign","NTuple with volume tried to be realigned","layerNum:modNum:volumeIDnum");
 419
 420
 421   if(!SelectFitter(fit))AliWarning("Incorrect fitter assignment!");
 422   if(!SelectMinimizer(minimizer))AliWarning("Incorrect minimizer assignment!");
 423   if(misalignmentFile=="")printf("NO FAKE MISALIGNMENT INTRODUCED \n");
 424   else {
 425     Bool_t misal=Misalign(misalignmentFile,"ITSAlignObjs");
 426     if(!misal)return;
 427   }
 428
 429
 430   TStopwatch *timer=new TStopwatch();
 431   timer->Start();
 432   BuildIndex();
 433   count=0;
 434   UShort_t volid;
 435
 436   if(startingfile=="")printf("Starting from default geometry \n");
 437   else {
 438     printf("Starting from AlignObjs file: %s",startingfile.Data());
 439     ReadAlignObjs(startingfile.Data(),"ITSAlignObjs");
 440   }
 441
 442   for(Int_t iter=0;iter<iterations;iter++){
 443     if(covUsed)SetCovIsUsed(kTRUE);
 444
 445
 446     //START HIERARCHY REALIGNMENT
 447
 448     if(layer==0&&(doGlobal==1||doGlobal==2)){
 449       for(Int_t siter=0;siter<5;siter++){
 450         fTrackFitter->SetMinNPoints(2);
 451         SetCovUpdate(kFALSE);
 452         AlignSPDHalfBarrelToSectorRef(kreferSect,3);
 453         //      AlignSPDBarrel(1);
 454         //      if(siter==0)SetCovUpdate(kFALSE);
 455         //      AlignSPDHalfBarrel(0,3);
 456         //      SetCovUpdate(kTRUE);
 457         AlignSPDHalfBarrelToHalfBarrel(1,3);
 458         //      AlignSPDHalfBarrelToSectorRef(kreferSect,3);
 459         for(Int_t sector=0;sector<10;sector++){
 460           SetMinNtracks(100);
 461           if(sector==kreferSect)continue;
 462           AlignSPDSectorWithSectors(sector,1);
 463         }
 464
 465
 466         for(Int_t lay=1;lay<=6;lay++){
 467           if(!AlignLayerToSector(lay,kreferSect,3))AlignLayerToSPDHalfBarrel(lay,0,3);
 468         }
 469         AlignSPDHalfBarrel(0,3);
 470
 471         Int_t layers[6]={2,2,1,0,0,0};
 472         fTrackFitter->SetMinNPoints(4);
 473         AlignLayersToLayers(layers,1);
 474
 475         fTrackFitter->SetMinNPoints(6);
 476         layers[2]=2;
 477         layers[3]=1;//{2,2,2,1,0,0};
 478         AlignLayersToLayers(layers,1);
 479
 480         fTrackFitter->SetMinNPoints(6);
 481         layers[3]=2;
 482         layers[4]=1;//{2,2,2,2,1,0};
 483         AlignLayersToLayers(layers,1);
 484
 485         fTrackFitter->SetMinNPoints(6);
 486         layers[4]=2;
 487         layers[5]=1;//{2,2,2,2,2,1};
 488         AlignLayersToLayers(layers,1);
 489
 490
 491         for(Int_t sector=0;sector<10;sector++){
 492           AlignSPDSectorToOuterLayers(sector,1);
 493         }
 494         WriteRealignObjArray("AfterGlobal.root",AliGeomManager::kSPD1,AliGeomManager::kSSD2);
 495       }
 496     }
 497
 498     if(doGlobal==2)return;
 499
 500     if(covUsed)SetCovUpdate(kTRUE);
 501     SetMinNtracks(1);
 502
 503
 504     // STARTS INDEPENDENT MOULES REALIGNMENT
 505
 506     fTrackFitter->SetMinNPoints(minTrackPoint);
 507     for(iLayerToAlign=(Int_t)AliGeomManager::kSPD1;iLayerToAlign<=(Int_t)AliGeomManager::kSSD2;iLayerToAlign++){
 508       if(layer!=0&&iLayerToAlign!=layer)continue;
 509       count=0;
 510       size=0;
 511       for(Int_t k=(Int_t)AliGeomManager::kSPD1;k<=(Int_t)AliGeomManager::kSSD2;k++){
 512         size+=AliGeomManager::LayerSize(k);
 513         printf("size: %d \n",size);
 514       }
 515
 516       for (Int_t iModule=0;iModule<AliGeomManager::LayerSize(iLayerToAlign);iModule++){
 517         count=0;
 518         if(minNtracks<0){
 519           if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNumPoint[iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer]*(-1*minNtracks/fact))continue;
 520         }
 521         else if(GetLastIndex(iLayerToAlign-(Int_t)AliGeomManager::kFirstLayer,iModule)<minNtracks)continue;
 522
 523         UShort_t volidAl = AliGeomManager::LayerToVolUID(iLayerToAlign,iModule);
 524
 525         volIDsFit.Reset(0);
 526         volIDsFit.Set(size-1);
 527         for (iLayer=(Int_t)AliGeomManager::kSPD1;iLayer<(Int_t)AliGeomManager::kTPC1;iLayer++){
 528           for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
 529
 530
 531             volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
 532
 533             if(AliGeomManager::LayerToVolUID(iLayer,iModule2)==volidAl)continue;
 534             volIDsFit.AddAt(volid,count);
 535             count++;
 536           }
 537         }
 538
 539         volIDs.AddAt((Int_t)volidAl,0);
 540         if(iter==iterations-1){
 541           volIDs2.AddAt(volidAl,lastVolid);
 542           lastVolid++;
 543         }
 544         volIDs2.AddAt(volidAl,lastVolid);
 545         AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iter1);
 546       }
 547     }
 548
 549     if((iter+1)%2==0||(iter+1)%5==0||iter==0||iter==1||iter==2||iter==3||iter==iterations-1){
 550       command="RealignObj";
 551       command+=(iter+1);
 552       command.Append(".root");
 553       WriteRealignObjArray(command.Data(),AliGeomManager::kSPD1,AliGeomManager::kSSD2);
 554     }
 555
 556   }
 557
 558   if(count==0){printf("count=0 \n");return;}
 559   for(Int_t k=0;k<volIDs2.GetSize();k++){
 560     if(volIDs2.At(k)==0)break;
 561     layerNum=AliGeomManager::VolUIDToLayer(volIDs2.At(k),modNum);
 562     ntVolumeAlign->Fill(layerNum,modNum,volIDs2.At(k));
 563   }
 564   printf("End of selecting modules cycle: %d modules selected \n",count);
 565   TFile *f=new TFile("RealignVolNt.root","RECREATE");
 566   f->cd();
 567   ntVolumeAlign->Write();
 568   f->Close();
 569
 570   timer->Stop();
 571   timer->Print();
 572   return;
 573
 574 }
 575
 576
 577 //______________________________________________________________________________
 578 void AliITSRealignTracks::InitAlignObjs()
 579 {
 580   // Initialize the alignment objects array
 581   TMatrixDSym c(6);
 582   Double_t cov[21];
 583   for(Int_t i=0;i<21;i++)cov[i]=0.;
 584   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
 585   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)
 586
 587   Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
 588   fAlignObjs = new AliAlignObj**[nLayers];
 589   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 590     fAlignObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
 591     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 592       UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
 593       fAlignObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
 594       fAlignObjs[iLayer][iModule]->SetCorrMatrix(cov);
 595       fAlignObjs[iLayer][iModule]->SetUniqueID(0);
 596     }
 597   }
 598 }
 599
 600 //______________________________________________________________________________
 601 void AliITSRealignTracks::ResetCorrModules(){
 602   // Initialize and reset to 0 the array with the information on correlation
 603   if(!fCorrModules){
 604     Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
 605     fCorrModules = new Double_t*[nLayers];
 606     for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 607       fCorrModules[iLayer] = new Double_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
 608       for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 609         fCorrModules[iLayer][iModule]=0.;
 610       }
 611     }
 612   }
 613   else{
 614     for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 615       for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 616         fCorrModules[iLayer][iModule]=0.;
 617       }
 618     }
 619   }
 620 }
 621
 622 //______________________________________________________________________________
 623 Bool_t AliITSRealignTracks::InitSurveyObjs(Bool_t infinite,Double_t factor,TString filename,TString arrayName){
 624   //Initialize the Survey Objects. There is the possibility to set them equal to external objects
 625   //   stored in file filename. Otherwuse they are set equal to 0 and with default values for the variances
 626   //infinite: set the cov matrix to extremly large values, so that the results of a minimization
 627   //   are never rejected by the comparison with the survey
 628   //factor: multiplication factor for the variances of the cov. matrix of the survey obj.
 629
 630   if(fSurveyObjs)DeleteSurveyObjs();
 631   Bool_t fromfile=kFALSE;
 632   TFile *surveyObj;
 633   TClonesArray *clnarray;
 634   if(!filename.IsNull()){
 635     //Initialize from file
 636     if(gSystem->AccessPathName(filename.Data(),kFileExists)){
 637       printf("Wrong Survey AlignObjs File Name \n");
 638       return kFALSE;
 639     }
 640     if(arrayName.IsNull()){
 641       printf("Null Survey Object Name! \n");
 642       return kFALSE;
 643     }
 644
 645     fromfile=kTRUE;
 646   }
 647
 648   // Initialize the alignment objects array with default values
 649   Double_t v=1.*factor;
 650   if(infinite)v*=100000.;
 651   TMatrixDSym c(6);
 652   Double_t cov[21];
 653   for(Int_t i=0;i<21;i++)cov[i]=0.;
 654   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
 655   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)
 656
 657   Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
 658   fSurveyObjs = new AliAlignObj**[nLayers];
 659   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 660     fSurveyObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
 661     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 662       fSurveyObjs[iLayer][iModule] = 0x0;
 663     }
 664   }
 665
 666   if(fromfile){
 667     surveyObj=TFile::Open(filename.Data());
 668     if (!surveyObj || !surveyObj->IsOpen()) {
 669       AliError(Form("Could not open SurveyObjs file: %s !",filename.Data()));
 670       return kFALSE;
 671     }
 672     printf("Getting TClonesArray \n");
 673     clnarray=(TClonesArray*)surveyObj->Get(arrayName);
 674     Int_t size=clnarray->GetSize();
 675     UShort_t volid;
 676     for(Int_t ivol=0;ivol<size;ivol++){
 677       AliAlignObjParams *a=(AliAlignObjParams*)clnarray->At(ivol);
 678       volid=a->GetVolUID();
 679       Int_t iModule;
 680       AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
 681       if(iLayer<=0)continue;
 682       if(a->GetUniqueID()==0)continue;
 683       printf("Updating survey for volume: %d ,layer: %d module: %d from file\n",volid,iLayer,iModule);
 684       fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule] = new AliAlignObjParams(*a);
 685       fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->SetUniqueID(a->GetUniqueID());
 686       fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->Print("");
 687     }
 688     delete clnarray;
 689     surveyObj->Close();
 690   }
 691
 692   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 693     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 694       UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
 695       if(!fSurveyObjs[iLayer][iModule]){
 696         printf("Updating survey for volume: %d ,layer: %d module: %d with default values \n",volid,iLayer,iModule);
 697         fSurveyObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
 698         fSurveyObjs[iLayer][iModule]->SetCorrMatrix(cov);
 699         fSurveyObjs[iLayer][iModule]->SetUniqueID(0);
 700       }
 701
 702     }
 703   }
 704
 705
 706   return kTRUE;
 707 }
 708
 709
 710 //______________________________________________________________________________
 711 Int_t AliITSRealignTracks::CheckWithSurvey(Double_t factor,const TArrayI *volids){
 712
 713   // Check the parameters of the alignment objects in volids (or of all objects if volids is null)
 714   // are into the boundaries set by the cov. matrix of the survey objs
 715   // Returns the number of objects out of boudaries
 716   AliAlignObj *alignObj;
 717   Int_t outofsurv=0;
 718   UShort_t volid;
 719   Double_t surveycov[21],transl[3],rot[3],survtransl[3],survrot[3];
 720   if(volids==0x0){
 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         volid=AliGeomManager::LayerToVolUIDSafe(iLayer+AliGeomManager::kFirstLayer,iModule);
 724         alignObj=GetAlignObj(volid);
 725         alignObj->GetPars(transl,rot);
 726         fSurveyObjs[iLayer][iModule]->GetCovMatrix(surveycov);
 727         fSurveyObjs[iLayer][iModule]->GetPars(survtransl,survrot);
 728         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])){
 729           printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
 730           //      *alignObj = *alignObjSurv;
 731           alignObj->SetPars(survtransl[0],survtransl[1],survtransl[2],survrot[0],survrot[1],survrot[2]);
 732           alignObj->SetUniqueID(0);
 733           if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
 734           outofsurv++;
 735         }
 736       }
 737     }
 738   }
 739   else{
 740     Int_t iLayer;
 741     Int_t iModule;
 742     for(Int_t j=0;j<volids->GetSize();j++){
 743       volid=volids->At(j);
 744       alignObj=GetAlignObj(volid);
 745       alignObj->GetPars(transl,rot);
 746       iLayer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volid,iModule)-(Int_t)AliGeomManager::kFirstLayer;
 747       fSurveyObjs[iLayer][iModule]->GetCovMatrix(surveycov);
 748       fSurveyObjs[iLayer][iModule]->GetPars(survtransl,survrot);
 749       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])){
 750         printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
 751         //        *alignObj = *alignObjSurv;
 752         alignObj->SetPars(survtransl[0],survtransl[1],survtransl[2],survrot[0],survrot[1],survrot[2]);
 753         alignObj->SetUniqueID(0);
 754         if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
 755         outofsurv++;
 756       }
 757     }
 758   }
 759   return outofsurv;
 760 }
 761
 762 //___________________________________________________________________
 763
 764 void AliITSRealignTracks::ResetAlignObjs(Bool_t all,TArrayI *volids)
 765 {
 766   // Reset the alignment objects in volids or all if all=kTRUE
 767   if(all){
 768     for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 769       for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++)
 770         fAlignObjs[iLayer][iModule]->SetPars(0,0,0,0,0,0);
 771     }
 772   }
 773   else{
 774     Int_t layer;
 775     Int_t mod;
 776     for(Int_t j=0;j<volids->GetSize();j++){
 777       layer=(Int_t)AliGeomManager::VolUIDToLayer(volids->At(j),mod)-(Int_t)AliGeomManager::kFirstLayer;
 778       fAlignObjs[layer][mod]->SetPars(0,0,0,0,0,0);
 779     }
 780   }
 781 }
 782
 783 //______________________________________________-
 784 void AliITSRealignTracks::DeleteSurveyObjs()
 785 {//destructor for the survey objs. array
 786
 787   if(!fSurveyObjs)return;
 788   // Delete the alignment objects array
 789   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 790     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++){
 791       if (fSurveyObjs[iLayer][iModule]) delete fSurveyObjs[iLayer][iModule];
 792     }
 793
 794     if(fSurveyObjs[iLayer])delete [] fSurveyObjs[iLayer];
 795   }
 796
 797   delete [] fSurveyObjs;
 798   fSurveyObjs = 0;
 799 }
 800
 801
 802 //______________________________________________________________________________
 803 Bool_t AliITSRealignTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName){
 804
 805   // Read alignment object from a file: update the alignobj already present with the one in the file
 806   // To be replaced by a call to CDB
 807
 808   if(gSystem->AccessPathName(alignObjFileName,kFileExists)){
 809     printf("Wrong AlignObjs File Name \n");
 810     return kFALSE;
 811   }
 812
 813   TFile *fRealign=TFile::Open(alignObjFileName);
 814   if (!fRealign || !fRealign->IsOpen()) {
 815     AliError(Form("Could not open Align Obj File file %s !",alignObjFileName));
 816     return kFALSE;
 817   }
 818   printf("Getting TClonesArray \n");
 819   TClonesArray *clnarray=(TClonesArray*)fRealign->Get(arrayName);
 820   Int_t size=clnarray->GetSize();
 821   UShort_t volid;
 822
 823   for(Int_t ivol=0;ivol<size;ivol++){
 824     AliAlignObjParams *a=(AliAlignObjParams*)clnarray->At(ivol);
 825     volid=a->GetVolUID();
 826     Int_t iModule;
 827     AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
 828     if(iLayer<AliGeomManager::kFirstLayer||iLayer>AliGeomManager::kSSD2)continue;
 829     printf("Updating volume: %d ,layer: %d module: %d \n",volid,iLayer,iModule);
 830     *fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule] *= *a;
 831   }
 832
 833   delete clnarray;
 834   fRealign->Close();
 835   return kTRUE;
 836 }
 837
 838 //_________________________________________
 839 Bool_t AliITSRealignTracks::FirstAlignmentLayers(const Bool_t *layers,Int_t minNtracks,Int_t iterations,Bool_t fitall,const TArrayI *volidsSet){
 840
 841   //Align all modules in the set of layers independently according to a sequence based on the number of tracks passing through a given module
 842
 843   BuildIndex();
 844   TString name="DrawFirstAlignment_Layers";
 845   UShort_t voluid;
 846   Int_t **lastIndex;
 847   Int_t laymax = 0;
 848   Int_t modmax = 0;
 849   Int_t maxntr=0,nMod=0,modAligned=0,size=0;
 850   for(Int_t i=0;i<6;i++){
 851     if(layers[i]==1){
 852       size+=AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);
 853       name+=i+1;
 854     }
 855   }
 856
 857   TArrayI *volFit=new TArrayI(size);
 858   TArrayI *volFit2=new TArrayI(size-1);
 859   TArrayI *sequence=new TArrayI(size);
 860   TArrayI *volIn=new TArrayI(1);
 861
 862   // Initialize the index arrays
 863   Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
 864   lastIndex = new Int_t*[nLayers];
 865   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 866     lastIndex[iLayer] = new Int_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
 867     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 868       lastIndex[iLayer][iModule] =  fLastIndex[iLayer][iModule];
 869       if(iLayer<=(AliGeomManager::kSSD2-AliGeomManager::kFirstLayer)&&layers[iLayer]==1){
 870         volFit->AddAt(AliGeomManager::LayerToVolUID(iLayer+AliGeomManager::kFirstLayer,iModule),maxntr);
 871         maxntr++;
 872       }
 873     }
 874   }
 875   Int_t found=0;
 876   maxntr=minNtracks+1;
 877   while (maxntr>minNtracks){
 878     maxntr=minNtracks;
 879     for (Int_t iLayer = 0; iLayer <= (AliGeomManager::kSSD2 - AliGeomManager::kFirstLayer); iLayer++) {
 880       if(layers[iLayer]==0)continue;
 881       for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 882         if(lastIndex[iLayer][iModule]>maxntr){
 883           maxntr=lastIndex[iLayer][iModule];
 884           laymax=iLayer;
 885           modmax=iModule;
 886         }
 887       }
 888     }
 889     if(maxntr>minNtracks){
 890       voluid=AliGeomManager::LayerToVolUID(laymax+AliGeomManager::kFirstLayer,modmax);
 891       sequence->AddAt(voluid,nMod);
 892       lastIndex[laymax][modmax]=0;
 893       nMod++;
 894     }
 895   }
 896
 897   sequence->Set(nMod);
 898
 899  Int_t ilayer,imod;
 900   for(Int_t iter=0;iter<iterations;iter++){
 901     if(iter>0&&fDraw)UpdateDraw(sequence,iter,iter);
 902     modAligned=0;
 903     for(Int_t k=0;k<nMod;k++){
 904       TArrayI *volFit3;
 905       voluid=sequence->At(k);
 906       ilayer=AliGeomManager::VolUIDToLayer(voluid,imod);
 907       volIn->AddAt(voluid,0);
 908       found=0;
 909       if(!fitall){
 910         for(Int_t j=0;j<nMod;j++){
 911           if(j==k){
 912             found=1;
 913             continue;
 914           }
 915           else volFit2->AddAt(sequence->At(j),j-found);
 916         }
 917         volFit2->Set(nMod-1);
 918       }
 919       else{
 920         for(Int_t j=0;j<volFit->GetSize();j++){
 921           if(volFit->At(j)!=volIn->At(0))volFit2->AddAt(volFit->At(j),j-found);
 922           else found=1;
 923         }
 924       }
 925
 926       if(volidsSet){
 927         volFit3=IntersectVolArray(volidsSet,volFit2);
 928       }
 929       else volFit3=new TArrayI(*volFit2);
 930
 931
 932       if(AlignVolumesITS(volIn,volFit3,AliGeomManager::kSPD1,AliGeomManager::kTPC1,2))modAligned++;
 933       delete volFit3;
 934
 935     }
 936   }
 937   Int_t noutofsurv=CheckWithSurvey(2.,sequence);
 938   printf("%d modules into the sequence \n %d modules re-aligned \n %d modules moved far away from survey (-> reset) \n",nMod,modAligned,noutofsurv);
 939   name.Append("_iter");
 940   name+=iterations;
 941   name.Append(".root");
 942   if(fDraw)WriteHists(name.Data());
 943   delete volFit;
 944   delete volFit2;
 945   delete sequence;
 946   for(Int_t m=0;m<nLayers;m++){
 947     delete [] lastIndex[m];
 948   }
 949   delete [] lastIndex;
 950   return kTRUE;
 951
 952 }
 953
 954 //__________________________________________
 955 Bool_t AliITSRealignTracks::FirstAlignmentSPD(Int_t minNtracks,Int_t iterations,Bool_t fitall,const TArrayI *volidsSet){
 956
 957   //OBSOLETE METHOD: perform a stand-alone realignment of the SPD modules
 958   //                 based on a sequence constructed accordingly to the number of tracks
 959   //                 passing through each module
 960
 961   BuildIndex();
 962
 963   UShort_t voluid;
 964   Int_t **lastIndex;
 965   Int_t laymax = 0;
 966   Int_t modmax = 0;
 967   Int_t maxntr=0,nMod=0,modAligned=0;
 968   TArrayI *volFit=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2));
 969   TArrayI *volFit2=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2)-1);
 970   TArrayI *sequence=new TArrayI(AliGeomManager::LayerSize(1)+AliGeomManager::LayerSize(2));
 971   TArrayI *volIn=new TArrayI(1);
 972
 973   // Initialize the index arrays
 974   Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
 975   lastIndex = new Int_t*[nLayers];
 976   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
 977     lastIndex[iLayer] = new Int_t[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
 978     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 979       lastIndex[iLayer][iModule] =  fLastIndex[iLayer][iModule];
 980       if(iLayer+AliGeomManager::kFirstLayer<=AliGeomManager::kSPD2){
 981         volFit->AddAt(AliGeomManager::LayerToVolUID(iLayer+AliGeomManager::kFirstLayer,iModule),maxntr);
 982         maxntr++;
 983       }
 984     }
 985   }
 986   Int_t found=0;
 987   maxntr=minNtracks+1;
 988   while (maxntr>minNtracks){
 989     maxntr=minNtracks;
 990     for (Int_t iLayer = 0; iLayer <= (AliGeomManager::kSPD2 - AliGeomManager::kFirstLayer); iLayer++) {
 991       for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
 992         if(lastIndex[iLayer][iModule]>maxntr){
 993           laymax=iLayer;
 994           modmax=iModule;
 995           maxntr=lastIndex[iLayer][iModule];
 996         }
 997       }
 998     }
 999     if(maxntr>minNtracks){
1000       voluid=AliGeomManager::LayerToVolUID(laymax+AliGeomManager::kFirstLayer,modmax);
1001       sequence->AddAt(voluid,nMod);
1002       lastIndex[laymax][modmax]=0;
1003       nMod++;
1004       volIn->AddAt(voluid,0);
1005     }
1006   }
1007   sequence->Set(nMod);
1008
1009   Int_t ilayer,imod;
1010   for(Int_t iter=0;iter<iterations;iter++){
1011     modAligned=0;
1012     for(Int_t k=0;k<nMod;k++){
1013       TArrayI *volFit3;
1014       voluid=sequence->At(k);
1015       ilayer=AliGeomManager::VolUIDToLayer(voluid,imod);
1016       volIn->AddAt(voluid,0);
1017       found=0;
1018       if(!fitall){
1019         for(Int_t j=0;j<nMod;j++){
1020           if(j==k){
1021             found=1;
1022             continue;
1023           }
1024           else volFit2->AddAt(sequence->At(j),j-found);
1025         }
1026         volFit2->Set(nMod-1);
1027       }
1028       else{
1029         for(Int_t j=0;j<volFit->GetSize();j++){
1030           if(volFit->At(j)!=volIn->At(0))volFit2->AddAt(volFit->At(j),j-found);
1031           else found=1;
1032         }
1033       }
1034
1035       if(volidsSet){
1036         volFit3=IntersectVolArray(volidsSet,volFit2);
1037       }
1038       else volFit3=new TArrayI(*volFit2);
1039
1040
1041       if(AlignVolumesITS(volIn,volFit3,AliGeomManager::kSPD1,AliGeomManager::kSDD1,2))modAligned++;
1042       delete volFit3;
1043       //      if(volidsSet)delete volFit3;
1044     }
1045   }
1046   Int_t noutofsurv=CheckWithSurvey(2.,sequence);
1047   printf("%d modules into the sequence \n %d modules re-aligned \n %d modules moved far away from survey (-> reset) \n",nMod,modAligned,noutofsurv);
1048   delete volFit;
1049   delete volFit2;
1050   delete sequence;
1051   for(Int_t m=0;m<nLayers;m++){
1052     delete [] lastIndex[m];
1053   }
1054   delete [] lastIndex;
1055
1056   return kTRUE;
1057 }
1058
1059
1060 //__________________________________
1061 Bool_t AliITSRealignTracks::SPDmodulesAlignToSSD(Int_t minNtracks,Int_t iterations){
1062   //Align each SPD module with at least minNtracks passing through it with respect to SSD
1063   //The selection based on the minimum number of tracks is a fast one:
1064   // the number considere here doesn't coincide with the tracks effectively used then in the
1065   // minimization, it's just the total number of tracks in the sample passing through the module
1066   // The procedure is iterated "iterations" times
1067   Int_t volSSD[6]={0,0,0,0,1,1};
1068   TArrayI *volOuter=GetLayersVolUID(volSSD);
1069   TArrayI *voluid=new TArrayI(1);
1070   for(Int_t iter=0;iter<iterations;iter++){
1071     //SPD1
1072     for(Int_t imod=0;imod<AliGeomManager::LayerSize(AliGeomManager::kSPD1);imod++){    if(GetLastIndex(AliGeomManager::kSPD1-AliGeomManager::kFirstLayer,imod)<minNtracks){
1073       printf("Not enough tracks for module: lay %d mod %d \n",1,imod );
1074       continue;
1075     }
1076     voluid->AddAt(AliGeomManager::LayerToVolUID(AliGeomManager::kSPD1,imod),0);
1077     AlignVolumesITS(voluid,volOuter,AliGeomManager::kSSD1,AliGeomManager::kSSD2,2);
1078     }
1079     //SPD2
1080     for(Int_t imod=0;imod<AliGeomManager::LayerSize(AliGeomManager::kSPD2);imod++){
1081       if(GetLastIndex(AliGeomManager::kSPD2-AliGeomManager::kFirstLayer,imod)<minNtracks){
1082         printf("Not enough tracks for module: lay %d mod %d \n",2,imod );
1083         continue;
1084       }
1085       voluid->AddAt(AliGeomManager::LayerToVolUID(AliGeomManager::kSPD2,imod),0);
1086       AlignVolumesITS(voluid,volOuter,AliGeomManager::kSSD1,AliGeomManager::kSSD2,2);
1087     }
1088   }
1089   return kTRUE;
1090 }
1091
1092 //______________________________________________________________________________
1093 Bool_t AliITSRealignTracks::AlignVolumesITS(const TArrayI *volids, const TArrayI *volidsfit,
1094                                      AliGeomManager::ELayerID layerRangeMin,
1095                                      AliGeomManager::ELayerID layerRangeMax,
1096                                      Int_t iterations){
1097
1098   // Align a set of detector volumes.
1099   // Tracks are fitted only within
1100   // the range defined by the user
1101   // (by layerRangeMin and layerRangeMax)
1102   // or within the set of volidsfit
1103   // Repeat the procedure 'iterations' times
1104
1105   Int_t nVolIds = volids->GetSize();
1106   if (nVolIds == 0) {
1107     AliError("Volume IDs array is empty!");
1108     return kFALSE;
1109   }
1110   Bool_t correlated=kFALSE;
1111   Double_t surveycov[21],transl[3],rot[3],survtransl[3],survrot[3];
1112   Double_t frac;
1113
1114   TGeoHMatrix hM;
1115   Double_t phiglob,smearing,rotorig[9],normplanevect[3]={0.,0.,0.},normplanevect2[3]={0.,0.,0.};
1116   Double_t *deltarot;
1117   TMatrixDSym covmatrx(6);
1118   AliAlignObj *modAlign;
1119
1120   // Load only the tracks with at least one
1121   // space point in the set of volume (volids)
1122   BuildIndex();
1123   AliTrackPointArray **points;
1124   Bool_t failed=kFALSE;
1125   Int_t pointsdim=0,skipped=0;
1126   // Start the iterations
1127   while (iterations > 0){
1128     normplanevect2[0]=0.;
1129     normplanevect2[1]=0.;
1130     normplanevect2[2]=0.;
1131     if(fLimitCorr>0.){
1132       ResetCorrModules();
1133       skipped=0;
1134     }
1135     Int_t nArrays = LoadPoints(volids, points,pointsdim);
1136
1137     if (nArrays < fmintracks||nArrays<=0){
1138       failed=kTRUE;
1139       printf("Not enough tracks to try minimization: volUID %d and following in volids \n", volids->At(0));
1140       UnloadPoints(pointsdim, points);
1141       break;
1142     }
1143     frac=1./(Double_t)nArrays;
1144     AliTrackResiduals *minimizer = CreateMinimizer();
1145     minimizer->SetNTracks(nArrays);
1146     minimizer->InitAlignObj();
1147     AliTrackFitter *fitter = CreateFitter();
1148
1149     //Here prepare to set the plane for GetPCArot
1150                                                        //    if(volids->GetSize()==1){//TEMPORARY: to be improved
1151     AliGeomManager::GetOrigRotation(volids->At(0),rotorig);
1152     if((Int_t)AliGeomManager::VolUIDToLayer(volids->At(0))==1){//TEMPORARY: to be improved
1153       normplanevect[0]=-rotorig[1];
1154       normplanevect[1]=-rotorig[4];
1155       normplanevect[2]=0.;
1156     }
1157     else{
1158           normplanevect[0]=rotorig[1];
1159           normplanevect[1]=rotorig[4];
1160           normplanevect[2]=0.;
1161     }
1162
1163     phiglob=TMath::ATan2(normplanevect[1],normplanevect[0]);
1164
1165     modAlign=GetAlignObj(volids->At(0));
1166     modAlign->GetMatrix(hM);
1167     deltarot=hM.GetRotationMatrix();
1168     for(Int_t j=0;j<3;j++){
1169       for(Int_t i=0;i<3;i++){
1170         normplanevect2[j]+=deltarot[j*3+i]*normplanevect[i];
1171       }
1172       // printf("Here the difference: norm1[%d]=%f  norm2[%d]=%f \n",j,normplanevect[j],j,normplanevect2[j]);
1173     }
1174
1175     if(fVarySigmaY){
1176       if(modAlign->GetUniqueID()==0)smearing=fsigmaY;
1177       else{
1178         modAlign->GetCovMatrix(covmatrx);
1179         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.;
1180         //This is a sort of average: the trace with the variances of the angles
1181         //weighted with 10 cm divided per 6 and the result multiplied per 25
1182         // (the sqrt would be 5 times a sort of "mean sigma" )
1183         //
1184           }
1185     }
1186     else smearing=fsigmaY;
1187     printf("This is the sigmaY value: %f \n",smearing);
1188     // the plane will be set into the loop on tracks
1189
1190
1191     for (Int_t iArray = 0; iArray < nArrays; iArray++) {
1192       if (!points[iArray]) continue;
1193       points[iArray]->Sort(kTRUE);
1194       fitter->SetTrackPointArray(points[iArray], kFALSE);
1195       // printf("Here normplane vect: %f \n",normplanevect2[1]); //TO BE REPLACED BY      fitter->SetNormPlaneVect(normplanevect2);
1196       if (fitter->Fit(volids,volidsfit,layerRangeMin,layerRangeMax) == kFALSE) continue;
1197
1198        if(fLimitCorr>0.){
1199         correlated=kFALSE;
1200         AliTrackPoint p;
1201         Int_t layer,module;
1202         TArrayI *volparray=new TArrayI(points[iArray]->GetNPoints());
1203         for(Int_t point=0;point<points[iArray]->GetNPoints();point++){
1204           points[iArray]->GetPoint(p,point);
1205           volparray->AddAt(p.GetVolumeID(),point);
1206         }
1207         TArrayI *volpArray=ExcludeVolidsFromVolidsArray(volids,volparray);
1208         for(Int_t point=0;point<volpArray->GetSize();point++){
1209           layer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volpArray->At(point),module);
1210           if(fCorrModules[layer-AliGeomManager::kFirstLayer][module]>fLimitCorr){
1211             correlated=kTRUE;
1212             //      printf("volid %d, iarray = %d : skipping %d for Volume: %d \n",volids->At(0),iArray,skipped,volpArray->At(point));
1213             skipped++;
1214             break;
1215           }
1216         }
1217         if(!correlated){
1218           for(Int_t point=0;point<volpArray->GetSize();point++){
1219             layer=(Int_t)AliGeomManager::VolUIDToLayerSafe(volpArray->At(point),module);
1220             //printf("Number of common tracks: %d \n",fCorrModules[layer-AliGeomManager::kFirstLayer][module]);
1221             fCorrModules[layer-AliGeomManager::kFirstLayer][module]+=frac;
1222             delete volparray;
1223             delete volpArray;
1224           }
1225         }
1226         else {
1227           delete volparray;
1228           delete volpArray;
1229           continue;
1230         }
1231        }
1232
1233       AliTrackPointArray *pVolId,*pTrack;
1234       fitter->GetTrackResiduals(pVolId,pTrack);
1235       minimizer->AddTrackPointArrays(pVolId,pTrack);
1236     }
1237
1238     printf("Number of tracks considered: %d \n",nArrays);
1239     frac=(Double_t)skipped/(Double_t)nArrays;
1240     printf("Number of tracks skipped cause of correlation: %d (fraction: %f )\n",skipped,frac);
1241
1242     Int_t ntracks=minimizer->GetNFilledTracks();
1243     frac=(Double_t)ntracks/(Double_t)nArrays;
1244     printf("Number of tracks into the minimizer: %d (fraction: %f )\n",ntracks,frac);
1245     if(ntracks<=fmintracks){
1246       printf("Not enough good tracks found: could not find parameter for volume %d (and following in volids)\n",volids->At(0));
1247       UnloadPoints(pointsdim, points);
1248       failed=kTRUE;
1249       break;
1250     }
1251
1252     failed=(!minimizer->Minimize());
1253
1254     // Update the alignment object(s)
1255     if (fDoUpdate) for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
1256       UShort_t volid = (*volids)[iVolId];
1257       if(!failed){
1258         Int_t iModule;
1259         AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
1260         AliAlignObj *alignObj = fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule];
1261
1262         //Check the last minimization is not too large
1263         minimizer->GetAlignObj()->GetPars(transl,rot);
1264         fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->GetCovMatrix(surveycov);
1265         fSurveyObjs[iLayer-AliGeomManager::kFirstLayer][iModule]->GetPars(survtransl,survrot);
1266         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])){
1267           printf("Results for module %d too large: can't update them \n",volid);
1268           alignObj->SetUniqueID(2);
1269           if(iterations==1){
1270             failed=kTRUE;
1271           }
1272         }
1273         else{
1274           if(fUpdateCov){
1275             *alignObj *= *minimizer->GetAlignObj();
1276             alignObj->SetUniqueID(1);
1277           }
1278           else{
1279             alignObj->GetCovMatrix(covmatrx);
1280             *alignObj *= *minimizer->GetAlignObj();
1281             alignObj->SetCorrMatrix(covmatrx);
1282             alignObj->SetUniqueID(1);
1283           }
1284
1285           /*alignObj->GetPars(transl,rot);
1286
1287           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])){
1288           printf("Results for module %d out of Survey: reinitializing it from survey \n",volid);
1289             //    *alignObj = *alignObjSurv;
1290             alignObj->SetPars(0.,0.,0.,0.,0.,0.);
1291             alignObj->SetUniqueID(0);
1292             if(fUpdateCov)alignObj->SetCorrMatrix(surveycov);
1293             if(iterations==1){
1294             failed=kTRUE;
1295             }
1296             }*/
1297         }
1298         if(iterations==1)alignObj->Print("");
1299       }
1300       else {
1301         printf("Minimization failed: cannot update AlignObj for volume: %d \n",volid);
1302       }
1303     }
1304     UnloadPoints(pointsdim,points);
1305     if(failed)break;
1306     minimizer->InitAlignObj();
1307     iterations--;
1308   }
1309
1310   printf("\n \n");
1311
1312   return (!failed);
1313 }
1314
1315
1316
1317 //______________________________________________
1318 Bool_t AliITSRealignTracks::AlignSPDBarrel(Int_t iterations){
1319   //Align the SPD barrel "iterations" times
1320
1321   Int_t size=0,size2=0;
1322   Int_t layers[6]={1,1,0,0,0,0};
1323   for(Int_t k=1;k<=2;k++){
1324     size+=AliGeomManager::LayerSize(k);
1325   }
1326   for(Int_t k=3;k<=6;k++){
1327     size2+=AliGeomManager::LayerSize(k);
1328     printf("size: %d \n",size2);
1329   }
1330
1331   printf("Aligning SPDBarrel: nmodules: %d \n",size);
1332   printf("Fitting modules: %d \n",size2);
1333
1334   TArrayI *volIDs=GetLayersVolUID(layers);
1335   layers[0]=0;
1336   layers[1]=0;
1337   layers[2]=1;
1338   layers[3]=1;
1339   layers[4]=1;
1340   layers[5]=1;
1341   TArrayI *volIDsFit=GetLayersVolUID(layers);
1342
1343   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSDD1,AliGeomManager::kTPC1,iterations);
1344
1345   return kTRUE;
1346 }
1347
1348 //______________________
1349 Bool_t AliITSRealignTracks::AlignSPDHalfBarrel(Int_t method,Int_t iterations){
1350   //Align a SPD Half barrel "iterations" times
1351   //method 0 : align SPDHalfBarrel Up without using the points on SPD Half Barrel down in the fits (only outer layers)
1352   //method 1 : align SPDHalfBarrel Down without using the points on SPD Half Barrel up in the fits (only outer layers)
1353   //method 10 : align SPDHalfBarrel Up using also the points on SPD Half Barrel down in the fits (and points on outer layers)
1354   //method 11 : align SPDHalfBarrel Down using also the points on SPD Half Barrel up in the fits (and points on outer layers)
1355
1356   Int_t size=0,size2=0;
1357   Int_t layers[6]={0,0,1,1,1,1};
1358   Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1359   Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1360
1361   TString updownstr;
1362   if(method==0)updownstr="UpNoDown";
1363   else if (method==1)updownstr="DownNoUp";
1364   else if (method==10)updownstr="UpWithDown";
1365   else if (method==11)updownstr="DownWithUp";
1366   else {
1367     AliWarning("Wrong AlignSPDHalfBarrel method selected ");
1368     return kFALSE;
1369   }
1370
1371   for(Int_t i=1;i<=2;i++){
1372     size+=AliGeomManager::LayerSize(i);
1373   }
1374
1375   for(Int_t i=3;i<=6;i++){
1376     size2+=AliGeomManager::LayerSize(i);
1377   }
1378
1379   size=size/2;
1380   if(method==10||method==11)size2+=size;
1381
1382   printf("Aligning  SPDHalfBarrel %s: nmodules: %d \n",updownstr.Data(),size);
1383   printf("Fitting modules: %d \n",size2);
1384   TArrayI *volIDsFit2;
1385   TArrayI *volids = NULL;
1386   TArrayI *volIDsFit=GetLayersVolUID(layers);
1387   if(method==0||method==10)volids=GetSPDSectorsVolids(sectorsUp);
1388   if(method==1||method==11)volids=GetSPDSectorsVolids(sectorsDown);
1389
1390   if(method==10)volIDsFit2=JoinVolArrays(GetSPDSectorsVolids(sectorsDown),volIDsFit);
1391   else if(method==11)volIDsFit2=JoinVolArrays(GetSPDSectorsVolids(sectorsUp),volIDsFit);
1392   else volIDsFit2=volIDsFit;
1393
1394   AlignVolumesITS(volids,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1395
1396   return kTRUE;
1397 }
1398
1399
1400 //______________________________________________________
1401 Bool_t AliITSRealignTracks::AlignLayer(Int_t layer,Int_t iterations){
1402   //Align the layer "layer" iterations times
1403
1404   Int_t size=0,size2=0;
1405   Int_t layers[6]={0,0,0,0,0,0};
1406   layers[layer-1]=1;
1407   TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1408   for(Int_t k=1;k<=6;k++){
1409     if(k!=layer)size2+=AliGeomManager::LayerSize(k);
1410   }
1411   size=AliGeomManager::LayerSize(layer);
1412
1413   printf("Aligning layer %s, nmodules %d ,fitted modules %d \n",layerstr[layer-1].Data(),size,size2);
1414
1415
1416   TArrayI *volIDs=GetLayersVolUID(layers);
1417   layers[0]=1;
1418   layers[1]=1;
1419   layers[2]=1;
1420   layers[3]=1;
1421   layers[4]=1;
1422   layers[5]=1;
1423   layers[layer]=0;
1424   TArrayI *volIDsFit=GetLayersVolUID(layers);
1425
1426   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSDD1,AliGeomManager::kSSD2,iterations);
1427
1428   return kTRUE;
1429 }
1430
1431 //___________________________________________
1432
1433 Bool_t AliITSRealignTracks::AlignLayersToLayers(const Int_t *layer,Int_t iterations){
1434
1435   //Align the set of layers A with respect to the set of layers B iterations time.
1436   //The two sets A and B are defined into *layer==layer[6] the following way:
1437   //   layer[i]=0 the layer is skipped both in the fits than in the minimization
1438   //   layer[i]=1 the layer is skipped in the fits and considered in the minimization
1439   //   layer[i]=2 the layer is considered in the fits and skipped in the minimization
1440   //   layer[i]=3 the layer is considered both in the fits and in the minimization
1441
1442   UShort_t volid;
1443   Int_t size=0,size2=0,j=0,k=0;
1444   Int_t iLayer;
1445   TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1446   TString command="",str;
1447   for(Int_t i=1;i<=6;i++){
1448     if(layer[i-1]==1||layer[i-1]==3){
1449       size+=AliGeomManager::LayerSize(i);
1450       command.Append(" ");
1451       command.Append(layerstr[i-1]);
1452     }
1453     if(layer[i-1]==2||layer[i-1]==3){
1454       size2+=AliGeomManager::LayerSize(i);
1455       str.Append(" ");
1456       str.Append(layerstr[i-1]);
1457     }
1458   }
1459
1460   printf("Aligning layers %s To layers %s, nmodules %d ,fitted modules %d \n",command.Data(),str.Data(),size,size2);
1461
1462
1463   TArrayI volIDs(size);
1464   TArrayI volIDsFit(size2);
1465
1466   for (iLayer=(Int_t)AliGeomManager::kSPD1;iLayer<(Int_t)AliGeomManager::kTPC1;iLayer++){
1467     if(layer[iLayer-AliGeomManager::kFirstLayer]==0)continue;
1468     if(layer[iLayer-AliGeomManager::kFirstLayer]==1){
1469       for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1470         volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1471         volIDs.AddAt(volid,j);
1472         j++;
1473       }
1474     }
1475     else if(layer[iLayer-AliGeomManager::kFirstLayer]==2){
1476       for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1477         volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1478         volIDsFit.AddAt(volid,k);
1479         k++;
1480       }
1481     }
1482     else if(layer[iLayer-AliGeomManager::kFirstLayer]==3){
1483       for (Int_t iModule2=0;iModule2<AliGeomManager::LayerSize(iLayer);iModule2++){
1484         volid = AliGeomManager::LayerToVolUID(iLayer,iModule2);
1485         volIDs.AddAt(volid,j);
1486         j++;
1487         volIDsFit.AddAt(volid,k);
1488         k++;
1489       }
1490     }
1491   }
1492
1493   AlignVolumesITS(&volIDs,&volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1494
1495   return kTRUE;
1496 }
1497
1498 //______________________________________________
1499
1500 Bool_t AliITSRealignTracks::AlignSPDSectorToOuterLayers(Int_t sector,Int_t iterations){
1501   //Align the SPD sector "sector" with respect to outer layers iterations times
1502
1503
1504   Int_t layers[6]={0,0,1,1,1,1};
1505   Bool_t spd=kFALSE;
1506   Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1507   Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1508
1509   if(sector<0){
1510     sector=-sector;
1511     spd=kTRUE;
1512   }
1513   sectorsIN[sector]=1;
1514   sectorsFit[sector]=0;
1515   TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1516   TArrayI *volIDsFit;
1517   if(spd){
1518     volIDsFit=JoinVolArrays(GetSPDSectorsVolids(sectorsFit),GetLayersVolUID(layers));
1519   }
1520   else volIDsFit=GetLayersVolUID(layers);
1521
1522   printf("Aligning SPD sector %d: nmodules: %d \n",sector,volIDs->GetSize());
1523   printf("Fitting modules: %d \n",volIDsFit->GetSize());
1524
1525   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1526
1527   return kTRUE;
1528 }
1529
1530 //______________________________________________
1531 Bool_t AliITSRealignTracks::AlignSPDSectorWithSectors(Int_t sector,Int_t iterations){
1532   //Align the SPD sector "sector" with respect to the other SPD sectors iterations times
1533
1534   Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1535   Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1536
1537   sectorsIN[sector]=1;
1538   sectorsFit[sector]=0;
1539   TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1540   TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);;
1541
1542   printf("Aligning SPD sector %d: nmodules: %d \n",sector,volIDs->GetSize());
1543   printf("Fitting modules: %d \n",volIDsFit->GetSize());
1544
1545
1546   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSDD1,iterations);
1547
1548   return kTRUE;
1549 }
1550
1551
1552
1553 //___________________________________________________
1554 Bool_t AliITSRealignTracks::AlignSPDSectorsWithSectors(const Int_t *sectorsIN,const Int_t *sectorsFit,Int_t iterations){
1555   //Align SPD sectors defined in "sectorsIN" with respect to
1556   //SPD sectors defined in "sectorsFit" iterations time
1557
1558   TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1559   TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);;
1560
1561   printf("Aligning SPD sectors: modules: %d \n",volIDs->GetSize());
1562   printf("Fitting modules: %d \n",volIDsFit->GetSize());
1563
1564
1565
1566   return AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSDD1,iterations);;
1567 }
1568
1569 //___________________________________________________
1570 Bool_t AliITSRealignTracks::AlignSPDStaves(const Int_t *staves,const Int_t *sectorsIN,const Int_t *sectorsFit,Int_t iterations){
1571   //Align SPD staves defined by staves and sectorsIN with respect to sectorsFit volumes iterations times
1572
1573   TArrayI *volIDs=GetSPDStavesVolids(sectorsIN,staves);
1574   TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1575
1576   if(volIDs->GetSize()==0){
1577     printf("EMPTY ARRAY !! \n");
1578     return kFALSE;
1579   }
1580   printf("Aligning SPD staves: modules: %d \n",volIDs->GetSize());
1581   printf("Fitting modules: %d \n",volIDsFit->GetSize());
1582
1583   TArrayI *volIDsFit2=ExcludeVolidsFromVolidsArray(volIDs,volIDsFit);
1584   return  AlignVolumesITS(volIDs,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSSD1,iterations);
1585 }
1586
1587
1588 //___________________________________________
1589
1590 Bool_t AliITSRealignTracks::AlignLayerToSPDHalfBarrel(Int_t layer,Int_t updown,Int_t iterations){
1591   //Align the layer "layer" with respect to SPD Half Barrel Up (updowon=0)
1592   //or Down (updown=1) iterations times
1593
1594
1595
1596   Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1597   Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1598   TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1599   TArrayI *volIDsFit;
1600   Int_t layers[6]={0,0,0,0,0,0};
1601   layers[layer-1]=1;
1602   Int_t size=AliGeomManager::LayerSize(layer);
1603   TArrayI *volIDs=GetLayersVolUID(layers);
1604
1605   if(updown==0){
1606     volIDsFit=GetSPDSectorsVolids(sectorsUp);
1607     printf("Aligning layer %s, nmodules %d ,to half barrel Up \n",layerstr[layer-1].Data(),size);
1608   }
1609   else if(updown==1){
1610     volIDsFit=GetSPDSectorsVolids(sectorsDown);
1611     printf("Aligning layer %s, nmodules %d ,to half barrel Down \n",layerstr[layer-1].Data(),size);
1612   }
1613   else {
1614     printf("Wrong Half Barrel selection! \n");
1615     return kFALSE;
1616   }
1617
1618   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1619
1620   return kTRUE;
1621 }
1622
1623 //___________________________________________
1624
1625 Bool_t AliITSRealignTracks::AlignLayerToSector(Int_t layer,Int_t sector,Int_t iterations){
1626   //Align the layer "layer" with respect to SPD sector "sector" iterations times
1627
1628   if(sector>9){
1629     printf("Wrong Sector selection! \n");
1630     return kFALSE;
1631   }
1632   Int_t sectors[10]={0,0,0,0,0,0,0,0,0,0};
1633   sectors[sector]=1;
1634   TString layerstr[6]={"SPD1","SPD2","SDD1","SDD2","SSD1","SSD2"};
1635   TArrayI *volIDsFit;
1636   Int_t layers[6]={0,0,0,0,0,0};
1637   layers[layer-1]=1;
1638   TArrayI *volIDs=GetLayersVolUID(layers);
1639   Int_t size=AliGeomManager::LayerSize(layer);
1640
1641
1642   volIDsFit=GetSPDSectorsVolids(sectors);
1643   printf("Aligning layer %s, nmodules %d ,to half barrel Up \n",layerstr[layer-1].Data(),size);
1644
1645   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSSD2,iterations);
1646
1647   return kTRUE;
1648 }
1649
1650 //_______________________________________________
1651
1652 Bool_t AliITSRealignTracks::AlignSPDHalfBarrelToHalfBarrel(Int_t updown,Int_t iterations){
1653   //Align the SPD Half Barrel Up[Down] with respect to HB Down[Up] iterations time if
1654   //updown=0[1]
1655
1656
1657   Int_t sectorsDown[10]={0,0,0,0,0,1,1,1,1,1};
1658   Int_t sectorsUp[10]={1,1,1,1,1,0,0,0,0,0};
1659
1660   TArrayI *volIDsUp=GetSPDSectorsVolids(sectorsUp);
1661   TArrayI *volIDsDown=GetSPDSectorsVolids(sectorsDown);
1662
1663   if(updown==0){
1664     printf("Aligning SPD HalfBarrel up to half Barrel down : nmodules: %d \n",volIDsUp->GetSize());
1665     printf("Fitting modules: %d \n",volIDsDown->GetSize());
1666     AlignVolumesITS(volIDsUp,volIDsDown,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1667   }
1668   else if(updown==1){
1669     printf("Aligning SPD HalfBarrel down to half Barrel Up : nmodules: %d \n",volIDsDown->GetSize());
1670     printf("Fitting modules: %d \n",volIDsUp->GetSize());
1671     AlignVolumesITS(volIDsDown,volIDsUp,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1672   }
1673   else {
1674     printf("Wrong Half Barrel selection! \n");
1675     return kFALSE;
1676   }
1677
1678   return kTRUE;
1679 }
1680
1681
1682 //_______________________
1683 Bool_t AliITSRealignTracks::AlignSPDHalfBarrelToSectorRef(Int_t sector,Int_t iterations){
1684   //Align the SPD Half Barrel Down with respect to sector "sector" iterations times
1685
1686   Int_t sectorsIN[10]={0,0,0,0,0,1,1,1,1,1};
1687   Int_t sectorsFit[10]={0,0,0,0,0,0,0,0,0,0};
1688
1689   sectorsFit[sector]=1;
1690
1691   TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1692   TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1693
1694   printf("Aligning SPD HalfBarrel to sector 0 %d: nmodules: %d \n",sector,volIDs->GetSize());
1695   printf("Fitting modules: %d \n",volIDsFit->GetSize());
1696
1697
1698   AlignVolumesITS(volIDs,volIDsFit,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1699
1700
1701   return kTRUE;
1702 }
1703 //_________________________________________
1704 Bool_t AliITSRealignTracks::AlignSPD1SectorRef(Int_t sector,Int_t iterations){
1705   //OBSOLETE METHOD: Align the SPD1 modules of sector "sector" with respect
1706   // to the other SPD volumes iterations times
1707
1708   Int_t sectorsIN[10]={0,0,0,0,0,0,0,0,0,0};
1709   Int_t sectorsFit[10]={1,1,1,1,1,1,1,1,1,1};
1710   sectorsIN[sector]=1;
1711   sectorsFit[sector]=0;
1712   TArrayI *volIDs=GetSPDSectorsVolids(sectorsIN);
1713   TArrayI *volIDsFit=GetSPDSectorsVolids(sectorsFit);
1714   Int_t size=volIDs->GetSize();
1715   Int_t size2=volIDsFit->GetSize();
1716   UShort_t volID;
1717   Int_t k=0;
1718
1719   TArrayI *volIDsSPD1=new TArrayI(size-8);
1720   TArrayI *volIDsFit2=new TArrayI(size2+8);
1721
1722   for(Int_t j=0;j<size;j++){
1723     volID=volIDs->At(j);
1724     if(AliGeomManager::VolUIDToLayer(volID)==AliGeomManager::kSPD1){
1725       volIDsSPD1->AddAt(volID,size2+k);
1726       k++;
1727     }
1728     else volIDsFit2->AddAt(volID,j-k);
1729   }
1730
1731
1732   for(Int_t j=0;j<size2;j++){
1733     volID=volIDsFit->At(j);
1734     volIDsFit2->AddAt(volID,size-k+j);
1735   }
1736
1737   printf("Aligning SPD Sector %d: nmodules: %d \n",sector,volIDsSPD1->GetSize());
1738   printf("Fitting modules: %d \n",volIDsFit2->GetSize());
1739
1740   AlignVolumesITS(volIDsSPD1,volIDsFit2,AliGeomManager::kSPD1,AliGeomManager::kSPD2,iterations);
1741
1742   return kTRUE;
1743 }
1744
1745 //_____________________________________________
1746
1747 AliAlignObjParams* AliITSRealignTracks::MediateAlignObj(const TArrayI *volIDs,Int_t lastVolid){
1748   //TEMPORARY METHOD: perform an average of the values of the parameters of the AlignObjs
1749   // defined by the array volIDs up to lastVolid position in this array
1750   //The aim of such a method is to look for collective movement of a given set of modules
1751
1752   UShort_t volid;
1753
1754   TGeoHMatrix hm;
1755   Double_t *rot,*transl;
1756   Double_t rotSum[9],translSum[3]={0.,0.,0.};
1757   for(Int_t k=0;k<8;k++)rotSum[k]=0.;
1758
1759
1760   for(Int_t ivol=0;ivol<lastVolid;ivol++){
1761     volid=volIDs->At(ivol);
1762
1763     GetAlignObj(volIDs->At(ivol))->GetMatrix(hm);
1764
1765     rot=hm.GetRotationMatrix();
1766     transl=hm.GetTranslation();
1767
1768     for(Int_t j=0;j<9;j++)rotSum[j]+=rot[j];
1769     for(Int_t jt=0;jt<3;jt++)translSum[jt]+=transl[jt];
1770   }
1771   if(lastVolid!=0){
1772     for(Int_t j=0;j<9;j++)rotSum[j]=rotSum[j]/lastVolid;
1773     for(Int_t jt=0;jt<3;jt++)translSum[jt]=translSum[jt]/lastVolid;
1774   }
1775   else printf("Try to mediate results for zero modules \n");
1776
1777   hm.SetRotation(rotSum);
1778   hm.SetTranslation(translSum);
1779
1780
1781
1782   AliAlignObjParams *alignObj=new AliAlignObjParams("average", 0,hm, kTRUE);
1783   return alignObj;
1784
1785 }
1786
1787
1788 //________________________________________________
1789 TArrayI* AliITSRealignTracks::GetSPDStavesVolids(const Int_t *sectors,const Int_t* staves){
1790
1791
1792   // This method gets the volID Array for the chosen staves into the
1793   // chosen sectors. You have to pass an array (10 dim) with a 1 for each
1794   // selected sector and an array (6 dim) with a 1 for each chosen stave.
1795   // The staves are numbered in this way: 0,1 for SPD1 and 2,3,4,5 for SPD2
1796   // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
1797   // staves[6]={0,1,1,0,0,1} -> Staves 1 on SPD1 and 0 and 3 on SPD2 selected
1798
1799   Int_t nSect=0,nStaves=0;
1800   Int_t last=0;
1801
1802
1803   for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1804     if(sectors[co]==1) nSect++;
1805   }
1806
1807   for(Int_t co=0;co<6;co++){ //counts the number of sectors chosen
1808     if(staves[co]==1) nStaves++;
1809   }
1810
1811   if(nSect<1||nStaves<1){ //if no sector chosen -> exit
1812     Printf("Error! No Sector/s or staves Selected!");
1813     return 0x0;
1814   }
1815
1816   TArrayI *volIDs = new TArrayI(nSect*nStaves*4);
1817   TString stave="/Stave",str,symn,laystr;
1818
1819   TArrayI *sectvol=GetSPDSectorsVolids(sectors);
1820   //SPD1
1821   laystr="SPD0";
1822   for(Int_t k=0;k<2;k++){
1823     if(staves[k]==1){
1824       str=stave;
1825       str+=k;
1826       for(Int_t i=0;i<sectvol->GetSize();i++){
1827         symn=AliGeomManager::SymName(sectvol->At(i));
1828         if(symn.Contains(str)&&symn.Contains(laystr)){
1829           //      printf("Adding: %s \n",symn.Data());
1830           volIDs->AddAt(sectvol->At(i),last);
1831           last++;
1832         }
1833       }
1834     }
1835   }
1836   //SPD1
1837   laystr="SPD1";
1838   for(Int_t k=2;k<6;k++){
1839     if(staves[k]==1){
1840       str=stave;
1841       str+=k-2;
1842       for(Int_t i=0;i<sectvol->GetSize();i++){
1843         symn=AliGeomManager::SymName(sectvol->At(i));
1844         if(symn.Contains(str)&&symn.Contains(laystr)){
1845           volIDs->AddAt(sectvol->At(i),last);
1846           printf("Adding: %s \n",symn.Data());
1847           last++;
1848         }
1849       }
1850     }
1851   }
1852
1853   volIDs->Set(last);
1854   return volIDs;
1855 }
1856
1857 //________________________________________________
1858 TArrayI* AliITSRealignTracks::GetSPDSectorsVolids(const Int_t *sectors)
1859 {
1860   //
1861   // This method gets the volID Array for the chosen sectors.
1862   // You have to pass an array with a 1 for each selected sector.
1863   // i.e. sectors[10] = {1,1,0,0,0,0,0,0,1,0} -> Sector 0, 1, 9 selected.
1864   //
1865
1866   Int_t nSect=0;
1867   Int_t iModule=0;
1868
1869
1870   for(Int_t co=0;co<10;co++){ //counts the number of sectors chosen
1871     if(sectors[co]==1) nSect++;
1872   }
1873
1874   if(nSect<1){ //if no sector chosen -> exit
1875     Printf("Error! No Sector/s Selected!");
1876     return 0x0;
1877   }
1878
1879   TArrayI *volIDs = new TArrayI(nSect*24);
1880
1881     if(sectors[0]==1){ //--->cSect = 0 <---
1882       volIDs->AddAt(2048,iModule); iModule++;
1883       volIDs->AddAt(2049,iModule); iModule++;
1884       volIDs->AddAt(2050,iModule); iModule++;
1885       volIDs->AddAt(2051,iModule); iModule++;
1886       volIDs->AddAt(2052,iModule); iModule++;
1887       volIDs->AddAt(2053,iModule); iModule++;
1888       volIDs->AddAt(2054,iModule); iModule++;
1889       volIDs->AddAt(2055,iModule); iModule++;
1890       volIDs->AddAt(4096,iModule); iModule++;
1891       volIDs->AddAt(4097,iModule); iModule++;
1892       volIDs->AddAt(4098,iModule); iModule++;
1893       volIDs->AddAt(4099,iModule); iModule++;
1894       volIDs->AddAt(4100,iModule); iModule++;
1895       volIDs->AddAt(4101,iModule); iModule++;
1896       volIDs->AddAt(4102,iModule); iModule++;
1897       volIDs->AddAt(4103,iModule); iModule++;
1898       volIDs->AddAt(4104,iModule); iModule++;
1899       volIDs->AddAt(4105,iModule); iModule++;
1900       volIDs->AddAt(4106,iModule); iModule++;
1901       volIDs->AddAt(4107,iModule); iModule++;
1902       volIDs->AddAt(4108,iModule); iModule++;
1903       volIDs->AddAt(4109,iModule); iModule++;
1904       volIDs->AddAt(4110,iModule); iModule++;
1905       volIDs->AddAt(4111,iModule); iModule++;
1906     }
1907     if(sectors[1]==1){ //--->cSect = 1 <//---
1908       volIDs->AddAt(2056,iModule); iModule++;
1909       volIDs->AddAt(2057,iModule); iModule++;
1910       volIDs->AddAt(2058,iModule); iModule++;
1911       volIDs->AddAt(2059,iModule); iModule++;
1912       volIDs->AddAt(2060,iModule); iModule++;
1913       volIDs->AddAt(2061,iModule); iModule++;
1914       volIDs->AddAt(2062,iModule); iModule++;
1915       volIDs->AddAt(2063,iModule); iModule++;
1916       volIDs->AddAt(4112,iModule); iModule++;
1917       volIDs->AddAt(4113,iModule); iModule++;
1918       volIDs->AddAt(4114,iModule); iModule++;
1919       volIDs->AddAt(4115,iModule); iModule++;
1920       volIDs->AddAt(4116,iModule); iModule++;
1921       volIDs->AddAt(4117,iModule); iModule++;
1922       volIDs->AddAt(4118,iModule); iModule++;
1923       volIDs->AddAt(4119,iModule); iModule++;
1924       volIDs->AddAt(4120,iModule); iModule++;
1925       volIDs->AddAt(4121,iModule); iModule++;
1926       volIDs->AddAt(4122,iModule); iModule++;
1927       volIDs->AddAt(4123,iModule); iModule++;
1928       volIDs->AddAt(4124,iModule); iModule++;
1929       volIDs->AddAt(4125,iModule); iModule++;
1930       volIDs->AddAt(4126,iModule); iModule++;
1931       volIDs->AddAt(4127,iModule); iModule++;
1932     }
1933     if(sectors[2]==1){//--->cSect = 2 <//---
1934       volIDs->AddAt(2064,iModule); iModule++;
1935       volIDs->AddAt(2065,iModule); iModule++;
1936       volIDs->AddAt(2066,iModule); iModule++;
1937       volIDs->AddAt(2067,iModule); iModule++;
1938       volIDs->AddAt(2068,iModule); iModule++;
1939       volIDs->AddAt(2069,iModule); iModule++;
1940       volIDs->AddAt(2070,iModule); iModule++;
1941       volIDs->AddAt(2071,iModule); iModule++;
1942       volIDs->AddAt(4128,iModule); iModule++;
1943       volIDs->AddAt(4129,iModule); iModule++;
1944       volIDs->AddAt(4130,iModule); iModule++;
1945       volIDs->AddAt(4131,iModule); iModule++;
1946       volIDs->AddAt(4132,iModule); iModule++;
1947       volIDs->AddAt(4133,iModule); iModule++;
1948       volIDs->AddAt(4134,iModule); iModule++;
1949       volIDs->AddAt(4135,iModule); iModule++;
1950       volIDs->AddAt(4136,iModule); iModule++;
1951       volIDs->AddAt(4137,iModule); iModule++;
1952       volIDs->AddAt(4138,iModule); iModule++;
1953       volIDs->AddAt(4139,iModule); iModule++;
1954       volIDs->AddAt(4140,iModule); iModule++;
1955       volIDs->AddAt(4141,iModule); iModule++;
1956       volIDs->AddAt(4142,iModule); iModule++;
1957       volIDs->AddAt(4143,iModule); iModule++;
1958     }
1959     if(sectors[3]==1){//--->cSect = 3 <//---
1960       volIDs->AddAt(2072,iModule); iModule++;
1961       volIDs->AddAt(2073,iModule); iModule++;
1962       volIDs->AddAt(2074,iModule); iModule++;
1963       volIDs->AddAt(2075,iModule); iModule++;
1964       volIDs->AddAt(2076,iModule); iModule++;
1965       volIDs->AddAt(2077,iModule); iModule++;
1966       volIDs->AddAt(2078,iModule); iModule++;
1967       volIDs->AddAt(2079,iModule); iModule++;
1968       volIDs->AddAt(4144,iModule); iModule++;
1969       volIDs->AddAt(4145,iModule); iModule++;
1970       volIDs->AddAt(4146,iModule); iModule++;
1971       volIDs->AddAt(4147,iModule); iModule++;
1972       volIDs->AddAt(4148,iModule); iModule++;
1973       volIDs->AddAt(4149,iModule); iModule++;
1974       volIDs->AddAt(4150,iModule); iModule++;
1975       volIDs->AddAt(4151,iModule); iModule++;
1976       volIDs->AddAt(4152,iModule); iModule++;
1977       volIDs->AddAt(4153,iModule); iModule++;
1978       volIDs->AddAt(4154,iModule); iModule++;
1979       volIDs->AddAt(4155,iModule); iModule++;
1980       volIDs->AddAt(4156,iModule); iModule++;
1981       volIDs->AddAt(4157,iModule); iModule++;
1982       volIDs->AddAt(4158,iModule); iModule++;
1983       volIDs->AddAt(4159,iModule); iModule++;
1984     }
1985     if(sectors[4]==1){//--->cSect = 4 <//---
1986       volIDs->AddAt(2080,iModule); iModule++;
1987       volIDs->AddAt(2081,iModule); iModule++;
1988       volIDs->AddAt(2082,iModule); iModule++;
1989       volIDs->AddAt(2083,iModule); iModule++;
1990       volIDs->AddAt(2084,iModule); iModule++;
1991       volIDs->AddAt(2085,iModule); iModule++;
1992       volIDs->AddAt(2086,iModule); iModule++;
1993       volIDs->AddAt(2087,iModule); iModule++;
1994       volIDs->AddAt(4160,iModule); iModule++;
1995       volIDs->AddAt(4161,iModule); iModule++;
1996       volIDs->AddAt(4162,iModule); iModule++;
1997       volIDs->AddAt(4163,iModule); iModule++;
1998       volIDs->AddAt(4164,iModule); iModule++;
1999       volIDs->AddAt(4165,iModule); iModule++;
2000       volIDs->AddAt(4166,iModule); iModule++;
2001       volIDs->AddAt(4167,iModule); iModule++;
2002       volIDs->AddAt(4168,iModule); iModule++;
2003       volIDs->AddAt(4169,iModule); iModule++;
2004       volIDs->AddAt(4170,iModule); iModule++;
2005       volIDs->AddAt(4171,iModule); iModule++;
2006       volIDs->AddAt(4172,iModule); iModule++;
2007       volIDs->AddAt(4173,iModule); iModule++;
2008       volIDs->AddAt(4174,iModule); iModule++;
2009       volIDs->AddAt(4175,iModule); iModule++;
2010     }
2011     if(sectors[5]==1){//--->cSect = 5 <//---
2012       volIDs->AddAt(2088,iModule); iModule++;
2013       volIDs->AddAt(2089,iModule); iModule++;
2014       volIDs->AddAt(2090,iModule); iModule++;
2015       volIDs->AddAt(2091,iModule); iModule++;
2016       volIDs->AddAt(2092,iModule); iModule++;
2017       volIDs->AddAt(2093,iModule); iModule++;
2018       volIDs->AddAt(2094,iModule); iModule++;
2019       volIDs->AddAt(2095,iModule); iModule++;
2020       volIDs->AddAt(4176,iModule); iModule++;
2021       volIDs->AddAt(4177,iModule); iModule++;
2022       volIDs->AddAt(4178,iModule); iModule++;
2023       volIDs->AddAt(4179,iModule); iModule++;
2024       volIDs->AddAt(4180,iModule); iModule++;
2025       volIDs->AddAt(4181,iModule); iModule++;
2026       volIDs->AddAt(4182,iModule); iModule++;
2027       volIDs->AddAt(4183,iModule); iModule++;
2028       volIDs->AddAt(4184,iModule); iModule++;
2029       volIDs->AddAt(4185,iModule); iModule++;
2030       volIDs->AddAt(4186,iModule); iModule++;
2031       volIDs->AddAt(4187,iModule); iModule++;
2032       volIDs->AddAt(4188,iModule); iModule++;
2033       volIDs->AddAt(4189,iModule); iModule++;
2034       volIDs->AddAt(4190,iModule); iModule++;
2035       volIDs->AddAt(4191,iModule); iModule++;
2036     }
2037     if(sectors[6]==1){//--->cSect = 6 <//---
2038       volIDs->AddAt(2096,iModule); iModule++;
2039       volIDs->AddAt(2097,iModule); iModule++;
2040       volIDs->AddAt(2098,iModule); iModule++;
2041       volIDs->AddAt(2099,iModule); iModule++;
2042       volIDs->AddAt(2100,iModule); iModule++;
2043       volIDs->AddAt(2101,iModule); iModule++;
2044       volIDs->AddAt(2102,iModule); iModule++;
2045       volIDs->AddAt(2103,iModule); iModule++;
2046       volIDs->AddAt(4192,iModule); iModule++;
2047       volIDs->AddAt(4193,iModule); iModule++;
2048       volIDs->AddAt(4194,iModule); iModule++;
2049       volIDs->AddAt(4195,iModule); iModule++;
2050       volIDs->AddAt(4196,iModule); iModule++;
2051       volIDs->AddAt(4197,iModule); iModule++;
2052       volIDs->AddAt(4198,iModule); iModule++;
2053       volIDs->AddAt(4199,iModule); iModule++;
2054       volIDs->AddAt(4200,iModule); iModule++;
2055       volIDs->AddAt(4201,iModule); iModule++;
2056       volIDs->AddAt(4202,iModule); iModule++;
2057       volIDs->AddAt(4203,iModule); iModule++;
2058       volIDs->AddAt(4204,iModule); iModule++;
2059       volIDs->AddAt(4205,iModule); iModule++;
2060       volIDs->AddAt(4206,iModule); iModule++;
2061       volIDs->AddAt(4207,iModule); iModule++;
2062     }
2063      if(sectors[7]==1){ //--->cSect = 7 <//---
2064        volIDs->AddAt(2104,iModule); iModule++;
2065        volIDs->AddAt(2105,iModule); iModule++;
2066        volIDs->AddAt(2106,iModule); iModule++;
2067        volIDs->AddAt(2107,iModule); iModule++;
2068        volIDs->AddAt(2108,iModule); iModule++;
2069        volIDs->AddAt(2109,iModule); iModule++;
2070        volIDs->AddAt(2110,iModule); iModule++;
2071        volIDs->AddAt(2111,iModule); iModule++;
2072        volIDs->AddAt(4208,iModule); iModule++;
2073        volIDs->AddAt(4209,iModule); iModule++;
2074        volIDs->AddAt(4210,iModule); iModule++;
2075        volIDs->AddAt(4211,iModule); iModule++;
2076        volIDs->AddAt(4212,iModule); iModule++;
2077        volIDs->AddAt(4213,iModule); iModule++;
2078        volIDs->AddAt(4214,iModule); iModule++;
2079        volIDs->AddAt(4215,iModule); iModule++;
2080        volIDs->AddAt(4216,iModule); iModule++;
2081        volIDs->AddAt(4217,iModule); iModule++;
2082        volIDs->AddAt(4218,iModule); iModule++;
2083        volIDs->AddAt(4219,iModule); iModule++;
2084        volIDs->AddAt(4220,iModule); iModule++;
2085        volIDs->AddAt(4221,iModule); iModule++;
2086        volIDs->AddAt(4222,iModule); iModule++;
2087        volIDs->AddAt(4223,iModule); iModule++;
2088      }
2089      if(sectors[8]==1){//--->cSect = 8 <//---
2090        volIDs->AddAt(2112,iModule); iModule++;
2091        volIDs->AddAt(2113,iModule); iModule++;
2092        volIDs->AddAt(2114,iModule); iModule++;
2093        volIDs->AddAt(2115,iModule); iModule++;
2094        volIDs->AddAt(2116,iModule); iModule++;
2095        volIDs->AddAt(2117,iModule); iModule++;
2096        volIDs->AddAt(2118,iModule); iModule++;
2097        volIDs->AddAt(2119,iModule); iModule++;
2098        volIDs->AddAt(4224,iModule); iModule++;
2099        volIDs->AddAt(4225,iModule); iModule++;
2100        volIDs->AddAt(4226,iModule); iModule++;
2101        volIDs->AddAt(4227,iModule); iModule++;
2102        volIDs->AddAt(4228,iModule); iModule++;
2103        volIDs->AddAt(4229,iModule); iModule++;
2104        volIDs->AddAt(4230,iModule); iModule++;
2105        volIDs->AddAt(4231,iModule); iModule++;
2106        volIDs->AddAt(4232,iModule); iModule++;
2107        volIDs->AddAt(4233,iModule); iModule++;
2108        volIDs->AddAt(4234,iModule); iModule++;
2109        volIDs->AddAt(4235,iModule); iModule++;
2110        volIDs->AddAt(4236,iModule); iModule++;
2111        volIDs->AddAt(4237,iModule); iModule++;
2112        volIDs->AddAt(4238,iModule); iModule++;
2113        volIDs->AddAt(4239,iModule); iModule++;
2114      }
2115      if(sectors[9]==1){//--->cSect = 9 <//---
2116        volIDs->AddAt(2120,iModule); iModule++;
2117        volIDs->AddAt(2121,iModule); iModule++;
2118        volIDs->AddAt(2122,iModule); iModule++;
2119        volIDs->AddAt(2123,iModule); iModule++;
2120        volIDs->AddAt(2124,iModule); iModule++;
2121        volIDs->AddAt(2125,iModule); iModule++;
2122        volIDs->AddAt(2126,iModule); iModule++;
2123        volIDs->AddAt(2127,iModule); iModule++;
2124        volIDs->AddAt(4240,iModule); iModule++;
2125        volIDs->AddAt(4241,iModule); iModule++;
2126        volIDs->AddAt(4242,iModule); iModule++;
2127        volIDs->AddAt(4243,iModule); iModule++;
2128        volIDs->AddAt(4244,iModule); iModule++;
2129        volIDs->AddAt(4245,iModule); iModule++;
2130        volIDs->AddAt(4246,iModule); iModule++;
2131        volIDs->AddAt(4247,iModule); iModule++;
2132        volIDs->AddAt(4248,iModule); iModule++;
2133        volIDs->AddAt(4249,iModule); iModule++;
2134        volIDs->AddAt(4250,iModule); iModule++;
2135        volIDs->AddAt(4251,iModule); iModule++;
2136        volIDs->AddAt(4252,iModule); iModule++;
2137        volIDs->AddAt(4253,iModule); iModule++;
2138        volIDs->AddAt(4254,iModule); iModule++;
2139        volIDs->AddAt(4255,iModule); iModule++;
2140      }
2141
2142   return volIDs;
2143 }
2144
2145 //___________________________________
2146 TArrayI* AliITSRealignTracks::GetLayersVolUID(const Int_t *layer){
2147
2148   //return a TArrayI with the volUIDs of the modules into the set of layers
2149   //defined by layer[6]
2150
2151   TArrayI *out=new TArrayI(2198);
2152   Int_t last=0;
2153   UShort_t voluid;
2154   for(Int_t i=0;i<6;i++){
2155     if(layer[i]==1){
2156       for(Int_t mod=0;mod<AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);mod++){
2157         voluid=AliGeomManager::LayerToVolUID(i+AliGeomManager::kFirstLayer,mod);
2158         out->AddAt(voluid,last);
2159         //      printf("voluid %d at position %d \n",out->At(last),last);
2160         last++;
2161       }
2162     }
2163   }
2164   out->Set(last);
2165   return out;
2166 }
2167
2168 //_________________
2169 TArrayI* AliITSRealignTracks::SelectLayerInVolids(const TArrayI *volidsIN,AliGeomManager::ELayerID layer){
2170   //Select between the modules specified by their volUIDs in volidsIN only those
2171   // of a given layer "layer"
2172
2173   Int_t size=volidsIN->GetSize();
2174   Int_t count=0;
2175   for(Int_t j=0;j<size;j++){
2176     if(AliGeomManager::VolUIDToLayer(volidsIN->At(j))==layer)count++;
2177   }
2178   TArrayI *volidsOUT=new TArrayI(count);
2179   count=0;
2180   for(Int_t j=0;j<size;j++){
2181     if(AliGeomManager::VolUIDToLayer(volidsIN->At(j))==layer){
2182       volidsOUT->AddAt(volidsIN->At(j),count);
2183       count++;
2184     }
2185   }
2186   return volidsOUT;
2187 }
2188
2189 //______________________________________________
2190
2191 TArrayI* AliITSRealignTracks::IntersectVolArray(const TArrayI *vol1,const TArrayI *vol2){
2192
2193   //Perform the intersection between the array vol1 and vol2
2194
2195   Int_t size1=vol1->GetSize();
2196   Int_t size2=vol2->GetSize();
2197   Int_t last=0,volid;
2198   Bool_t found;
2199   TArrayI *volidOut=new TArrayI(size1+size2);
2200
2201   for(Int_t k=0;k<size1;k++){
2202     found=kFALSE;
2203     volid=vol1->At(k);
2204     for(Int_t j=0;j<size2;j++){
2205       if(vol2->At(j)==volid)found=kTRUE;
2206     }
2207     if(found){
2208       volidOut->AddAt(volid,last);
2209       last++;
2210     }
2211   }
2212   volidOut->Set(last);
2213   return volidOut;
2214 }
2215 //_________________________________________
2216
2217 TArrayI* AliITSRealignTracks::JoinVolArrays(const TArrayI *vol1,const TArrayI *vol2){
2218   //!BE CAREFUL: If an index is repeated into vol1 or into vol2 will be repeated also in the final array
2219
2220   Int_t size1=vol1->GetSize();
2221   Int_t size2=vol2->GetSize();
2222   Int_t count=0;
2223   UShort_t volid;
2224   Bool_t found;
2225   TArrayI *volidOut=new TArrayI(size1+size2);
2226
2227   for(Int_t k=0;k<size1;k++){
2228     volid=vol1->At(k);
2229     volidOut->AddAt(volid,k);
2230   }
2231
2232   for(Int_t k=0;k<size2;k++){
2233     found=kFALSE;
2234     volid=vol2->At(k);
2235     for(Int_t j=0;j<size1;j++){
2236       if(volidOut->At(j)==volid)found=kTRUE;
2237     }
2238     if(!found){
2239       volidOut->AddAt(volid,size1+count);
2240       count++;
2241     }
2242   }
2243   volidOut->Set(size1+count);
2244   return volidOut;
2245 }
2246
2247 //______________________________________
2248
2249 TArrayI* AliITSRealignTracks::ExcludeVolidsFromVolidsArray(const TArrayI *volidsToExclude,const TArrayI *volStart){
2250   //Excludes the modules defined by their volUID in the array volidsToExclude from the array volStart
2251
2252   Int_t size1=volidsToExclude->GetSize();
2253   Int_t size2=volStart->GetSize();
2254   Int_t last=0;
2255   UShort_t volid;
2256   Bool_t found;
2257   TArrayI *volidOut=new TArrayI(size2);
2258
2259   for(Int_t k=0;k<size2;k++){
2260     found=kFALSE;
2261     volid=volStart->At(k);
2262     for(Int_t j=0;j<size1;j++){
2263       if(volidsToExclude->At(j)==volid){
2264         found=kTRUE;
2265         break;
2266       }
2267     }
2268     if(!found){
2269       volidOut->AddAt(volid,last);
2270       last++;
2271     }
2272   }
2273   volidOut->Set(last);
2274   return volidOut;
2275 }
2276
2277
2278 //________________________________________
2279
2280 TArrayI* AliITSRealignTracks::GetLayerVolumes(const Int_t *layer){
2281   //returns a TArrayI with the volUIDs of the modules of the layers
2282   //specified into *layer
2283
2284   TArrayI *out=new TArrayI(2198);
2285   Int_t last=0;
2286   UShort_t voluid;
2287   for(Int_t i=0;i<6;i++){
2288     if(layer[i]==1){
2289       for(Int_t mod=0;mod<AliGeomManager::LayerSize(i+AliGeomManager::kFirstLayer);mod++){
2290         voluid=AliGeomManager::LayerToVolUID(i+AliGeomManager::kFirstLayer,mod);
2291         out->AddAt(voluid,last);
2292         //      printf("voluid %d at position %d \n",out->At(last),last);
2293         last++;
2294       }
2295     }
2296   }
2297   out->Set(last);
2298   return out;
2299 }
2300
2301
2302
2303 //______________________________
2304 TArrayI* AliITSRealignTracks::GetAlignedVolumes(char *filename){
2305   //Open the file "filename" which is expected to contain
2306   //a TClonesArray named "ITSAlignObjs" with stored a set of AlignObjs
2307   //returns an array with the volumes UID of the modules considered realigned
2308
2309   if(gSystem->AccessPathName(filename)){
2310     printf("Wrong Realignment file name \n");
2311     return 0x0;
2312   }
2313   TFile *f=TFile::Open(filename,"READ");
2314   TClonesArray *array=(TClonesArray*)f->Get("ITSAlignObjs");
2315   AliAlignObjParams *a;
2316   Int_t last=0;
2317   TArrayI *volidOut=new TArrayI(2200);
2318   for(Int_t j=0;j<array->GetSize();j++){
2319     a=(AliAlignObjParams*)array->At(j);
2320     if(a->GetUniqueID()==0)continue;
2321
2322     else {
2323       volidOut->AddAt(a->GetVolUID(),last);
2324       last++;
2325     }
2326   }
2327   volidOut->Set(last);
2328   f->Close();
2329   return volidOut;
2330 }
2331
2332
2333 //________________________________________
2334 void AliITSRealignTracks::SetDraw(Bool_t draw,Bool_t refresh){
2335   //TEPMORARY METHOD: method to switch on/off the drawing of histograms
2336   // if refresh=kTRUE deletes the old histos and constructs new ones
2337
2338   if(refresh){
2339     // WriteHists();
2340     if(fAlignDrawObjs)DeleteDrawHists();
2341     InitDrawHists();
2342   }
2343   fDraw=draw;
2344   return;
2345 }
2346
2347 void AliITSRealignTracks::DeleteDrawHists(){
2348   //Delete the pointers to the histograms
2349
2350   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
2351     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
2352       delete fAlignDrawObjs[iLayer][iModule];
2353     }
2354     if(fAlignDrawObjs[iLayer])delete [] fAlignDrawObjs[iLayer];
2355   }
2356
2357   delete [] fAlignDrawObjs;
2358   fAlignDrawObjs = 0;
2359
2360
2361   delete fCanvPar;
2362   delete fCanvGr;
2363   delete fgrIterMeanX;
2364   delete fgrIterRMSX;
2365   delete fgrIterMeanY;
2366   delete fgrIterRMSY;
2367   delete fgrIterMeanZ;
2368   delete fgrIterRMSZ;
2369   delete fgrIterMeanPsi;
2370   delete fgrIterRMSPsi;
2371   delete fgrIterMeanTheta;
2372   delete fgrIterRMSTheta;
2373   delete fgrIterMeanPhi;
2374   delete fgrIterRMSPhi;
2375
2376 }
2377
2378 void AliITSRealignTracks::InitDrawHists(){
2379   //Initialize the histograms to monitor the results
2380
2381   Int_t nLayers = AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer;
2382   fAlignDrawObjs = new AliAlignObj**[nLayers];
2383   for (Int_t iLayer = 0; iLayer < (AliGeomManager::kLastLayer - AliGeomManager::kFirstLayer); iLayer++) {
2384     fAlignDrawObjs[iLayer] = new AliAlignObj*[AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer)];
2385     for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
2386       UShort_t volid = AliGeomManager::LayerToVolUID(iLayer+ AliGeomManager::kFirstLayer,iModule);
2387       fAlignDrawObjs[iLayer][iModule] = new AliAlignObjParams(AliGeomManager::SymName(volid),volid,0,0,0,0,0,0,kTRUE);
2388       fAlignDrawObjs[iLayer][iModule]->SetUniqueID(1);
2389     }
2390   }
2391
2392
2393   TH1F *hX=new TH1F("hX","hX",1000,-10000.,10000.);
2394   TH1F *hY=new TH1F("hY","hY",1000,-10000.,10000.);
2395   TH1F *hZ=new TH1F("hZ","hZ",1000,-10000.,10000.);
2396   TH1F *hPsi=new TH1F("hPsi","hPsi",1000,-5000.,5000.);
2397   TH1F *hTheta=new TH1F("hTheta","hTheta",1000,-5000.,5000.);
2398   TH1F *hPhi=new TH1F("hPhi","hPhi",1000,-5000.,5000.);
2399
2400   fCanvPar=new TCanvas("fCanvPar","Parameters trend during iterations: Convergence \n");
2401   fCanvPar->Divide(3,2);
2402   fCanvPar->cd(1);
2403   hX->Draw();
2404   hX->SetXTitle("#mum");
2405   fCanvPar->cd(2);
2406   hY->Draw();
2407   hY->SetXTitle("#mum");
2408   fCanvPar->cd(3);
2409   hZ->SetXTitle("#mum");
2410   hZ->Draw();
2411   fCanvPar->cd(4);
2412   hPsi->SetXTitle("mdeg");
2413   hPsi->Draw();
2414   fCanvPar->cd(5);
2415   hTheta->SetXTitle("mdeg");
2416   hTheta->Draw();
2417   fCanvPar->cd(6);
2418   hPhi->SetXTitle("mdeg");
2419   hPhi->Draw();
2420   fCanvPar->Update();
2421
2422
2423   fCanvGr=new TCanvas("fCanvGr","Parameters trend during iterations: Convergence \n");
2424   fCanvGr->Divide(3,2);
2425
2426   fCanvGr->cd(1);
2427   fgrIterMeanX=new TGraph(1);
2428   fgrIterRMSX=new TGraph(1);
2429   fgrIterRMSX->GetYaxis()->SetRangeUser(-1000.,1000.);
2430   fgrIterRMSX->SetName("fgrIterRMSX");
2431   fgrIterRMSX->SetLineColor(2);
2432   fgrIterMeanX->SetName("fgrIterMeanX");
2433   fgrIterMeanX->SetTitle("Convergence of #deltaX \n");
2434   fgrIterMeanX->GetXaxis()->SetTitle("#mum");
2435   fgrIterRMSX->Draw("acp");
2436   fgrIterMeanX->Draw("cp");
2437
2438   fCanvGr->cd(2);
2439   fgrIterMeanY=new TGraph(1);
2440   fgrIterRMSY=new TGraph(1);
2441   fgrIterRMSY->GetYaxis()->SetRangeUser(-1000.,1000.);
2442   fgrIterRMSY->SetName("fgrIterRMSY");
2443   fgrIterRMSY->SetLineColor(2);
2444   fgrIterMeanY->SetName("fgrIterMeanY");
2445   fgrIterMeanY->SetTitle("Convergence of #deltaY \n");
2446   fgrIterMeanY->GetXaxis()->SetTitle("#mum");
2447   fgrIterRMSY->Draw("acp");
2448   fgrIterMeanY->Draw("cp");
2449
2450   fCanvGr->cd(3);
2451   fgrIterMeanZ=new TGraph(1);
2452   fgrIterRMSZ=new TGraph(1);
2453   fgrIterRMSZ->GetYaxis()->SetRangeUser(-1000.,1000.);
2454   fgrIterRMSZ->SetName("fgrIterRMSZ");
2455   fgrIterRMSZ->SetLineColor(2);
2456   fgrIterMeanZ->SetName("fgrIterMeanZ");
2457   fgrIterMeanZ->SetTitle("Convergence of #deltaZ \n");
2458   fgrIterMeanZ->GetXaxis()->SetTitle("#mum");
2459   fgrIterRMSZ->Draw("acp");
2460   fgrIterMeanZ->Draw("cp");
2461
2462   fCanvGr->cd(4);
2463   fgrIterMeanPsi=new TGraph(1);
2464   fgrIterRMSPsi=new TGraph(1);
2465   fgrIterRMSPsi->GetYaxis()->SetRangeUser(-1000.,1000.);
2466   fgrIterRMSPsi->SetName("fgrIterRMSPsi");
2467   fgrIterRMSPsi->SetLineColor(2);
2468   fgrIterMeanPsi->SetName("fgrIterMeanPsi");
2469   fgrIterMeanPsi->SetTitle("Convergence of #deltaPsi \n");
2470   fgrIterMeanPsi->GetXaxis()->SetTitle("mdeg");
2471   fgrIterRMSPsi->Draw("acp");
2472   fgrIterMeanPsi->Draw("cp");
2473
2474   fCanvGr->cd(5);
2475   fgrIterMeanTheta=new TGraph(1);
2476   fgrIterRMSTheta=new TGraph(1);
2477   fgrIterRMSTheta->GetYaxis()->SetRangeUser(-1000.,1000.);
2478   fgrIterRMSTheta->SetName("fgrIterRMSTheta");
2479   fgrIterRMSTheta->SetLineColor(2);
2480   fgrIterMeanTheta->SetName("fgrIterMeanTheta");
2481   fgrIterMeanTheta->SetTitle("Convergence of #deltaTheta \n");
2482   fgrIterMeanTheta->GetXaxis()->SetTitle("mdeg");
2483   fgrIterRMSTheta->Draw("acp");
2484   fgrIterMeanTheta->Draw("cp");
2485
2486   fCanvGr->cd(6);
2487   fgrIterMeanPhi=new TGraph(1);
2488   fgrIterRMSPhi=new TGraph(1);
2489   fgrIterRMSPhi->GetYaxis()->SetRangeUser(-1000.,1000.);
2490   fgrIterRMSPhi->SetName("fgrIterRMSPhi");
2491   fgrIterRMSPhi->SetLineColor(2);
2492   fgrIterMeanPhi->SetName("fgrIterMeanPhi");
2493   fgrIterMeanPhi->SetTitle("Convergence of #deltaPhi \n");
2494   fgrIterMeanPhi->GetXaxis()->SetTitle("mdeg");
2495   fgrIterRMSPhi->Draw("acp");
2496   fgrIterMeanPhi->Draw("cp");
2497
2498
2499
2500 }
2501
2502 void AliITSRealignTracks::UpdateDraw(TArrayI *volids,Int_t iter,Int_t color){
2503   //Updates the histograms to monitor the results. Only the histograms
2504   //of the volumes specified in *volids will be updated.
2505   // iter is just a flag for the names of the histo
2506   // color specifies the color of the lines of the histograms for this update
2507
2508   TString name="hX_";
2509   name+=iter;
2510   name.Append("iter");
2511   TH1F *hX=new TH1F("hX",name.Data(),1000,-10000.,10000.);
2512
2513   name="hY_";
2514   name+=iter;
2515   name.Append("iter");
2516   TH1F *hY=new TH1F("hY",name.Data(),1000,-10000.,10000.);
2517
2518   name="hZ_";
2519   name+=iter;
2520   name.Append("iter");
2521   TH1F *hZ=new TH1F("hZ",name.Data(),1000,-10000.,10000.);
2522
2523   name="hPsi_";
2524   name+=iter;
2525   name.Append("iter");
2526   TH1F *hPsi=new TH1F("hPsi",name.Data(),1000,-5000.,5000.);
2527
2528   name="hTheta_";
2529   name+=iter;
2530   name.Append("iter");
2531   TH1F *hTheta=new TH1F("hTheta",name.Data(),1000,-5000.,5000.);
2532
2533   name="hPhi_";
2534   name+=iter;
2535   name.Append("iter");
2536   TH1F *hPhi=new TH1F("hPhi",name.Data(),1000,-5000.,5000.);
2537
2538   Int_t layer,mod;
2539   Double_t transl[3],rot[3],transldr[3],rotdr[3];
2540
2541   for(Int_t i=0;i<volids->GetSize();i++){
2542     layer=AliGeomManager::VolUIDToLayer(volids->At(i),mod);
2543     fAlignObjs[layer-AliGeomManager::kFirstLayer][mod]->GetPars(transl,rot);
2544     fAlignDrawObjs[layer-AliGeomManager::kFirstLayer][mod]->GetPars(transldr,rotdr);
2545
2546     hX->Fill(10000.*(transl[0]-transldr[0]));
2547     hY->Fill(10000.*(transl[1]-transldr[1]));
2548     hZ->Fill(10000.*(transl[2]-transldr[2]));
2549     hPsi->Fill(1000.*(rot[0]-rotdr[0]));
2550     hTheta->Fill(1000.*(rot[1]-rotdr[1]));
2551     hPhi->Fill(1000.*(rot[1]-rotdr[2]));
2552     //Update the pars of the draw object
2553     fAlignDrawObjs[layer-AliGeomManager::kFirstLayer][mod]->SetPars(transl[0],transl[1],transl[2],rot[0],rot[1],rot[2]);
2554   }
2555
2556   hX->SetLineColor(color);
2557   hY->SetLineColor(color);
2558   hZ->SetLineColor(color);
2559   hPsi->SetLineColor(color);
2560   hTheta->SetLineColor(color);
2561   hPhi->SetLineColor(color);
2562
2563
2564   fCanvPar->cd(1);
2565   hX->Draw("Same");
2566   fCanvPar->cd(2);
2567   hY->Draw("Same");
2568   fCanvPar->cd(3);
2569   hZ->Draw("Same");
2570   fCanvPar->cd(4);
2571   hPsi->Draw("Same");
2572   fCanvPar->cd(5);
2573   hTheta->Draw("Same");
2574   fCanvPar->cd(6);
2575   hPhi->Draw("Same");
2576   gPad->Modified();
2577   fCanvPar->Update();
2578   fCanvPar->Modified();
2579
2580   fgrIterMeanX->SetPoint(fgrIterMeanX->GetN()+1,iter,hX->GetMean());
2581   fgrIterRMSX->SetPoint(fgrIterRMSX->GetN()+1,iter,hX->GetRMS());
2582   fgrIterMeanY->SetPoint(fgrIterMeanY->GetN()+1,iter,hY->GetMean());
2583   fgrIterRMSY->SetPoint(fgrIterRMSY->GetN()+1,iter,hY->GetRMS());
2584   fgrIterMeanZ->SetPoint(fgrIterMeanZ->GetN()+1,iter,hZ->GetMean());
2585   fgrIterRMSZ->SetPoint(fgrIterRMSZ->GetN()+1,iter,hZ->GetRMS());
2586   fgrIterMeanPsi->SetPoint(fgrIterMeanPsi->GetN()+1,iter,hPsi->GetMean());
2587   fgrIterRMSPsi->SetPoint(fgrIterRMSPsi->GetN()+1,iter,hPsi->GetRMS());
2588   fgrIterMeanTheta->SetPoint(fgrIterMeanTheta->GetN()+1,iter,hTheta->GetMean());
2589   fgrIterRMSTheta->SetPoint(fgrIterRMSTheta->GetN()+1,iter,hTheta->GetRMS());
2590   fgrIterMeanPhi->SetPoint(fgrIterMeanPhi->GetN()+1,iter,hPhi->GetMean());
2591   fgrIterRMSPhi->SetPoint(fgrIterRMSPhi->GetN()+1,iter,hPhi->GetRMS());
2592
2593   gPad->Modified();
2594   fCanvGr->Update();
2595   fCanvGr->Update();
2596 }
2597
2598 void AliITSRealignTracks::WriteHists(const char *outfile){
2599   //Writes the histograms for the monitoring of the results
2600   // in a file named "outfile"
2601
2602   TFile *f=new TFile(outfile,"RECREATE");
2603   f->cd();
2604   fCanvPar->Write();
2605   fCanvGr->Write();
2606   fgrIterMeanX->Write();
2607   fgrIterRMSX->Write();
2608   fgrIterMeanY->Write();
2609   fgrIterRMSY->Write();
2610   fgrIterMeanZ->Write();
2611   fgrIterRMSZ->Write();
2612   fgrIterMeanPsi->Write();
2613   fgrIterRMSPsi->Write();
2614   fgrIterMeanTheta->Write();
2615   fgrIterRMSTheta->Write();
2616   fgrIterMeanPhi->Write();
2617   fgrIterRMSPhi->Write();
2618
2619   f->Close();
2620   return;
2621
2622 }