]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliAlignmentTracks.cxx
Setter added.
[u/mrichter/AliRoot.git] / STEER / AliAlignmentTracks.cxx
index cd8efdac45e422d7b73b1ae2c24d7464f69e5734..e734e375bc5db3e31e3de9e684e9a62a2abb1e90 100644 (file)
@@ -25,6 +25,7 @@
 #include <TChain.h>
 #include <TFile.h>
 #include <TVector3.h>
+#include <TSystem.h>
 
 #include "AliAlignmentTracks.h"
 #include "AliTrackPointArray.h"
@@ -49,7 +50,8 @@ AliAlignmentTracks::AliAlignmentTracks():
   fMisalignObjs(0),
   fTrackFitter(0),
   fMinimizer(0),
-  fDoUpdate(kTRUE)
+  fDoUpdate(kTRUE),
+  fCovIsUsed(kFALSE)
 {
   // Default constructor
   InitIndex();
@@ -69,7 +71,8 @@ AliAlignmentTracks::AliAlignmentTracks(TChain *esdchain):
   fMisalignObjs(0),
   fTrackFitter(0),
   fMinimizer(0),
-  fDoUpdate(kTRUE)
+  fDoUpdate(kTRUE),
+  fCovIsUsed(kFALSE)
 {
   // Constructor in the case
   // the user provides an already
@@ -92,7 +95,8 @@ AliAlignmentTracks::AliAlignmentTracks(const char *esdfilename, const char *esdt
   fMisalignObjs(0),
   fTrackFitter(0),
   fMinimizer(0),
-  fDoUpdate(kTRUE)
+  fDoUpdate(kTRUE),
+  fCovIsUsed(kFALSE)
 {
   // Constructor in the case
   // the user provides a single ESD file
@@ -143,6 +147,7 @@ void AliAlignmentTracks::AddESD(const char *esdfilename, const char *esdtreename
   }
 }
 
+
 //________________________________________________________________________
 void AliAlignmentTracks::ProcessESD(Bool_t onlyITS,
                                    Int_t minITSpts,
@@ -165,7 +170,7 @@ void AliAlignmentTracks::ProcessESD(Bool_t onlyITS,
   fESDChain->SetBranchAddress("ESDfriend.",&esdf);
 
   // Open the output file
-  if (fPointsFilename.Data() == "") {
+  if (fPointsFilename.IsNull()) {
     AliWarning("Incorrect output filename!");
     return;
   }
@@ -209,31 +214,41 @@ void AliAlignmentTracks::ProcessESD(Bool_t onlyITS,
       array = track->GetTrackPointArray();
 
       if(onlyITS) {
-       array2 = new AliTrackPointArray(track->GetNcls(0));
-       Bool_t smallAngleWrtITSModulePlanes=kFALSE;
+       Bool_t layerOK[6]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
        Int_t ipt,volId,modId,layerId;
        Int_t jpt=0;
        for(ipt=0; ipt<array->GetNPoints(); ipt++) {
          array->GetPoint(point,ipt);
          volId = point.GetVolumeID();
          layerId = AliGeomManager::VolUIDToLayer(volId,modId);
-         if(layerId<7) {
-           array2->AddPoint(jpt,&point);
-           jpt++;
-         }
+         if(layerId>6) continue;
          // check minAngleWrtITSModulePlanes
          if(cuts) {
            Double_t p[3]; track->GetDirection(p);
            TVector3 pvec(p[0],p[1],p[2]);
            Double_t rot[9]; AliGeomManager::GetOrigRotation(volId,rot);
-           TVector3 normvec(rot[1],-rot[0],0.);
+           TVector3 normvec(rot[1],rot[4],rot[7]);
            Double_t angle = pvec.Angle(normvec);
            if(angle>0.5*TMath::Pi()) angle = TMath::Pi()-angle;
            angle = 0.5*TMath::Pi()-angle;
-           if(angle<minAngleWrtITSModulePlanes) smallAngleWrtITSModulePlanes=kTRUE;
+           if(angle<minAngleWrtITSModulePlanes) {
+             layerOK[layerId-1]=kFALSE;
+             continue;
+           }
          }
+         jpt++;
+       }
+       if(jpt < minITSpts) continue;      
+       array2 = new AliTrackPointArray(jpt);
+       jpt=0;
+       for(ipt=0; ipt<array->GetNPoints(); ipt++) {
+         array->GetPoint(point,ipt);
+         volId = point.GetVolumeID();
+         layerId = AliGeomManager::VolUIDToLayer(volId,modId);
+         if(layerId>6 || !layerOK[layerId-1]) continue;
+         array2->AddPoint(jpt,&point);
+         jpt++;
        }
-       if(smallAngleWrtITSModulePlanes) continue;
       } // end if(onlyITS)
  
       pointsTree->Fill();
@@ -273,7 +288,7 @@ void AliAlignmentTracks::ProcessESDCosmics(Bool_t onlyITS,
   fESDChain->SetBranchAddress("ESDfriend.",&esdf);
 
   // Open the output file
-  if (fPointsFilename.Data() == "") {
+  if (fPointsFilename.IsNull()) {
     AliWarning("Incorrect output filename!");
     return;
   }
@@ -355,18 +370,11 @@ void AliAlignmentTracks::ProcessESDCosmics(Bool_t onlyITS,
     AliESDtrack * track1 = esd->GetTrack(good1);
     AliESDtrack * track2 = esd->GetTrack(good2);
 
-    Int_t ntotpts;
-    if(onlyITS) {
-      ntotpts = track1->GetNcls(0)+track2->GetNcls(0);
-    } else {
-      ntotpts = track1->GetTrackPointArray()->GetNPoints()+
-               track2->GetTrackPointArray()->GetNPoints();
-    }
-    array2 = new AliTrackPointArray(ntotpts);
     AliTrackPoint point;
     Int_t ipt,volId,modId,layerId;
     Int_t jpt=0;
-    Bool_t smallAngleWrtITSModulePlanes=kFALSE;
+    Bool_t layerOK[6][2]; 
+    for(Int_t l1=0;l1<6;l1++) for(Int_t l2=0;l2<2;l2++) layerOK[l1][l2]=kTRUE; 
     array = track1->GetTrackPointArray();
     for(ipt=0; ipt<array->GetNPoints(); ipt++) {
       array->GetPoint(point,ipt);
@@ -379,14 +387,16 @@ void AliAlignmentTracks::ProcessESDCosmics(Bool_t onlyITS,
          Double_t p[3]; track1->GetDirection(p);
          TVector3 pvec(p[0],p[1],p[2]);
          Double_t rot[9]; AliGeomManager::GetOrigRotation(volId,rot);
-         TVector3 normvec(rot[1],-rot[0],0.);
+         TVector3 normvec(rot[1],rot[4],rot[7]);
          Double_t angle = pvec.Angle(normvec);
          if(angle>0.5*TMath::Pi()) angle = TMath::Pi()-angle;
          angle = 0.5*TMath::Pi()-angle;
-         if(angle<minAngleWrtITSModulePlanes) smallAngleWrtITSModulePlanes=kTRUE;
+         if(angle<minAngleWrtITSModulePlanes) {
+           layerOK[layerId-1][0]=kFALSE;
+           continue;
+         }
        }
       }
-      array2->AddPoint(jpt,&point);
       jpt++;
     }
     array = track2->GetTrackPointArray();
@@ -401,18 +411,45 @@ void AliAlignmentTracks::ProcessESDCosmics(Bool_t onlyITS,
          Double_t p[3]; track2->GetDirection(p);
          TVector3 pvec(p[0],p[1],p[2]);
          Double_t rot[9]; AliGeomManager::GetOrigRotation(volId,rot);
-         TVector3 normvec(rot[1],-rot[0],0.);
+         TVector3 normvec(rot[1],rot[4],rot[7]);
          Double_t angle = pvec.Angle(normvec);
          if(angle>0.5*TMath::Pi()) angle = TMath::Pi()-angle;
          angle = 0.5*TMath::Pi()-angle;
-         if(angle<minAngleWrtITSModulePlanes) smallAngleWrtITSModulePlanes=kTRUE;
+         if(angle<minAngleWrtITSModulePlanes) {
+           layerOK[layerId-1][0]=kFALSE;
+           continue;
+         }
        }
       }
+      jpt++;
+    }
+
+    if(jpt < 2*minITSpts) continue;
+    array2 = new AliTrackPointArray(jpt);
+    jpt=0;
+    array = track1->GetTrackPointArray();
+    for(ipt=0; ipt<array->GetNPoints(); ipt++) {
+      array->GetPoint(point,ipt);
+      if(onlyITS) {
+       volId = point.GetVolumeID();
+       layerId = AliGeomManager::VolUIDToLayer(volId,modId);
+       if(layerId>6 || !layerOK[layerId-1][0]) continue;
+      }
+      array2->AddPoint(jpt,&point);
+      jpt++;
+    }
+    array = track2->GetTrackPointArray();
+    for(ipt=0; ipt<array->GetNPoints(); ipt++) {
+      array->GetPoint(point,ipt);
+      if(onlyITS) {
+       volId = point.GetVolumeID();
+       layerId = AliGeomManager::VolUIDToLayer(volId,modId);
+       if(layerId>6 || !layerOK[layerId-1][1]) continue;
+      }
       array2->AddPoint(jpt,&point);
       jpt++;
     }
 
-    if(smallAngleWrtITSModulePlanes) continue;
     pointsTree->Fill();
   }
 
@@ -556,11 +593,37 @@ void AliAlignmentTracks::DeleteIndex()
 //______________________________________________________________________________
 Bool_t AliAlignmentTracks::ReadAlignObjs(const char *alignObjFileName, const char* arrayName)
 {
-  // Read alignment object from a file
+  // Read alignment object from a file: update the alignobj already present with the one in the file
   // To be replaced by a call to CDB
-  AliWarning(Form("Method not yet implemented (%s in %s) !",arrayName,alignObjFileName));
+  
+  if(gSystem->AccessPathName(alignObjFileName,kFileExists)){
+    printf("Wrong AlignObjs File Name \n");
+    return kFALSE;
+  } 
 
-  return kFALSE;
+  TFile *fRealign=TFile::Open(alignObjFileName);
+  if (!fRealign || !fRealign->IsOpen()) {
+    AliError(Form("Could not open Align Obj File file %s !",alignObjFileName));
+    return kFALSE;
+  }  
+  printf("Getting TClonesArray \n");
+  TClonesArray *clnarray=(TClonesArray*)fRealign->Get(arrayName);
+  Int_t size=clnarray->GetSize();
+  UShort_t volid;
+
+  for(Int_t ivol=0;ivol<size;ivol++){
+    AliAlignObjParams *a=(AliAlignObjParams*)clnarray->At(ivol);
+    volid=a->GetVolUID();
+    Int_t iModule;
+    AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
+    if(iLayer<AliGeomManager::kFirstLayer||iLayer>AliGeomManager::kSSD2)continue;
+    printf("Updating volume: %d ,layer: %d module: %d \n",volid,iLayer,iModule);
+    *fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule] *= *a;
+  }
+  delete clnarray;
+  fRealign->Close();
+  return kTRUE;
 }
 
 //______________________________________________________________________________
@@ -602,11 +665,11 @@ void AliAlignmentTracks::DeleteAlignObjs()
   fAlignObjs = 0;
 }
 
-void AliAlignmentTracks::AlignDetector(AliGeomManager::ELayerID firstLayer,
-                                      AliGeomManager::ELayerID lastLayer,
-                                      AliGeomManager::ELayerID layerRangeMin,
-                                      AliGeomManager::ELayerID layerRangeMax,
-                                      Int_t iterations)
+Bool_t AliAlignmentTracks::AlignDetector(AliGeomManager::ELayerID firstLayer,
+                                        AliGeomManager::ELayerID lastLayer,
+                                        AliGeomManager::ELayerID layerRangeMin,
+                                        AliGeomManager::ELayerID layerRangeMax,
+                                        Int_t iterations)
 {
   // Align detector volumes within
   // a given layer range
@@ -627,17 +690,19 @@ void AliAlignmentTracks::AlignDetector(AliGeomManager::ELayerID firstLayer,
     }
   }
 
+  Bool_t result = kFALSE;
   while (iterations > 0) {
-    AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
+    if (!(result = AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax))) break;
     iterations--;
   }
+  return result;
 }
 
 //______________________________________________________________________________
-void AliAlignmentTracks::AlignLayer(AliGeomManager::ELayerID layer,
-                                   AliGeomManager::ELayerID layerRangeMin,
-                                   AliGeomManager::ELayerID layerRangeMax,
-                                   Int_t iterations)
+Bool_t AliAlignmentTracks::AlignLayer(AliGeomManager::ELayerID layer,
+                                     AliGeomManager::ELayerID layerRangeMin,
+                                     AliGeomManager::ELayerID layerRangeMax,
+                                     Int_t iterations)
 {
   // Align detector volumes within
   // a given layer.
@@ -650,14 +715,16 @@ void AliAlignmentTracks::AlignLayer(AliGeomManager::ELayerID layer,
     volIds.AddAt(volId,iModule);
   }
 
+  Bool_t result = kFALSE;
   while (iterations > 0) {
-    AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax);
+    if (!(result = AlignVolumes(&volIds,0x0,layerRangeMin,layerRangeMax))) break;
     iterations--;
   }
+  return result;
 }
 
 //______________________________________________________________________________
-void AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit,
+Bool_t AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit,
                                     Int_t iterations)
 {
   // Align single detector volume to
@@ -669,14 +736,16 @@ void AliAlignmentTracks::AlignVolume(UShort_t volId, UShort_t volIdFit,
   TArrayI volIdsFit(1);
   volIdsFit.AddAt(volIdFit,0);
 
+  Bool_t result = kFALSE;
   while (iterations > 0) {
-    AlignVolumes(&volIds,&volIdsFit);
+    if (!(result = AlignVolumes(&volIds,&volIdsFit))) break;
     iterations--;
   }
+  return result;
 }
 
 //______________________________________________________________________________
-void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *volidsfit,
+Bool_t AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *volidsfit,
                                     AliGeomManager::ELayerID layerRangeMin,
                                     AliGeomManager::ELayerID layerRangeMax,
                                     Int_t iterations)
