OCDB details added in runSimulation.C and runReconstruction.C
authorauras <auras@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 3 Apr 2012 00:42:18 +0000 (00:42 +0000)
committerauras <auras@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 3 Apr 2012 00:42:18 +0000 (00:42 +0000)
15 files changed:
MFT/AliMFT.cxx
MFT/AliMFTClusterQA.cxx
MFT/AliMFTClusterQA.h
MFT/AliMFTConstants.h
MFT/AliMFTDigit.h
MFT/AliMFTPlane.cxx
MFT/AliMFTPlane.h
MFT/AliMFTSegmentation.cxx
MFT/AliMFTSegmentation.h
MFT/AliMuonForwardTrackAnalysis.C
MFT/AliMuonForwardTrackPair.cxx
MFT/AliMuonForwardTrackPair.h
MFT/SetMFTGeometry.C
MFT/runReconstruction.C
MFT/runSimulation.C

index f4b9a4c..1a6ae16 100644 (file)
@@ -499,13 +499,62 @@ void AliMFT::Hits2SDigitsLocal(TClonesArray *hits, const TObjArray *pSDig, Int_t
     }
 
     for (Int_t iSideDigit=0; iSideDigit<fSideDigits->GetEntries(); iSideDigit++) {
-      AliMFTDigit *newDigit = (AliMFTDigit*) fSideDigits->At(iSideDigit);
-      new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(*newDigit);
+      AliMFTDigit *newSDigit = (AliMFTDigit*) fSideDigits->At(iSideDigit);
+      new ((*pSDigList[sDigit.GetPlane()])[pSDigList[sDigit.GetPlane()]->GetEntries()]) AliMFTDigit(*newSDigit);
     }
 
-    fSideDigits->Clear();    
+    fSideDigits->Delete();  
 
   }
+
+  // ------------ In case we should simulate a rectangular pattern of pixel...
+  
+  for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) { 
+    if (fSegmentation->GetPlane(iPlane)->HasPixelRectangularPatternAlongY()) {
+      Int_t nSDigits = pSDigList[iPlane]->GetEntries();
+      for (Int_t iSDigit=0; iSDigit<nSDigits; iSDigit++) {
+       AliMFTDigit *mySDig = (AliMFTDigit*) (pSDigList[iPlane]->At(iSDigit));
+       if (mySDig->GetPixelX()%2 == mySDig->GetPixelY()%2) {   // both pair or both odd
+         Int_t xPixelNew = mySDig->GetPixelX();
+         Int_t yPixelNew = mySDig->GetPixelY()+1;
+         if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
+           AliMFTDigit newSDigit;
+           newSDigit.SetEloss(0.);
+           newSDigit.SetDetElemID(mySDig->GetDetElemID());
+           newSDigit.SetPlane(iPlane);
+           newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
+           newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()), 
+                                 fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
+                                 fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));  
+           newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew), 
+                                  fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
+                                  fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0)); 
+           new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
+         }
+       }
+       else {   // pair-odd
+         Int_t xPixelNew = mySDig->GetPixelX();
+         Int_t yPixelNew = mySDig->GetPixelY()-1;
+         if (fSegmentation->DoesPixelExist(mySDig->GetDetElemID(), xPixelNew, yPixelNew)) {
+           AliMFTDigit newSDigit;
+           newSDigit.SetEloss(0.);
+           newSDigit.SetDetElemID(mySDig->GetDetElemID());
+           newSDigit.SetPlane(iPlane);
+           newSDigit.SetPixID(xPixelNew, yPixelNew, 0);
+           newSDigit.SetPixWidth(fSegmentation->GetPixelSizeX(newSDigit.GetDetElemID()), 
+                                 fSegmentation->GetPixelSizeY(newSDigit.GetDetElemID()),
+                                 fSegmentation->GetPixelSizeZ(newSDigit.GetDetElemID()));  
+           newSDigit.SetPixCenter(fSegmentation->GetPixelCenterX(newSDigit.GetDetElemID(), xPixelNew), 
+                                  fSegmentation->GetPixelCenterY(newSDigit.GetDetElemID(), yPixelNew),
+                                  fSegmentation->GetPixelCenterZ(newSDigit.GetDetElemID(), 0)); 
+           new ((*pSDigList[iPlane])[pSDigList[iPlane]->GetEntries()]) AliMFTDigit(newSDigit);
+         }
+       }
+      }
+    }
+  }
+
+  //------------------------------------------------------------------------
   
   AliDebug(1,"Stop Hits2SDigitsLocal");
 
