]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOSTrackSegmentMakerv1.cxx
added cdb stuff
[u/mrichter/AliRoot.git] / PHOS / AliPHOSTrackSegmentMakerv1.cxx
index a4d319e551ecea76491c39f204116bc8d21d5610..203abda60f7c51fa5731b9b7d629c57afbcd8e3d 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 /* $Id$ */
+
+/* History of cvs commits:
+ *
+ * $Log$
+ * Revision 1.79  2006/04/25 12:41:15  hristov
+ * Moving non-persistent data to AliESDfriend (Yu.Belikov)
+ *
+ * Revision 1.78  2005/11/18 13:04:51  hristov
+ * Bug fix
+ *
+ * Revision 1.77  2005/11/17 23:34:36  hristov
+ * Corrected logics
+ *
+ * Revision 1.76  2005/11/17 22:29:12  hristov
+ * Faster version, no attempt to match tracks outside the PHOS acceptance
+ *
+ * Revision 1.75  2005/11/17 12:35:27  hristov
+ * Use references instead of objects. Avoid to create objects when they are not really needed
+ *
+ * Revision 1.74  2005/07/08 14:01:36  hristov
+ * Tracking in non-uniform nmagnetic field (Yu.Belikov)
+ *
+ * Revision 1.73  2005/05/28 14:19:05  schutz
+ * Compilation warnings fixed by T.P.
+ *
+ */
+
 //_________________________________________________________________________
 // Implementation version 1 of algorithm class to construct PHOS track segments
 // Track segment for PHOS is list of 
@@ -61,26 +88,70 @@ ClassImp( AliPHOSTrackSegmentMakerv1)
 
 
 //____________________________________________________________________________
-  AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1() : AliPHOSTrackSegmentMaker()
+AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1() :
+  AliPHOSTrackSegmentMaker(),
+  fDefaultInit(kTRUE),
+  fWrite(kFALSE),
+  fNTrackSegments(0),
+  fRcpv(0.f),
+  fRtpc(0.f),
+  fLinkUpArray(0),
+  fEmcFirst(0),
+  fEmcLast(0),
+  fCpvFirst(0),
+  fCpvLast(0),
+  fModule(0),
+  fTrackSegmentsInRun(0)
+  
 {
   // default ctor (to be used mainly by Streamer)
-
   InitParameters() ; 
-  fDefaultInit = kTRUE ; 
 }
 
 //____________________________________________________________________________
- AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(const TString alirunFileName, const TString eventFolderName)
-   :AliPHOSTrackSegmentMaker(alirunFileName, eventFolderName)
+AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(const TString & alirunFileName, const TString & eventFolderName) :
+  AliPHOSTrackSegmentMaker(alirunFileName, eventFolderName),
+  fDefaultInit(kFALSE),
+  fWrite(kFALSE),
+  fNTrackSegments(0),
+  fRcpv(0.f),
+  fRtpc(0.f),
+  fLinkUpArray(0),
+  fEmcFirst(0),
+  fEmcLast(0),
+  fCpvFirst(0),
+  fCpvLast(0),
+  fModule(0),
+  fTrackSegmentsInRun(0)
 {
   // ctor
-
   InitParameters() ; 
   Init() ;
-  fDefaultInit = kFALSE ; 
   fESD = 0;
 }
 