@@ -691,17 +760,19 @@ void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *voli
   Int_t nVolIds = volids->GetSize();
   if (nVolIds == 0) {
     AliError("Volume IDs array is empty!");
-    return;
+    return kFALSE;
   }
-
   // Load only the tracks with at least one
   // space point in the set of volume (volids)
   BuildIndex();
   AliTrackPointArray **points;
+  Int_t pointsdim;
   // Start the iterations
+  Bool_t result = kFALSE;
   while (iterations > 0) {
-    Int_t nArrays = LoadPoints(volids, points);
-    if (nArrays == 0) return;
+    Int_t nArrays = LoadPoints(volids, points,pointsdim);
+    if (nArrays == 0) return kFALSE;
 
     AliTrackResiduals *minimizer = CreateMinimizer();
     minimizer->SetNTracks(nArrays);
@@ -715,7 +786,7 @@ void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *voli
       fitter->GetTrackResiduals(pVolId,pTrack);
       minimizer->AddTrackPointArrays(pVolId,pTrack);
     }
-    minimizer->Minimize();
+    if (!(result = minimizer->Minimize())) break;
 
     // Update the alignment object(s)
     if (fDoUpdate) for (Int_t iVolId = 0; iVolId < nVolIds; iVolId++) {
@@ -724,17 +795,18 @@ void AliAlignmentTracks::AlignVolumes(const TArrayI *volids, const TArrayI *voli
       AliGeomManager::ELayerID iLayer = AliGeomManager::VolUIDToLayer(volid,iModule);
       AliAlignObj *alignObj = fAlignObjs[iLayer-AliGeomManager::kFirstLayer][iModule];      
       *alignObj *= *minimizer->GetAlignObj();
-      alignObj->Print("");
+      if(iterations==1)alignObj->Print("");
     }
 
-    UnloadPoints(nArrays, points);
+    UnloadPoints(pointsdim, points);
     
     iterations--;
   }