index 26bb99a..6168b97 100644 (file)
@@ -42,11 +42,12 @@ AliMFTClusterQA::AliMFTClusterQA():
   // default constructor
 
   for (Int_t iPlane=0; iPlane<fNMaxPlanes; iPlane++) {
-    fHistNClustersPerEvent[iPlane] = 0;
-    fHistNPixelsPerCluster[iPlane] = 0;
-    fHistClusterSizeX[iPlane] = 0; 
-    fHistClusterSizeY[iPlane] = 0;
-    fClusterScatterPlotXY[iPlane] = 0;
+    fHistNClustersPerEvent[iPlane]     = 0;
+    fHistNPixelsPerCluster[iPlane]     = 0;
+    fHistClusterSizeX[iPlane]          = 0; 
+    fHistClusterSizeY[iPlane]          = 0;
+    fHistClusterRadialPosition[iPlane] = 0;
+    fClusterScatterPlotXY[iPlane]      = 0;
   }
 
 }
@@ -96,7 +97,8 @@ Bool_t AliMFTClusterQA::LoadNextEvent() {
       fHistNPixelsPerCluster[iPlane] -> Fill(cluster->GetSize());
       fHistClusterSizeX[iPlane]      -> Fill(cluster->GetErrX()*1.e4);   // converted in microns
       fHistClusterSizeY[iPlane]      -> Fill(cluster->GetErrY()*1.e4);   // converted in microns
-      fClusterScatterPlotXY[iPlane]  -> Fill(cluster->GetX(), cluster->GetY());
+      fHistClusterRadialPosition[iPlane] -> Fill(TMath::Sqrt(cluster->GetX()*cluster->GetX()+cluster->GetY()*cluster->GetY()));
+      fClusterScatterPlotXY[iPlane]      -> Fill(cluster->GetX(), cluster->GetY());
     }
   }
 
@@ -140,12 +142,22 @@ void AliMFTClusterQA::BookHistos() {
 
     //------------------------------------------------------------
 
-    Int_t rMax = Int_t(10.*(fMFT->GetSegmentation()->GetPlane(iPlane)->GetRMaxSupport()));
+    Int_t rMax = Int_t(10.*(fMFT->GetSegmentation()->GetPlane(iPlane)->GetRMaxSupport()));  // expressed in mm
+
+    fHistClusterRadialPosition[iPlane] = new TH1D(Form("fHistClusterRadialPosition_Pl%02d",iPlane),
+                                                 Form("Cluster radial position (Plane%02d)",iPlane),
+                                                 rMax, 0, Double_t(rMax)/10.);
+
     fClusterScatterPlotXY[iPlane] = new TH2D(Form("fClusterScatterPlotXY_Pl%02d",iPlane),
                                             Form("Cluster scatter plot (Plane%02d)",iPlane),
                                             2*rMax+1, (-rMax-0.5)/10., (rMax+0.5)/10., 2*rMax+1, (-rMax-0.5)/10., (rMax+0.5)/10.);
     
-    fClusterScatterPlotXY[iPlane] -> Sumw2();
+    fHistClusterRadialPosition[iPlane] -> SetXTitle("R  [cm]");
+    fClusterScatterPlotXY[iPlane]      -> SetXTitle("X  [cm]");
+    fClusterScatterPlotXY[iPlane]      -> SetYTitle("Y  [cm]");
+
+    fHistClusterRadialPosition[iPlane] -> Sumw2();
+    fClusterScatterPlotXY[iPlane]      -> Sumw2();
     
   }
   
@@ -157,14 +169,31 @@ void AliMFTClusterQA::Terminate() {
 
   AliInfo("Writing QA histos...");
 
+  // ---- equalize radial clusters distribution ----------------------
+
+  for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
+    for (Int_t iBin=1; iBin<=fHistClusterRadialPosition[iPlane]->GetNbinsX(); iBin++) {
+      Double_t rMin = fHistClusterRadialPosition[iPlane]->GetBinLowEdge(iBin);        // in cm
+      Double_t rMax = fHistClusterRadialPosition[iPlane]->GetBinWidth(iBin) + rMin;   // in cm
+      Double_t area = 100.*TMath::Pi()*(rMax*rMax - rMin*rMin);                       // in mm^2
+      Double_t density = fHistClusterRadialPosition[iPlane]->GetBinContent(iBin)/area;
+      fHistClusterRadialPosition[iPlane]->SetBinContent(iBin, density);
+      fHistClusterRadialPosition[iPlane]->SetBinError(iBin, fHistClusterRadialPosition[iPlane]->GetBinError(iBin)/area);
+    }
+    fHistClusterRadialPosition[iPlane] -> SetBinContent(1, fEv);      // "scaler" bin
+  }
+
+  // -----------------------------------------------------------------
+
   fFileOut->cd();
 
   for (Int_t iPlane=0; iPlane<fNPlanes; iPlane++) {
-    fHistNClustersPerEvent[iPlane] -> Write();
-    fHistNPixelsPerCluster[iPlane] -> Write();
-    fHistClusterSizeX[iPlane]      -> Write();
-    fHistClusterSizeY[iPlane]      -> Write();
-    fClusterScatterPlotXY[iPlane]  -> Write();
+    fHistNClustersPerEvent[iPlane]     -> Write();
+    fHistNPixelsPerCluster[iPlane]     -> Write();
+    fHistClusterSizeX[iPlane]          -> Write();
+    fHistClusterSizeY[iPlane]          -> Write();
+    fHistClusterRadialPosition[iPlane] -> Write();
+    fClusterScatterPlotXY[iPlane]      -> Write();
   }
 
   fFileOut -> Close();