+
+AliPHOSTrackSegmentMakerv1::AliPHOSTrackSegmentMakerv1(const AliPHOSTrackSegmentMakerv1 & tsm) :
+  AliPHOSTrackSegmentMaker(tsm),
+  fDefaultInit(kFALSE),
+  fWrite(kFALSE),
+  fNTrackSegments(0),
+  fRcpv(0.f),
+  fRtpc(0.f),
+  fLinkUpArray(0),
+  fEmcFirst(0),
+  fEmcLast(0),
+  fCpvFirst(0),
+  fCpvLast(0),
+  fModule(0),
+  fTrackSegmentsInRun(0)
+{
+  // cpy ctor: no implementation yet
+  // requested by the Coding Convention
+  Fatal("cpy ctor", "not implemented") ;
+}
+
+
 //____________________________________________________________________________
  AliPHOSTrackSegmentMakerv1::~AliPHOSTrackSegmentMakerv1()
 { 
@@ -130,10 +201,10 @@ Float_t  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint *
   // Calculates the distance between the EMC RecPoint and the CPV RecPoint
   // Clusters are sorted in "rows" and "columns" of width 1 cm
 
-  Float_t delta = 1 ;  // Width of the rows in sorting of RecPoints (in cm)
+  //Float_t delta = 1 ;  // Width of the rows in sorting of RecPoints (in cm)
                        // if you change this value, change it as well in xxxRecPoint::Compare()
   Float_t distance2Cpv   = fRcpv ;
-  Float_t distance2Track = fRcpv ; 
+  Float_t distance2Track = fRtpc ; 
 
   trackindex = -1 ; // closest track within fRCpv 
 
@@ -151,21 +222,59 @@ Float_t  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint *
     distance2Cpv = (vecCpv - vecEmc).Mag() ;
     
     if (fESD != 0x0) {
+      AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
+      const AliPHOSGeometry * geom = gime->PHOSGeometry() ; 
+
+      Double_t rPHOS = geom->GetIPtoCrystalSurface();
+
+      //PH Acceptance boundaries for each PHOS module
+      Int_t nModules = geom->GetNModules();
+      Double_t * thmin = new Double_t[nModules];// theta min
+      Double_t * thmax = new Double_t[nModules];// theta max
+      Double_t * phmin = new Double_t[nModules];// phi min
+      Double_t * phmax = new Double_t[nModules];// phi max
+      
+      for (Int_t imod=0; imod<nModules; imod++) {
+       geom->EmcModuleCoverage(imod,
+                               thmin[imod],thmax[imod],
+                               phmin[imod],phmax[imod]);
+      }
+
       // Extrapolate the global track direction if any to CPV and find the closest track
       Int_t nTracks = fESD->GetNumberOfTracks();
       Int_t iClosestTrack = -1;
       Double_t minDistance = 1e6;
       Double_t pxyz[3], xyz[3];
+      TVector3 inPHOS; //PH Used to calculate theta and phi
+
+      //PH Loop on tracks
       AliESDtrack *track;
       for (Int_t iTrack=0; iTrack<nTracks; iTrack++) {
        track = fESD->GetTrack(iTrack);
-       if (track->IsPHOS()) 
-         continue ; 
-       track->GetOuterXYZ(xyz);     // track coord on the cylinder of PHOS radius
+       if (!track->GetXYZAt(rPHOS, fESD->GetMagneticField(), xyz))
+           continue; //track coord on the cylinder of PHOS radius
        if ((TMath::Abs(xyz[0])+TMath::Abs(xyz[1])+TMath::Abs(xyz[2]))<=0)
-         continue;
-       track->GetOuterPxPyPz(pxyz); // track momentum ibid.
-       vecDist = PropagateToPlane(xyz,pxyz,"CPV",cpvClu->GetPHOSMod());
+          continue;
+       //PH Here one has to cut out the tracks which are not inside the PHOS
+       //PH acceptance
+       inPHOS.SetXYZ(xyz[0],xyz[1],xyz[2]);
+       Double_t inPhi = inPHOS.Phi();
+       Double_t inTheta = inPHOS.Theta();
+
+       Bool_t skip = kTRUE;
+       for (Int_t imod=0; imod<nModules; imod++) {
+         //PH Loop on modules to check if the track enters in the acceptance 
+         if (thmin[imod] < inTheta && thmax[imod] > inTheta && 
+             phmin[imod] < inPhi   && phmax[imod] > inPhi) {
+           skip = kFALSE;
+           break;
+         }
+       }
+       if (skip) continue; //PH Skip, if not in the PHOS acceptance
+
+       if (!track->GetPxPyPzAt(rPHOS, fESD->GetMagneticField(), pxyz))
+           continue; // track momentum ibid.
+       PropagateToPlane(vecDist,xyz,pxyz,"CPV",cpvClu->GetPHOSMod());
        //      Info("GetDistanceInPHOSPlane","Track %d propagation to CPV = (%f,%f,%f)",
        //     iTrack,vecDist.X(),vecDist.Y(),vecDist.Z());
        vecDist -= vecCpv;
@@ -177,26 +286,31 @@ Float_t  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint *
        }
       }
 
+      delete [] thmin;
+      delete [] thmax;
+      delete [] phmin;
+      delete [] phmax;
+
       if (iClosestTrack != -1) {
        track = fESD->GetTrack(iClosestTrack);
-       track->GetOuterPxPyPz(pxyz); // track momentum ibid.
+       if (track->GetPxPyPzAt(rPHOS, fESD->GetMagneticField(), pxyz)) { // track momentum ibid.
        TVector3 vecCpvGlobal; // Global position of the CPV recpoint
-       AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
-       const AliPHOSGeometry * geom = gime->PHOSGeometry() ; 
        geom->GetGlobal((AliRecPoint*)cpvClu,vecCpvGlobal);
        for (Int_t ixyz=0; ixyz<3; ixyz++)
          xyz[ixyz] = vecCpvGlobal[ixyz];
-       vecDist = PropagateToPlane(xyz,pxyz,"EMC",cpvClu->GetPHOSMod());
+       PropagateToPlane(vecDist,xyz,pxyz,"EMC",cpvClu->GetPHOSMod());
 //     Info("GetDistanceInPHOSPlane","Track %d propagation to EMC = (%f,%f,%f)",
 //          iClosestTrack,vecDist.X(),vecDist.Y(),vecDist.Z());
        vecDist -= vecEmc;
        distance2Track = TMath::Sqrt(vecDist.X()*vecDist.X() + vecDist.Z()*vecDist.Z());
+       }
       }
 //     } else {
 //       // If no ESD exists, than simply find EMC-CPV distance
 //       distance = (vecCpv - vecEmc).Mag() ;
     
-      if(distance2Track < fRcpv + 2*delta )
+      //if(distance2Track < fRcpv + 2*delta )
+      if(distance2Track < fRtpc )
        trackindex = iClosestTrack ; 
       //      toofar = kFALSE ;
     }