+  return result;
 }
   
 //______________________________________________________________________________
-Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray** &points)
+Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray** &points,Int_t &pointsdim)
 {
   // Load track point arrays with at least
   // one space point in a given set of detector
@@ -771,7 +843,7 @@ Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray**
 
   if (nArrays == 0) {
     AliError("There are no space-points belonging to all of the volumes which are to be aligned!");
-    points = 0;
+    points = 0x0;
     return 0;
   }
 
@@ -779,6 +851,7 @@ Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray**
   fPointsTree->SetBranchAddress("SP", &array);
 
   // Allocate the pointer to the space-point arrays
+  pointsdim=nArrays;
   points = new AliTrackPointArray*[nArrays];
   for (Int_t i = 0; i < nArrays; i++) points[i] = 0x0;
 
@@ -802,7 +875,10 @@ Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray**
 
       // Get tree entry
       Int_t entry = (*index)[iArrayId];
-      if (indexUsed[entry] == kTRUE) continue;
+      if (indexUsed[entry] == kTRUE) {
+       nArrays--;
+       continue;
+      }
       fPointsTree->GetEvent(entry);
       if (!array) {
        AliWarning("Wrong space point array index!");
@@ -841,8 +917,21 @@ Int_t AliAlignmentTracks::LoadPoints(const TArrayI *volids, AliTrackPointArray**
        }
        // End of misalignment
 
+
        AliAlignObj *alignObj = fAlignObjs[layer-AliGeomManager::kFirstLayer][modnum];
-       alignObj->Transform(p);
+       UShort_t volp=p.GetVolumeID();
+       Bool_t found=kFALSE;
+       if(fCovIsUsed){
+         for (Int_t iVol = 0; iVol < nVolIds; iVol++) {
+           UShort_t vol = (*volids)[iVol];
+           if(volp==vol){
+             alignObj->Transform(p,kFALSE);
+             found=kTRUE;
+             break;
+           }
+         }
+       }
+       if(!found)alignObj->Transform(p,fCovIsUsed);
        points[iArray]->AddPoint(iPoint,&p);
       }
       iArray++;
@@ -932,7 +1021,38 @@ Bool_t AliAlignmentTracks::Misalign(const char *misalignObjFileName, const char*
     {
       AliAlignObj* alObj = (AliAlignObj*)array->UncheckedAt(i);
       alObj->GetVolUID(layerId,modId);
+      if(layerId<AliGeomManager::kFirstLayer) {
+       AliWarning(Form("Alignment object is ignored: %s",alObj->GetSymName()));
+       continue;
+      }
       fMisalignObjs[layerId-AliGeomManager::kFirstLayer][modId] = alObj;
     }
   return kTRUE;
 }
+
+
+//________________________________________________
+void AliAlignmentTracks::WriteRealignObjArray(TString outfilename,AliGeomManager::ELayerID layerRangeMin,AliGeomManager::ELayerID layerRangeMax){
+  
+  Int_t last=0;
+  TClonesArray *clonesArray=new TClonesArray("AliAlignObjParams",2200);
+  TClonesArray &alo=*clonesArray;
+  for (Int_t iLayer = layerRangeMin-AliGeomManager::kFirstLayer; iLayer <= (layerRangeMax - AliGeomManager::kFirstLayer);iLayer++) {
+  
+    for (Int_t iModule = 0; iModule < AliGeomManager::LayerSize(iLayer + AliGeomManager::kFirstLayer); iModule++) {
+     
+      AliAlignObj *alignObj = fAlignObjs[iLayer][iModule]; 
+      new(alo[last])AliAlignObjParams(*alignObj);
+      last++;
+    }
+  }
+  TFile *file=new TFile(outfilename.Data(),"RECREATE");
+  file->cd();
+  alo.Write("ITSAlignObjs",TObject::kSingleKey);
+  file->Close();    
+     
+  delete clonesArray;
+  return;
+}