index b1cb3a5..cc0f9fa 100644 (file)
@@ -45,7 +45,7 @@ protected:
   static const Int_t fNMaxPlanes = AliMFTConstants::fNMaxPlanes;
 
   TH1D *fHistNClustersPerEvent[fNMaxPlanes], *fHistNPixelsPerCluster[fNMaxPlanes];
-  TH1D *fHistClusterSizeX[fNMaxPlanes], *fHistClusterSizeY[fNMaxPlanes];
+  TH1D *fHistClusterSizeX[fNMaxPlanes], *fHistClusterSizeY[fNMaxPlanes], *fHistClusterRadialPosition[fNMaxPlanes];
   TH2D *fClusterScatterPlotXY[fNMaxPlanes];
 
   AliLoader *fMFTLoader;
@@ -62,7 +62,7 @@ protected:
 
 };
 
-//======================================================================================================
+//====================================================================================================================================================
  
 #endif
 
index 524c381..91dea82 100644 (file)
@@ -20,7 +20,7 @@ public:
   
   static const Int_t fNMaxPlanes = 20; 
 
-  static const Int_t fNMaxDigitsPerCluster = 12;  ///< max number of digits per cluster
+  static const Int_t fNMaxDigitsPerCluster = 50;  ///< max number of digits per cluster
   static const Double_t fCutForAvailableDigits;   ///<
   static const Double_t fCutForAttachingDigits;   ///<
 
index 715e45a..26db393 100644 (file)
@@ -60,7 +60,7 @@ public:
   Double_t  GetPixelWidthX()   const { return fPixelWidthX; }
   Double_t  GetPixelWidthY()   const { return fPixelWidthY; }
   Double_t  GetPixelWidthZ()   const { return fPixelWidthZ; }
-    
+
 protected:
     
   static const Double_t fElossPerElectron;
index 2212f81..4197447 100644 (file)
@@ -63,7 +63,8 @@ AliMFTPlane::AliMFTPlane():
   fEquivalentSiliconBeforeBack(0),
   fActiveElements(0),
   fReadoutElements(0),
