]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - MUON/MUONAlignment.C
Protection if the compiler is not gcc
[u/mrichter/AliRoot.git] / MUON / MUONAlignment.C
index 96bd50e9f149cdbf09cf28045ae46fc024c3a540..40af52d6ca8e52f5325b5eb3a2a829071d7dd451 100644 (file)
 // Authors: B. Becker and J. Castillo
 // ---
 
-void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TString fileList = "")
+#if !defined(__CINT__) || defined(__MAKECINT__)
+
+#include "AliMUONAlignment.h"
+#include "AliMUONTrack.h"
+#include "AliMUONTrackExtrap.h"
+#include "AliMUONTrackParam.h"
+#include "AliMUONGeometryTransformer.h"
+#include "AliMUONDataInterface.h"
+#include "AliMUONVTrackStore.h"
+
+#include "AliMagFMaps.h"
+#include "AliTracker.h"
+#include "AliCDBManager.h"
+#include "AliCDBMetaData.h"
+#include "AliCDBId.h"
+#include "AliGeomManager.h"
+
+#include <TString.h>
+#include <TError.h>
+#include <TH1.h>
+#include <TGraphErrors.h>
+#include <TFile.h>
+#include <TClonesArray.h>
+#include <Riostream.h>
+
+#include <fstream>
+
+#endif
+
+void MUONAlignment(Int_t nEvents = 100000, char* geoFilename = "geometry.root", TString fileName = "galice.root", TString fileList = "")
 {
  
+  // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
+  if ( ! AliGeomManager::GetGeometry() ) {
+    AliGeomManager::LoadGeometry(geoFilename);
+    if (! AliGeomManager::GetGeometry() ) {
+      Error("MUONAlignment", "getting geometry from file %s failed", geoFilename);
+      return;
+    }
+  }
+  
+  // set  mag field 
+  // waiting for mag field in CDB 
+  printf("Loading field map...\n");
+  AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 1, 1., 10., AliMagFMaps::k5kG);
+  AliTracker::SetFieldMap(field, kFALSE);
+  // set the magnetic field for track extrapolations
+  AliMUONTrackExtrap::SetField(AliTracker::GetFieldMap());
+
   Double_t parameters[3*156];
   Double_t errors[3*156];
   Double_t pulls[3*156];
@@ -48,10 +94,10 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
 
   Double_t trackParams[8] = {0.,0.,0.,0.,0.,0.,0.,0.};
 
-  Int_t iPar = 0;
   // Set initial values here, good guess may help convergence
   // St 1 
-  //   parameters[iPar++] =  0.010300 ;  parameters[iPar++] =  0.010600 ;  parameters[iPar++] =  0.000396 ;  
+  //  Int_t iPar = 0;
+  //  parameters[iPar++] =  0.010300 ;  parameters[iPar++] =  0.010600 ;  parameters[iPar++] =  0.000396 ;  
 
   bool bLoop = kFALSE;
   ifstream sFileList;
@@ -71,16 +117,43 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
   AliMUONAlignment* alig = new AliMUONAlignment();
   alig->InitGlobalParameters(parameters);
 
-  AliMUONGeometryTransformer *transform = new AliMUONGeometryTransformer(true);
-  transform->ReadGeometryData("volpath.dat", "transform.dat");
+  AliMUONGeometryTransformer *transform = new AliMUONGeometryTransformer();
+  transform->LoadGeometryData();
   alig->SetGeometryTransformer(transform);
 
-  char cFileName[100];  
-  AliMUONDataInterface amdi;
+  // Do alignment with magnetic field off
+  alig->SetBFieldOn(kFALSE);
+  
+  // Set tracking station to use
+  Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
 
+  // Fix parameters or add constraints here
+  for (Int_t iSt=0; iSt<5; iSt++)
+    if (!bStOnOff[iSt]) alig->FixStation(iSt+1);
+
+  // Left and right sides of the detector are independent, one can choose to align 
+  // only one side
+  Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE};
+  alig->FixHalfSpectrometer(bStOnOff,bSpecLROnOff);
+
+  // Set predifined global constrains: X, Y, P, XvsZ, YvsZ, PvsZ, XvsY, YvsY, PvsY
+  Bool_t bVarXYT[9] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
+  Bool_t bDetTLBR[4] = {kFALSE,kTRUE,kFALSE,kTRUE};
+  alig->AddConstraints(bStOnOff,bVarXYT,bDetTLBR,bSpecLROnOff);
+
+
+  char cFileName[100];  
+  AliMUONDataInterface *amdi =0x0;
+  AliMUONVTrackStore *trackStore =0x0;
+  AliMUONTrack* track =0x0;  
+    
   Int_t lMaxFile = 1000;
   Int_t iFile = 0;