@@ -208,8 +322,11 @@ Float_t  AliPHOSTrackSegmentMakerv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint *
 }
 
 //____________________________________________________________________________
-TVector3  AliPHOSTrackSegmentMakerv1::PropagateToPlane(Double_t *x, Double_t *p,
-                                                      char *det, Int_t moduleNumber) const
+void AliPHOSTrackSegmentMakerv1::PropagateToPlane(TVector3& globalIntersection,
+                                                 Double_t *x,
+                                                 Double_t *p,
+                                                 const char *det,
+                                                 Int_t moduleNumber) const
 {
   // Propagate a straight-line track from the origin point x
   // along the direction p to the CPV or EMC module moduleNumber
@@ -217,17 +334,18 @@ TVector3  AliPHOSTrackSegmentMakerv1::PropagateToPlane(Double_t *x, Double_t *p,
 
   AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ; 
-  TVector3 moduleCenter = geom->GetModuleCenter(det,moduleNumber);
-  TVector3 vertex(x);
-  TVector3 direction(p);
+  TVector3 moduleCenter;
+  geom->GetModuleCenter(moduleCenter,det,moduleNumber);
+  TVector3 vertex; vertex.SetXYZ(x[0],x[1],x[2]);
+  TVector3 direction; direction.SetXYZ(p[0],p[1],p[2]);
 
 //   Info("PropagateToCPV","Center of the %s module %d is (%f,%f,%f)",
 //        det,moduleNumber,moduleCenter[0],moduleCenter[1],moduleCenter[2]);
 
   Double_t time = (moduleCenter.Mag2() - vertex.Dot(moduleCenter)) /
     (direction.Dot(moduleCenter));
-  TVector3 globalIntersection = vertex + direction*time;
-  return geom->Global2Local(globalIntersection,moduleNumber);
+  vertex += direction*time;
+  geom->Global2Local(globalIntersection,vertex,moduleNumber);
 }
 
 //____________________________________________________________________________
@@ -235,7 +353,9 @@ void  AliPHOSTrackSegmentMakerv1::Init()
 {
   // Make all memory allocations that are not possible in default constructor
   
-  AliPHOSGetter* gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
+  AliPHOSGetter* gime = AliPHOSGetter::Instance();
+  if(!gime)
+    gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
   
   fLinkUpArray  = new TClonesArray("AliPHOSLink", 1000); 
   if ( !gime->TrackSegmentMaker() ) {
@@ -247,12 +367,14 @@ void  AliPHOSTrackSegmentMakerv1::Init()
 void  AliPHOSTrackSegmentMakerv1::InitParameters()
 {
   //Initializes parameters
-  fRcpv      = 10. ;   
+  fRcpv      = 10. ;
+  fRtpc      = 4. ;
   fEmcFirst  = 0 ;    
   fEmcLast   = 0 ;   
   fCpvFirst  = 0 ;   
   fCpvLast   = 0 ;   
   fLinkUpArray = 0 ;
+  fWrite                   = kTRUE ;
   fTrackSegmentsInRun       = 0 ; 
   SetEventRange(0,-1) ;
 }
@@ -347,7 +469,6 @@ void  AliPHOSTrackSegmentMakerv1::MakePairs()
         
        (dynamic_cast<AliPHOSTrackSegment *>(trackSegments->At(fNTrackSegments)))->SetIndexInList(fNTrackSegments);
        fNTrackSegments++ ;
-       
        emcExist[linkUp->GetEmc()-fEmcFirst] = -1 ; //Mark emc  that Cpv was found 
        //mark CPV recpoint as already used 
        cpvExist[linkUp->GetCpv()-fCpvFirst] = kFALSE ;
@@ -388,18 +509,20 @@ void  AliPHOSTrackSegmentMakerv1::Exec(Option_t *option)
     return ; 
   }
   
-  AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;  
+  AliPHOSGetter * gime = AliPHOSGetter::Instance() ;  
  
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ; 
 
-  if (fLastEvent == -1) fLastEvent = gime->MaxEvent() - 1 ;
-  else fLastEvent = TMath::Min(fFirstEvent,gime->MaxEvent());
+  if (fLastEvent == -1) 
+    fLastEvent = gime->MaxEvent() - 1 ;
+  else 
+    fLastEvent = TMath::Min(fFirstEvent,gime->MaxEvent());
   Int_t nEvents   = fLastEvent - fFirstEvent + 1;
 
   Int_t ievent ; 
   for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
     gime->Event(ievent,"R") ;
-    //Make some initializations 
+   //Make some initializations 
     fNTrackSegments = 0 ;
     fEmcFirst = 0 ;    
     fEmcLast  = 0 ;   
@@ -423,7 +546,6 @@ void  AliPHOSTrackSegmentMakerv1::Exec(Option_t *option)
     
     //increment the total number of track segments per run 
     fTrackSegmentsInRun += gime->TrackSegments()->GetEntriesFast() ; 
-
   }
   
   if(strstr(option,"tim")){
@@ -432,7 +554,8 @@ void  AliPHOSTrackSegmentMakerv1::Exec(Option_t *option)
           gBenchmark->GetCpuTime("PHOSTSMaker"), 
           gBenchmark->GetCpuTime("PHOSTSMaker")/nEvents) ;
    }