-  fSupportElements(0)
+  fSupportElements(0),
+  fHasPixelRectangularPatternAlongY(kFALSE)
 {
 
   // default constructor
@@ -91,7 +92,8 @@ AliMFTPlane::AliMFTPlane(const Char_t *name, const Char_t *title):
   fEquivalentSiliconBeforeBack(0),
   fActiveElements(new TClonesArray("THnSparseC")),
   fReadoutElements(new TClonesArray("THnSparseC")),
-  fSupportElements(new TClonesArray("THnSparseC"))
+  fSupportElements(new TClonesArray("THnSparseC")),
+  fHasPixelRectangularPatternAlongY(kFALSE)
 {
 
   // constructor
@@ -119,7 +121,8 @@ AliMFTPlane::AliMFTPlane(const AliMFTPlane& plane):
   fEquivalentSiliconBeforeBack(plane.fEquivalentSiliconBeforeBack),
   fActiveElements(new TClonesArray("THnSparseC")),
   fReadoutElements(new TClonesArray("THnSparseC")),
-  fSupportElements(new TClonesArray("THnSparseC"))
+  fSupportElements(new TClonesArray("THnSparseC")),
+  fHasPixelRectangularPatternAlongY(plane.fHasPixelRectangularPatternAlongY)
 {
 
   // copy constructor
@@ -142,24 +145,25 @@ AliMFTPlane& AliMFTPlane::operator=(const AliMFTPlane& plane) {
     // base class assignement
     TNamed::operator=(plane);
     
-    fPlaneNumber                  = plane.fPlaneNumber;
-    fZCenter                      = plane.fZCenter; 
-    fRMinSupport                  = plane.fRMinSupport; 
-    fRMax                         = plane.fRMax;
-    fRMaxSupport                  = plane.fRMaxSupport;
-    fPixelSizeX                   = plane.fPixelSizeX;
-    fPixelSizeY                   = plane.fPixelSizeY; 
-    fThicknessActive              = plane.fThicknessActive; 
-    fThicknessSupport             = plane.fThicknessSupport; 
-    fThicknessReadout             = plane.fThicknessReadout;
-    fZCenterActiveFront           = plane.fZCenterActiveFront;
-    fZCenterActiveBack            = plane.fZCenterActiveBack;
-    fEquivalentSilicon            = plane.fEquivalentSilicon;
-    fEquivalentSiliconBeforeFront = plane.fEquivalentSiliconBeforeFront;
-    fEquivalentSiliconBeforeBack  = plane.fEquivalentSiliconBeforeBack;
-    *fActiveElements              = *plane.fActiveElements;
-    *fReadoutElements             = *plane.fReadoutElements;
-    *fSupportElements             = *plane.fSupportElements;
+    fPlaneNumber                      = plane.fPlaneNumber;
+    fZCenter                          = plane.fZCenter; 
+    fRMinSupport                      = plane.fRMinSupport; 
+    fRMax                             = plane.fRMax;
+    fRMaxSupport                      = plane.fRMaxSupport;
+    fPixelSizeX                       = plane.fPixelSizeX;
+    fPixelSizeY                       = plane.fPixelSizeY; 
+    fThicknessActive                  = plane.fThicknessActive; 
+    fThicknessSupport                 = plane.fThicknessSupport; 
+    fThicknessReadout                 = plane.fThicknessReadout;
+    fZCenterActiveFront               = plane.fZCenterActiveFront;
+    fZCenterActiveBack                = plane.fZCenterActiveBack;
+    fEquivalentSilicon                = plane.fEquivalentSilicon;
+    fEquivalentSiliconBeforeFront     = plane.fEquivalentSiliconBeforeFront;
+    fEquivalentSiliconBeforeBack      = plane.fEquivalentSiliconBeforeBack;
+    *fActiveElements                  = *plane.fActiveElements;
+    *fReadoutElements                 = *plane.fReadoutElements;
+    *fSupportElements                 = *plane.fSupportElements;
+    fHasPixelRectangularPatternAlongY = plane.fHasPixelRectangularPatternAlongY;
   }
 
   return *this;
@@ -176,7 +180,8 @@ Bool_t AliMFTPlane::Init(Int_t    planeNumber,
                         Double_t pixelSizeY, 
                         Double_t thicknessActive, 
                         Double_t thicknessSupport, 
-                        Double_t thicknessReadout) {
+                        Double_t thicknessReadout,
+                        Bool_t   hasPixelRectangularPatternAlongY) {
 
   AliDebug(1, Form("Initializing Plane Structure for Plane %s", GetName()));
 
@@ -190,6 +195,8 @@ Bool_t AliMFTPlane::Init(Int_t    planeNumber,
   fThicknessSupport = thicknessSupport;
   fThicknessReadout = thicknessReadout;
 
+  fHasPixelRectangularPatternAlongY = hasPixelRectangularPatternAlongY;
+
   fZCenterActiveFront = fZCenter - 0.5*fThicknessSupport - 0.5*fThicknessActive;
   fZCenterActiveBack  = fZCenter + 0.5*fThicknessSupport + 0.5*fThicknessActive;
 
index 274f945..fc2384d 100644 (file)
@@ -44,7 +44,8 @@ public:
              Double_t pixelSizeY, 
              Double_t thicknessActive, 
              Double_t thicknessSupport, 
-             Double_t thicknessReadout);
+             Double_t thicknessReadout,
+             Bool_t   hasPixelRectangularPatternAlongY);
   
   Bool_t CreateStructure();
 
@@ -80,6 +81,7 @@ public:
   Double_t GetEquivalentSiliconBeforeBack()  const { return fEquivalentSiliconBeforeBack; }
 
   Int_t GetNumberOfChips(Option_t *opt);
+  Bool_t HasPixelRectangularPatternAlongY() { return fHasPixelRectangularPatternAlongY; }
   
 private:
 
@@ -101,6 +103,8 @@ private:
 
   TClonesArray *fActiveElements, *fReadoutElements, *fSupportElements;
 
+  Bool_t fHasPixelRectangularPatternAlongY;
+
   ClassDef(AliMFTPlane, 1)
 
 };
index fc89866..675fbc5 100644 (file)
@@ -51,7 +51,7 @@ AliMFTSegmentation::AliMFTSegmentation(const Char_t *nameGeomFile):
   // constructor
 
   Float_t zCenter, rMin, rMax, pixelSizeX, pixelSizeY, thicknessActive, thicknessSupport, thicknessReadout;
-  Float_t equivalentSilicon, equivalentSiliconBeforeFront, equivalentSiliconBeforeBack;
+  Float_t equivalentSilicon, equivalentSiliconBeforeFront, equivalentSiliconBeforeBack, hasPixelRectangularPatternAlongY;
 
   TFile *geomFile = new TFile(nameGeomFile);
   TNtuple *geomNtuple = (TNtuple*) geomFile->Get("AliMFTGeometry");
@@ -67,6 +67,10 @@ AliMFTSegmentation::AliMFTSegmentation(const Char_t *nameGeomFile):
   geomNtuple -> SetBranchAddress("equivalentSilicon",            &equivalentSilicon);
   geomNtuple -> SetBranchAddress("equivalentSiliconBeforeFront", &equivalentSiliconBeforeFront);
   geomNtuple -> SetBranchAddress("equivalentSiliconBeforeBack",  &equivalentSiliconBeforeBack);
+  if (geomNtuple -> GetBranch("hasPixelRectangularPatternAlongY")) {
+    geomNtuple -> SetBranchAddress("hasPixelRectangularPatternAlongY", &hasPixelRectangularPatternAlongY);
+  }
+  else hasPixelRectangularPatternAlongY = 0.;
   
   Int_t nPlanes = geomNtuple->GetEntries();
 
@@ -80,7 +84,18 @@ AliMFTSegmentation::AliMFTSegmentation(const Char_t *nameGeomFile):
     zCenter = TMath::Abs(zCenter);
 
     AliMFTPlane *plane = new AliMFTPlane(Form("MFTPlane_%02d", iPlane), Form("MFTPlane_%02d", iPlane));
-    plane -> Init(iPlane, zCenter, rMin, rMax, pixelSizeX, pixelSizeY, thicknessActive, thicknessSupport, thicknessReadout);
+
+    plane -> Init(iPlane, 
+                 zCenter, 
+                 rMin, 
+                 rMax, 
+                 pixelSizeX, 
+                 pixelSizeY, 
+                 thicknessActive, 
+                 thicknessSupport, 
+                 thicknessReadout, 
+                 (hasPixelRectangularPatternAlongY>0.5));
+
     plane -> SetEquivalentSilicon(equivalentSilicon);
     plane -> SetEquivalentSiliconBeforeFront(equivalentSiliconBeforeFront);
     plane -> SetEquivalentSiliconBeforeBack(equivalentSiliconBeforeBack);
@@ -115,6 +130,8 @@ THnSparseC* AliMFTSegmentation::GetDetElem(Int_t detElemID) const {
 
 Bool_t AliMFTSegmentation::Hit2PixelID(Double_t xHit, Double_t yHit, Int_t detElemID, Int_t &xPixel, Int_t &yPixel) {
 
+  // xPixel and yPixel start from 0
+
   THnSparseC *detElem = GetDetElem(detElemID);
 
   if ( xHit<detElem->GetAxis(0)->GetXmin() ||
@@ -131,3 +148,14 @@ Bool_t AliMFTSegmentation::Hit2PixelID(Double_t xHit, Double_t yHit, Int_t detEl
 
 //====================================================================================================================================================
 
+Bool_t AliMFTSegmentation::DoesPixelExist(Int_t detElemID, Int_t xPixel, Int_t yPixel) {
+
+  THnSparseC *detElem = GetDetElem(detElemID);
+
+  if (xPixel>=0 && xPixel<detElem->GetAxis(0)->GetNbins() && yPixel>=0 && yPixel<detElem->GetAxis(1)->GetNbins()) return kTRUE;
+  else return kFALSE;
+
+}
+
+//====================================================================================================================================================
+   
index 9600a3d..5072bae 100644 (file)
@@ -48,6 +48,8 @@ public:
   Int_t GetNPlanes() const { return fMFTPlanes->GetEntries(); }
 
   AliMFTPlane* GetPlane(Int_t iPlane) const { if (iPlane>=0 && iPlane<fMFTPlanes->GetEntries()) return (AliMFTPlane*) fMFTPlanes->At(iPlane); else return NULL; }
+
+  Bool_t DoesPixelExist(Int_t detElemID, Int_t xPixel, Int_t yPixel);
  
 protected:
 
index b812175..74be606 100644 (file)
@@ -7,13 +7,14 @@ void AliMuonForwardTrackAnalysis(const Char_t *readDir= ".",
                                 Int_t nMassBin = 100, 
                                 Double_t massMin = 0.,
                                 Double_t massMax = 10.,
+                                Bool_t useCutOnOffsetChi2 = kFALSE,
+                                Int_t maxNWrongClusters = 999,
                                 const Char_t *outDir = ".",
                                 Bool_t singleMuonAnalysis = kTRUE,
                                 Bool_t muonPairAnalysis = kTRUE,
                                 Int_t firstEvent = -1,
                                 Int_t lastEvent = -1, 
                                 Int_t myRandom = 0,
-                                Int_t maxNWrongClusters = 999,
                                 Double_t ptMinSingleMuons = 0.0) {
   
   gROOT -> LoadMacro("./AliMuonForwardTrackAnalysis.cxx+");
@@ -29,16 +30,15 @@ void AliMuonForwardTrackAnalysis(const Char_t *readDir= ".",
   myAnalysis->SetInputDir(readDir);
   myAnalysis->SetOutputDir(outDir);
   myAnalysis->SetMassRange(nMassBin, massMin, massMax);
-  myAnalysis->SetPtDimuRange(10, 0., 5.);
   myAnalysis->SetSingleMuonAnalysis(singleMuonAnalysis);
   myAnalysis->SetMuonPairAnalysis(muonPairAnalysis);
   myAnalysis->SetOption(option);
-  myAnalysis->SetMatchTrigger(kTRUE);
   myAnalysis->SetMaxNWrongClustersMC(maxNWrongClusters);
   myAnalysis->SetPtMinSingleMuons(ptMinSingleMuons);
+  myAnalysis->UseCutOnOffsetChi2(useCutOnOffsetChi2);
 
-  myAnalysis->UseCutOnOffsetChi2(kFALSE);        // cut on the single muons
-
+  myAnalysis->SetPtDimuRange(10, 0., 5.);
+  myAnalysis->SetMatchTrigger(kTRUE);
   myAnalysis->UseBransonForCut(kFALSE);
   myAnalysis->UseBransonForKinematics(kFALSE);
 
index 474c52c..b56cf66 100644 (file)
@@ -39,7 +39,9 @@ ClassImp(AliMuonForwardTrackPair)
 AliMuonForwardTrackPair::AliMuonForwardTrackPair():
   TObject(),
   fMuonForwardTracks(0),
-  fKinemMC(0,0,0,0)
+  fKinemMC(0,0,0,0),
+  fKinem(0,0,0,0),
+  fIsKinemSet(kFALSE)
 {
 
   // default constructor
@@ -53,7 +55,9 @@ AliMuonForwardTrackPair::AliMuonForwardTrackPair():
 AliMuonForwardTrackPair::AliMuonForwardTrackPair(AliMuonForwardTrack *track0, AliMuonForwardTrack *track1):
   TObject(),
   fMuonForwardTracks(0),
-  fKinemMC(0,0,0,0)
+  fKinemMC(0,0,0,0),
+  fKinem(0,0,0,0),
+  fIsKinemSet(kFALSE)
 {
 
   fMuonForwardTracks = new TClonesArray("AliMuonForwardTrack", 2);
@@ -70,7 +74,9 @@ AliMuonForwardTrackPair::AliMuonForwardTrackPair(AliMuonForwardTrack *track0, Al
 AliMuonForwardTrackPair::AliMuonForwardTrackPair(const AliMuonForwardTrackPair& trackPair): 
   TObject(trackPair),
   fMuonForwardTracks(trackPair.fMuonForwardTracks),
-  fKinemMC(trackPair.fKinemMC)
+  fKinemMC(trackPair.fKinemMC),
+  fKinem(trackPair.fKinem),
+  fIsKinemSet(trackPair.fIsKinemSet)
 {
 
   // copy constructor
@@ -93,6 +99,9 @@ AliMuonForwardTrackPair& AliMuonForwardTrackPair::operator=(const AliMuonForward
   Clear();
   
   fMuonForwardTracks = trackPair.fMuonForwardTracks;
+  fKinemMC = trackPair.fKinemMC;
+  fKinem = trackPair.fKinem;
+  fIsKinemSet = trackPair.fIsKinemSet;
 
   return *this;
 
@@ -125,55 +134,6 @@ Double_t AliMuonForwardTrackPair::GetWeightedOffset(Double_t x, Double_t y, Doub
 
 //====================================================================================================================================================
 
-Double_t AliMuonForwardTrackPair::GetMass(Double_t z, Int_t nClusters) {
-
-  Int_t idCluster[2] = {0};
-  if (nClusters>0) {
-    idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMFTClusters() - nClusters;
-    idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMFTClusters() - nClusters;
-  }
-  if (idCluster[0]<0) idCluster[0] = 0;
-  if (idCluster[1]<0) idCluster[1] = 0;
-
-  Double_t momentum[2] = {0};
-  
-  AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMFTCluster(idCluster[0]);
-  AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMFTCluster(idCluster[1]);
-
-  AliDebug(2, Form("MFT before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
-                  param0->Px(), param0->Py(), param0->Pz(), 
-                  param1->Px(), param1->Py(), param1->Pz()));
-
-  if (TMath::Abs(z)<1e6) {
-    AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
-    AliMUONTrackExtrap::ExtrapToZCov(param0, z);
-    AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
-    AliMUONTrackExtrap::ExtrapToZCov(param1, z);
-  }
-
-  AliDebug(2, Form("MFT after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
-                  param0->Px(), param0->Py(), param0->Pz(), 
-                  param1->Px(), param1->Py(), param1->Pz()));
-
-  momentum[0] = (param0->P());
-  momentum[1] = (param1->P());
-
-  Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
-
-  TLorentzVector dimu;
-
-  dimu.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
-
-  dimu.SetPx(param0->Px() + param1->Px());
-  dimu.SetPy(param0->Py() + param1->Py());
-  dimu.SetPz(param0->Pz() + param1->Pz());
-
-  return dimu.M();
-
-}
-
-//====================================================================================================================================================
-
 Double_t AliMuonForwardTrackPair::GetMassWithoutMFT(Double_t x, Double_t y, Double_t z, Int_t nClusters) {
 
   Int_t idCluster[2] = {0};
@@ -223,6 +183,9 @@ Double_t AliMuonForwardTrackPair::GetMassWithoutMFT(Double_t x, Double_t y, Doub
 
 void AliMuonForwardTrackPair::SetKinemMC() {
 
+  if ( !(((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()) || 
+       !(((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetMCTrackRef()) ) return;
+
   AliDebug(2, Form("MC: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
                   ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Px(),
                   ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetMCTrackRef()->Py(),
@@ -247,6 +210,57 @@ void AliMuonForwardTrackPair::SetKinemMC() {
 
 //====================================================================================================================================================
 
+void AliMuonForwardTrackPair::SetKinem(Double_t z, Int_t nClusters) {
+
+//   if (!fMuonForwardTracks) return kFALSE;
+//   if (!fMuonForwardTracks->At(0) || !fMuonForwardTracks->At(1)) return kFALSE;
+
+  Int_t idCluster[2] = {0};
+  if (nClusters>0) {
+    idCluster[0] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetNMFTClusters() - nClusters;
+    idCluster[1] = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetNMFTClusters() - nClusters;
+  }
+  if (idCluster[0]<0) idCluster[0] = 0;
+  if (idCluster[1]<0) idCluster[1] = 0;
+
+  Double_t momentum[2] = {0};
+  
+  AliMUONTrackParam *param0 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(0))->GetTrackParamAtMFTCluster(idCluster[0]);
+  AliMUONTrackParam *param1 = ((AliMuonForwardTrack*) fMuonForwardTracks->At(1))->GetTrackParamAtMFTCluster(idCluster[1]);
+
+  AliDebug(2, Form("MFT before extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
+                  param0->Px(), param0->Py(), param0->Pz(), 
+                  param1->Px(), param1->Py(), param1->Pz()));
+
+  if (TMath::Abs(z)<1e6) {
+    AliDebug(2, Form("Extrapolating 1st muon from z = %f to z = %f", param0->GetZ(), z));
+    AliMUONTrackExtrap::ExtrapToZCov(param0, z);
+    AliDebug(2, Form("Extrapolating 2nd muon from z = %f to z = %f", param1->GetZ(), z));
+    AliMUONTrackExtrap::ExtrapToZCov(param1, z);
+  }
+
+  AliDebug(2, Form("MFT after extrap: 1st muon = (%f, %f, %f) 2nd muon = (%f, %f, %f)", 
+                  param0->Px(), param0->Py(), param0->Pz(), 
+                  param1->Px(), param1->Py(), param1->Pz()));
+
+  momentum[0] = (param0->P());
+  momentum[1] = (param1->P());
+
+  Double_t mMu = TDatabasePDG::Instance()->GetParticle("mu-")->Mass();
+
+  fKinem.SetE(TMath::Sqrt(mMu*mMu + momentum[0]*momentum[0]) + TMath::Sqrt(mMu*mMu + momentum[1]*momentum[1]));
+  fKinem.SetPx(param0->Px() + param1->Px());
+  fKinem.SetPy(param0->Py() + param1->Py());
+  fKinem.SetPz(param0->Pz() + param1->Pz());
+
+  fIsKinemSet = kTRUE;
+
+  //  return fKinem.M();
+
+}
+
+//====================================================================================================================================================
+
 Bool_t AliMuonForwardTrackPair::IsResonance() {
 
   Bool_t result = kFALSE;
index 9e39a5f..7a69e63 100644 (file)
@@ -42,20 +42,25 @@ public:
   }
 
   void SetKinemMC();
+  void SetKinem(Double_t z, Int_t nClusters=-1);
+  Bool_t IsKinemSet() { return fIsKinemSet; }
 
   Double_t GetWeightedOffset(Double_t x, Double_t y, Double_t z);
-  Double_t GetMass(Double_t z, Int_t nClusters=-1);
   Double_t GetMassWithoutMFT(Double_t x, Double_t y, Double_t z, Int_t nClusters=-1);
   Double_t GetMassMC()     { return fKinemMC.M(); }
   Double_t GetRapidityMC() { return fKinemMC.Rapidity(); }
   Double_t GetPtMC()       { return fKinemMC.Pt(); }
+  Double_t GetMass()     { return fKinem.M(); }
+  Double_t GetRapidity() { return fKinem.Rapidity(); }
+  Double_t GetPt()       { return fKinem.Pt(); }
 
   Bool_t IsResonance();
 
 protected:
 
   TClonesArray *fMuonForwardTracks;
-  TLorentzVector fKinemMC;
+  TLorentzVector fKinemMC, fKinem;
+  Bool_t fIsKinemSet;
 
   ClassDef(AliMuonForwardTrackPair,1)
     
index a86e6f9..15ff883 100644 (file)
@@ -26,7 +26,9 @@ void SetMFTGeometry() {
   const Float_t equivalentSiliconBeforeFront[nPlanes] = {   0.e-4,   0.e-4,   0.e-4,   0.e-4,   0.e-4};    // expressed in cm
   const Float_t equivalentSiliconBeforeBack[nPlanes]  = { 250.e-4, 250.e-4, 250.e-4, 250.e-4, 250.e-4};    // expressed in cm
 
-  TNtuple *geomMFT = new TNtuple("AliMFTGeometry", "ALICE MFT Geometry", "zCenter:rMin:rMax:pixelSizeX:pixelSizeY:thicknessActive:thicknessSupport:thicknessReadout:equivalentSilicon:equivalentSiliconBeforeFront:equivalentSiliconBeforeBack");
+  const Float_t hasPixelRectangularPatternAlongY[nPlanes] = {0., 0., 0., 0., 0.};
+
+  TNtuple *geomMFT = new TNtuple("AliMFTGeometry", "ALICE MFT Geometry", "zCenter:rMin:rMax:pixelSizeX:pixelSizeY:thicknessActive:thicknessSupport:thicknessReadout:equivalentSilicon:equivalentSiliconBeforeFront:equivalentSiliconBeforeBack:hasPixelRectangularPatternAlongY");
 
   for (Int_t iPlane=0; iPlane<nPlanes; iPlane++) geomMFT -> Fill(zCenter[iPlane],
                                                                 rMin[iPlane],
@@ -38,7 +40,8 @@ void SetMFTGeometry() {
                                                                 thicknessReadout[iPlane],
                                                                 equivalentSilicon[iPlane],
                                                                 equivalentSiliconBeforeFront[iPlane],
-                                                                equivalentSiliconBeforeBack[iPlane]);
+                                                                equivalentSiliconBeforeBack[iPlane],
+                                                                hasPixelRectangularPatternAlongY[iPlane]);
 
   TFile *fileGeomMFT = new TFile("AliMFTGeometry.root", "recreate");
   geomMFT -> Write();
index 6c52848..f8bf6f1 100644 (file)
@@ -10,8 +10,9 @@ void runReconstruction(Int_t seed, const Char_t *recOptions) {
   // GRP from local OCDB\r
   reco->SetSpecificStorage("GRP/GRP/Data",Form("local://%s",gSystem->pwd()));\r
   \r
-  // MUON Tracker\r
-  // reco->SetSpecificStorage("MUON/Align/Data","local:///$OCDB/simulation/2008/v4-15-Release/Residual");\r
+  // MUON Tracker -> local:///$OCDB should reflect the content of alien://folder=/alice\r
+  reco->SetSpecificStorage("MUON/Align/Data",      "local:///$OCDB/simulation/2008/v4-15-Release/Residual");\r
+  reco->SetSpecificStorage("MUON/Calib/RecoParam", "local:///$OCDB/simulation/2008/v4-15-Release/Full");\r
 \r
   reco->SetRunReconstruction("MUON MFT");\r
   reco->SetRunLocalReconstruction("MUON MFT");\r
index 46eb45b..70eab84 100644 (file)
@@ -13,8 +13,9 @@ void runSimulation(Int_t seed,
   simulator->SetRunQA("ALL");\r
   simulator->SetRunHLT("");\r
 \r
-  // MUON Tracker\r
-  // simulator->SetSpecificStorage("MUON/Align/Data", "local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
+  // MUON Tracker -> local:///$OCDB should reflect the content of alien://folder=/alice\r
+  simulator->SetSpecificStorage("MUON/Align/Data", "local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
+  simulator->SetSpecificStorage("MUON/Calib/Gains","local:///$OCDB/simulation/2008/v4-15-Release/Ideal");\r
 \r
   // The rest\r
   TStopwatch timer;\r