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