-  Unload();
+  if(fWrite) //do not unload in "on flight" mode
+    Unload();
 }
 
 //____________________________________________________________________________
@@ -445,7 +568,7 @@ void AliPHOSTrackSegmentMakerv1::Unload()
 }
 
 //____________________________________________________________________________
-void AliPHOSTrackSegmentMakerv1::Print()const
+void AliPHOSTrackSegmentMakerv1::Print(const Option_t *)const
 {
   //  Print TrackSegmentMaker parameters
 
@@ -478,15 +601,17 @@ void AliPHOSTrackSegmentMakerv1::WriteTrackSegments()
   TClonesArray * trackSegments = gime->TrackSegments() ; 
   trackSegments->Expand(trackSegments->GetEntriesFast()) ;
 
-  TTree * treeT = gime->TreeT();
-  //First TS
-  Int_t bufferSize = 32000 ; 
-  TBranch * tsBranch = treeT->Branch("PHOSTS",&trackSegments,bufferSize);
-  tsBranch->Fill() ;  
-
-  gime->WriteTracks("OVERWRITE");
-  gime->WriteTrackSegmentMaker("OVERWRITE");
+  if(fWrite){ //We write TreeT
+    TTree * treeT = gime->TreeT();
+    
+    //First TS
+    Int_t bufferSize = 32000 ; 
+    TBranch * tsBranch = treeT->Branch("PHOSTS",&trackSegments,bufferSize);
+    tsBranch->Fill() ;  
+    
+    gime->WriteTracks("OVERWRITE");
+    gime->WriteTrackSegmentMaker("OVERWRITE");
+  }
 }