+  Int_t iEvent = 0;
   bool bKeepLoop = kTRUE;
+  Int_t iTrackTot=0;
+  Int_t iTrackOk=0;
+
   while(bKeepLoop && iFile<lMaxFile){
     iFile++;
     if (bLoop) {
@@ -93,29 +166,39 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
     }
     if (!strstr(cFileName,"galice.root")) continue;      
     cout << "Using file: " << cFileName << endl;
-    amdi.SetFile(cFileName);
-    Int_t nevents = amdi.NumberOfEvents();
+
+    if (!amdi) amdi = new AliMUONDataInterface(cFileName);
+    else amdi->Open(cFileName);
+
+    Int_t nevents = amdi->NumberOfEvents();
+    cout << "... with " << nevents << endl;
     for(Int_t event = 0; event < nevents; event++) {
-      amdi.GetEvent(event);
-      Int_t ntracks = amdi.NumberOfRecTracks();
-      cout << " there are " << ntracks << " tracks in event " << event << endl;
-      Int_t iTrack=0;
-      Int_t iTrackOk=0;
-      AliMUONTrack* track = amdi.RecTrack(iTrack);
-      while(track) {   
-       Double_t invBenMom = track->GetTrackParamAtVertex()->GetInverseBendingMomentum();
+      if (iEvent >= nEvents){
+       bKeepLoop = kFALSE;
+       break;
+      }
+      iEvent++;
+
+      trackStore = amdi->TrackStore(event);
+
+      Int_t ntracks = trackStore->GetSize();
+      if (!event%100)
+       cout << " there are " << ntracks << " tracks in event " << event << endl;
+      TIter next(trackStore->CreateIterator());
+      while ( ( track = static_cast<AliMUONTrack*>(next()) ) ) {
+       AliMUONTrackParam trackParam(*((AliMUONTrackParam*)(track->GetTrackParamAtHit()->First())));
+       AliMUONTrackExtrap::ExtrapToVertex(&trackParam,0.,0.,0.);
+       Double_t invBenMom = trackParam.GetInverseBendingMomentum();
        fInvBenMom->Fill(invBenMom);
        fBenMom->Fill(1./invBenMom);
        if (TMath::Abs(invBenMom)<=1.04) {
-         cout << "Track " << iTrack << endl;
          alig->ProcessTrack(track);
-         cout << "Calling LocalFit" << endl;
          alig->LocalFit(iTrackOk++,trackParams,0);
        }
-       iTrack++;
-       track = amdi.RecTrack(iTrack);
+       iTrackTot++;
       }
     }
+    cout << "Processed " << iTrackTot << " Tracks so far." << endl;
   }
   alig->GlobalFit(parameters,errors,pulls);
 
@@ -134,7 +217,7 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
   Double_t MSDEytErr[156] = {0};
   Int_t lNDetElem = 4*2+4*2+18*2+26*2+26*2;
   Int_t lNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
-  Int_t lSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
+  // Int_t lSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
   Int_t idOffset = 0; // 400
   Int_t lSDetElemCh = 0;
   for(Int_t iDE=0; iDE<lNDetElem; iDE++){
@@ -186,7 +269,7 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
   newTransform->WriteTransformations("transform2ReAlign.dat");
   
   // Generate realigned data in local cdb
-  TClonesArray* array = newTransform->GetMisAlignmentData();
+  const TClonesArray* array = newTransform->GetMisAlignmentData();
    
   // CDB manager
   AliCDBManager* cdbManager = AliCDBManager::Instance();
@@ -196,6 +279,6 @@ void MUONAlignment(Int_t nEvents = 100000, TString fileName = "galice.root", TSt
   cdbData->SetResponsible("Dimuon Offline project");
   cdbData->SetComment("MUON alignment objects with residual misalignment");
   AliCDBId id("MUON/Align/Data", 0, 0); 
-  cdbManager->Put(array, id, cdbData);
+  cdbManager->Put(const_cast<TClonesArray*>(array), id, cdbData);
 
 }