- Many changes in the code to make combined cluster / track finder

  working again
- Using TObjArray instead of TClonesArray to store hits for rec.; two
  function prototypes have 2 extra parameters
(Sasha)

diff --git a/MUON/AliMUONDetElement.cxx b/MUON/AliMUONDetElement.cxx
index 2b0c087aab18e3eccc64a35215e5633ea482bd79..a37781e1b578630930cc39937a04ce43eb529073 100644 (file)
--- a/MUON/AliMUONDetElement.cxx
+++ b/MUON/AliMUONDetElement.cxx
@@ -19,27 +19,34 @@
 // Class AliMUONDetElement
 // --------------------------
 // Detection element object containing information for combined 
-// cluster / track finder in MUON arm 
+// cluster / track finder in the MUON arm 
 // Author: Alexander Zinchenko, JINR Dubna
  
-
-#include <TObjArray.h>
-#include <TClonesArray.h>
 #include "AliMUONDetElement.h"
 #include "AliMUON.h"
 #include "AliMUONSegmentation.h"
 #include "AliMUONDigit.h"
 #include "AliMUONHitMapA1.h"
+#include "AliMUONData.h"
 #include "AliMUONRawCluster.h"
 #include "AliMUONHitForRec.h"
 #include "AliMUONClusterInput.h"
 #include "AliMUONClusterFinderAZ.h"
+#include "AliMUONGeometryModuleTransformer.h" //AZ
+#include "AliMUONVGeometryDESegmentation.h" //AZ
+#include "AliMpVSegmentation.h" //AZ
+
 #include "AliRun.h"
 #include "AliLog.h"
 
+#include <TObjArray.h>
+#include <TClonesArray.h>
+#include <TVector2.h> //AZ
+
 /// \cond CLASSIMP
 ClassImp(AliMUONDetElement) // Class implementation in ROOT context
 /// \endcond
+  FILE *lun = 0x0; //fopen("hitmap.dat","w");
 
 //_________________________________________________________________________
 AliMUONDetElement::AliMUONDetElement()
@@ -53,7 +60,7 @@ AliMUONDetElement::AliMUONDetElement()
     fHitsForRec(0x0),
     fRecModel(0x0)
 {
-/// Default constructor
+  /// Default constructor
   for (Int_t i = 0; i < 2; i++) {
     fHitMap[i] = NULL;
     fDigits[i] = NULL;
@@ -70,43 +77,55 @@ AliMUONDetElement::AliMUONDetElement(Int_t idDE, AliMUONDigit *dig, AliMUONClust
     fZ(0.),
     fNHitsForRec(0),
     fRawClus(0x0),
-    fHitsForRec(0x0),
+    fHitsForRec(new TObjArray(10)),
     fRecModel(recModel)
 {
-/// Constructor
+  /// Constructor
   fDigits[0] = new TObjArray(10);
   fDigits[1] = new TObjArray(10);
   fRawClus = new TObjArray(10);
-  fHitsForRec = new TClonesArray("AliMUONHitForRec",10);
   AliMUON *pMUON = (AliMUON*) gAlice->GetModule("MUON");
   AliMUONSegmentation *pSegmentation = pMUON->GetSegmentation();
-  fSeg[0] = pSegmentation->GetModuleSegmentationByDEId(fidDE, 0);
-  fSeg[1] = pSegmentation->GetModuleSegmentationByDEId(fidDE, 1);
+  fSeg[0] = pSegmentation->GetModuleSegmentationByDEId(idDE, 0);
+  fSeg[1] = pSegmentation->GetModuleSegmentationByDEId(idDE, 1);
+  /*
   Float_t x, y, z;
   fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
   fZ = z;
+  */
   AddDigit(dig);
 }
 
 //_________________________________________________________________________
 AliMUONDetElement::~AliMUONDetElement()
 {
-/// Destructor
+  /// Destructor
+
   for (Int_t i = 0; i < 2; i++) {
     delete fHitMap[i]; fHitMap[i] = NULL; 
     delete fDigits[i]; fDigits[i] = NULL; 
+    //fDigits[i]->Delete(); delete fDigits[i]; fDigits[i] = NULL; 
   }
+  //fHitsForRec->Clear(); delete fHitsForRec; fHitsForRec = 0;
   if (fRawClus) { fRawClus->Delete(); delete fRawClus; fRawClus = 0; }
   //if (fRawClus) { delete fRawClus; fRawClus = 0; }
+  for (Int_t i = 0; i < fNHitsForRec; i++) {
+    Int_t ntracks = ((AliMUONHitForRec*)fHitsForRec->UncheckedAt(i))->GetNTrackHits();
+    //cout << ntracks << endl;
+    if (ntracks) fHitsForRec->RemoveAt(i);
+  }
+  //fHitsForRec->Delete(); 
   delete fHitsForRec; fHitsForRec = 0;
+  //fHitsForRec->SetOwner(kFALSE);
+  //fHitsForRec->Clear(); cout << " here " << fHitsForRec->GetEntriesFast() << endl; delete fHitsForRec; fHitsForRec = 0;
 }
 
 //_________________________________________________________________________
 Int_t AliMUONDetElement::Compare(const TObject* detElem) const
 {
-/// "Compare" function to sort in Z (towards interaction point)
-/// Returns -1 (0, +1) if charge of current pixel
-/// is greater than (equal to, less than) charge of pixel
+  /// "Compare" function to sort in Z (towards interaction point)
+  /// Returns -1 (0, +1) if charge of current pixel
+  /// is greater than (equal to, less than) charge of pixel
   if (fZ > ((AliMUONDetElement*)detElem)->Z()) return(+1);
   else if (fZ == ((AliMUONDetElement*)detElem)->Z()) return( 0);
   else return(-1);
@@ -115,7 +134,7 @@ Int_t AliMUONDetElement::Compare(const TObject* detElem) const
 //_________________________________________________________________________
 void AliMUONDetElement::Fill(AliMUONData */*data*/)
 {
-/// Fill hit maps
+  /// Fill hit maps
   fLeft[0] = fDigits[0]->GetEntriesFast();
   fLeft[1] = fDigits[1]->GetEntriesFast();
 
@@ -150,34 +169,76 @@ void AliMUONDetElement::Fill(AliMUONData */*data*/)
 //_________________________________________________________________________
 void AliMUONDetElement::AddDigit(AliMUONDigit *dig)
 {
-/// Add digit
+  /// Add digit
 
-  fDigits[dig->Cathode()]->Add(dig);
+  //fDigits[dig->Cathode()]->Add(dig);
+  fDigits[dig->Cathode()]->Add(new AliMUONDigit(*dig));
+  Float_t x, y, z;
+  fSeg[dig->Cathode()]->GetPadC(fidDE, dig->PadX(), dig->PadY(), x, y, z);
+  fZ += z;
+  /*
+  Int_t ndig = fDigits[dig->Cathode()]->GetEntriesFast();
+  new((*fDigits[dig->Cathode()])[ndig]) AliMUONDigit(*dig);
+  */
 }
 
 //_________________________________________________________________________
-Bool_t AliMUONDetElement::Inside(Double_t x, Double_t y, Double_t z) const
+Bool_t AliMUONDetElement::Inside(Double_t x, Double_t y, Double_t z, Double_t dx, Double_t dy) const
 {
-/// Check if point is inside detection element
+  /// Check if point is inside detection element
 
+  dx = TMath::Max (dx, 1.);
+  dy = TMath::Max (dy, 1.);
   Int_t ix, iy;
+  ix = iy = 0;
+  Double_t xl, yl, zl;
   for (Int_t i = 0; i < 2; i++) {
     if (!fSeg[i]) continue;
+    if (i == 0 || fSeg[0] == 0x0) {
+      fSeg[i]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
+      //cout << xl << " " << yl << " " << zl << endl;
+      if (TMath::Abs (zl) > 0.5) return kFALSE;
+      //if (fChamber < 4 && (xl < 0 || yl < 0)) return kFALSE;
+    }
     fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy);
-    //cout << x << " " << y << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << endl; 
-    if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
+    //cout << x << " " << y << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z) << endl; 
+    //if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
+    if (fSeg[i]->HasPad(fidDE, (Float_t)x, (Float_t)y, (Float_t)z)) return kTRUE;
   }
-  // Check for edge effect (extrapolated track "right outside" det. elem. boundaries (+- 1cm in X and Y)
+  //Double_t xl, yl, zl;
+  //fSeg[0]->GetTransformer()->Global2Local(fidDE, x, y, z, xl, yl, zl);
+  //cout << xl << " " << yl << " " << zl << endl;
+
+  Int_t cath = fSeg[0] ? 0 : 1;
+  const AliMpVSegmentation *mpSeg = fSeg[cath]->GetDESegmentation(fidDE)->GetMpSegmentation();
+  TVector2 dim = mpSeg->Dimensions(); 
+  if (fChamber < 4) dim *= 2.; // quadrants
+  //cout << fChamber << " " << dim.X() << " " << dim.Y() << endl;
+  if (TMath::Abs (xl) - dim.X() > dx) return kFALSE;
+  if (TMath::Abs (yl) - dim.Y() > dy) return kFALSE;
+
+  Int_t flagX = 0, flagY = 0;
+  if (TMath::Abs (xl) - dim.X() < -2) flagX = 1;
+  if (TMath::Abs (yl) - dim.Y() < -2) flagY = 1;
+  if (flagX && flagY) flagX = flagY = 0; // both coordinates "inside" det. elem. limits
+
+  // Check for edge effect (extrapolated track "right outside" det. elem. boundaries 
+  // (+- dx * cath in X and Y)
+  dx /= 2;
+  dy /= 2;
   for (Int_t i = 0; i < 2; i++) {
     if (!fSeg[i]) continue;
     for (Int_t idx = -1; idx < 2; idx++) {
-      Double_t x1 = x + 1. * idx;
+      if (flagX && idx) continue;
+      Double_t x1 = x + dx * idx * (i + 1);
       for (Int_t idy = -1; idy < 2; idy++) {
        if (idx == 0 && idy == 0) continue;
-       Double_t y1 = y + 1. * idy;
+       if (flagY && idy) continue;
+       Double_t y1 = y + dy * idy * (i + 1);
        fSeg[i]->GetPadI(fidDE, x1, y1, z, ix, iy);
-       //cout << x1 << " " << y1 << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << endl; 
-       if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
+       //cout << x1 << " " << y1 << " " << z << " " << fChamber << " " << ix << " " << iy << " " << fSeg[i]->Npx(fidDE) << " " << fSeg[i]->Npy(fidDE) /*<< " " << fSeg[i]->GetPadI(fidDE, x, y, z, ix, iy)*/ << " " << fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z) << endl; 
+       //if (ix > 0 && iy > 0 && ix <= fSeg[i]->Npx(fidDE) && iy <= fSeg[i]->Npy(fidDE)) return kTRUE;
+       if (fSeg[i]->HasPad(fidDE, (Float_t)x1, (Float_t)y1, (Float_t)z)) return kTRUE;
       }
     }
   }
@@ -185,23 +246,23 @@ Bool_t AliMUONDetElement::Inside(Double_t x, Double_t y, Double_t z) const
 }
 
 //_________________________________________________________________________
-void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
+void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack, Double_t dx, Double_t dy)
 {
-/// Run cluster reconstruction around point (x,y)
+  /// Run cluster reconstruction around point (x,y) inside window (dx,dy)
 
   if (fLeft[0] == 0 && fLeft[1] == 0) return; // all digits have been used 
-  Float_t dx, dy;
-  dx = dy = 5; // 5 cm for now 
   AliMUONClusterInput::Instance()->SetDigits(fChamber, fidDE, 
              (TClonesArray*)fDigits[0], (TClonesArray*)fDigits[1]);
 
   // Mark used pads
   for (Int_t cath = 0; cath < 2; cath++) {
-    if (fDigits[cath]->GetEntriesFast() == 0) continue; // empty cathode
-    for (Int_t i = 0; i < fDigits[cath]->GetEntriesFast(); i++) {
+    Int_t ndig = fDigits[cath]->GetEntriesFast();
+    if (ndig == 0) continue; // empty cathode
+    //cout << fidDE << " " << cath << " " << ndig << endl; 
+    for (Int_t i = 0; i < ndig; i++) {
       if (fLeft[cath] == 0) { fRecModel->SetUsed(cath,i); continue; }
       AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath]->UncheckedAt(i);
-      //cout << i << " " << dig->PadX() << " " << dig->PadY() << " " << fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) << endl;
+      //cout << i << " " << dig << " " /*<< dig->PadX() << " " << dig->PadY() << " " << fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) */<< endl;
       if (fHitMap[cath]->TestHit(dig->PadX(), dig->PadY()) == kUsed) fRecModel->SetUsed(cath,i);
       else fRecModel->SetUnused(cath,i);
     }
@@ -210,7 +271,8 @@ void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
   fRecModel->ResetRawClusters();
 
   for (Int_t cath = 0; cath < 2; cath++) {
-    if (fDigits[cath]->GetEntriesFast() == 0) continue; // empty cathode
+    Int_t ndig = fDigits[cath]->GetEntriesFast();
+    if (ndig == 0) continue; // empty cathode
     // Loop over pads
     for (fSeg[cath]->FirstPad(fidDE, xTrack, yTrack, fZ, dx, dy);
         fSeg[cath]->MorePads(fidDE);
@@ -220,14 +282,17 @@ void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
       if (fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kEmpty ||
          fHitMap[cath]->TestHit(fSeg[cath]->Ix(), fSeg[cath]->Iy()) == kUsed) continue;
 
+      Int_t iOK = 0;
       // Set starting pad
-      for (Int_t j = 0; j < fDigits[cath]->GetEntriesFast(); j++) {
+      for (Int_t j = 0; j < ndig; j++) {
        AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath]->UncheckedAt(j);
        if (dig->PadX() != fSeg[cath]->Ix() || dig->PadY() != fSeg[cath]->Iy()) continue;
-       //cout << fidDE << " " << j << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << endl;
+       //cout << fidDE << " " << j << " " << fDigits[cath]->GetEntriesFast() << " " << fSeg[cath]->Ix() << " " << fSeg[cath]->Iy() << endl;
        fRecModel->SetStart(cath, j);
+       iOK = 1;
        break;
       }
+      if (!iOK) continue;
 
       fRecModel->FindRawClusters();
       Int_t nClusEnd = fRecModel->GetRawClusters()->GetEntriesFast();
@@ -239,7 +304,8 @@ void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
       }
       // Mark used pads
       for (Int_t cath1 = 0; cath1 < 2; cath1++) {
-       for (Int_t j = 0; j < fDigits[cath1]->GetEntriesFast(); j++) {
+       Int_t ndig1 = fDigits[cath1]->GetEntriesFast();
+       for (Int_t j = 0; j < ndig1; j++) {
          if (fLeft[cath1] == 0) break;
          AliMUONDigit *dig = (AliMUONDigit*) fDigits[cath1]->UncheckedAt(j);
          Float_t x, y, z;
@@ -248,7 +314,7 @@ void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
          if (!fRecModel->GetUsed(cath1,j)) continue;
          if (fHitMap[cath1]->TestHit(dig->PadX(), dig->PadY()) == kUsed) continue;
          fHitMap[cath1]->FlagHit(dig->PadX(), dig->PadY());
-         AliDebug(10,Form(" %d %d %d %d \n", cath1, fidDE, dig->PadX(), dig->PadY()));
+         if (lun) fprintf(lun," %d %d %d %d \n", cath1, fidDE, dig->PadX(), dig->PadY());
          fLeft[cath1]--;
        }
       }
@@ -259,16 +325,21 @@ void AliMUONDetElement::ClusterReco(Double_t xTrack, Double_t yTrack)
 //_________________________________________________________________________
 void AliMUONDetElement::AddHitForRec(AliMUONRawCluster *clus)
 {
-/// Make HitForRec from raw cluster (rec. point)
+  /// Make HitForRec from raw cluster (rec. point)
 
   fRawClus->Add(new AliMUONRawCluster(*clus));
-  AliMUONHitForRec *hitForRec = 
-    new ((*fHitsForRec)[fNHitsForRec++]) AliMUONHitForRec(clus);
+  //AliMUONHitForRec *hitForRec = 
+  //new ((*fHitsForRec)[fNHitsForRec++]) AliMUONHitForRec(clus);
+  AliMUONHitForRec *hitForRec = new AliMUONHitForRec(clus);
+  fHitsForRec->Add(hitForRec);
+  fNHitsForRec++;
 
   // more information into HitForRec
   //  resolution: info should be already in raw cluster and taken from it ????
-  hitForRec->SetBendingReso2(-1); //fBendingResolution * fBendingResolution);
-  hitForRec->SetNonBendingReso2(-1); //fNonBendingResolution * fNonBendingResolution);
+  //hitForRec->SetBendingReso2(-1); //fBendingResolution * fBendingResolution);
+  //hitForRec->SetNonBendingReso2(-1); //fNonBendingResolution * fNonBendingResolution);
+  hitForRec->SetBendingReso2(clus->GetErrY() * clus->GetErrY());
+  hitForRec->SetNonBendingReso2(clus->GetErrX() * clus->GetErrX());
   //  original raw cluster
   hitForRec->SetChamberNumber(fChamber);
   hitForRec->SetZ(clus->GetZ(0));

diff --git a/MUON/AliMUONDetElement.h b/MUON/AliMUONDetElement.h
index ac97f9ed86febeba049877ad3e63ea7687e87095..09739114177129b4bf6afb9509ebe90294278b01 100644 (file)
--- a/MUON/AliMUONDetElement.h
+++ b/MUON/AliMUONDetElement.h
@@ -8,7 +8,7 @@
 /// \ingroup rec
 /// \class AliMUONDetElement
 /// \brief Detection element object containing information for combined 
-/// cluster / track finder in MUON arm 
+/// cluster / track finder in the MUON arm 
 ///
 /// \author Alexander Zinchenko, JINR Dubna
  
@@ -37,9 +37,10 @@ class AliMUONDetElement : public TObject
   //Int_t GetMapElem(AliMUONDigit *digit) const; // get map element
   TObjArray *Digits(Int_t cath) const { return fDigits[cath]; } // array of digits
   TObjArray *RawClusters() const { return fRawClus; } // array of raw clusters
-  TClonesArray *HitsForRec() const { return fHitsForRec; } // hits for rec.
+  TObjArray *HitsForRec() const { return fHitsForRec; } // hits for rec.
   Int_t NHitsForRec() const { return fNHitsForRec; } // No. of hits for rec.
-  Bool_t Inside(Double_t x, Double_t y, Double_t z) const; // check if point inside DE
+  //Bool_t Inside(Double_t x, Double_t y, Double_t z) const; // check if point inside DE
+  Bool_t Inside(Double_t x, Double_t y, Double_t z, Double_t dx, Double_t dy) const; // check if point inside DE
   
   void SetID(Int_t idDE) { fidDE = idDE; } // set det. elem. ID
   void SetIndex(Int_t index) { fIndex = index; } // set position index
@@ -48,7 +49,7 @@ class AliMUONDetElement : public TObject
   //void SetMapElem(const AliMUONDigit *digit, Int_t flag); // set map element 
   void AddDigit(AliMUONDigit *dig); // add digit
   void Fill(AliMUONData *data); // fill hit maps 
-  void ClusterReco(Double_t x, Double_t y); // run cluster reco around (x,y)
+  void ClusterReco(Double_t x, Double_t y, Double_t dx, Double_t dy); // run cluster reco around (x,y)
   void AddHitForRec(AliMUONRawCluster *clus); // make HitForRec
   // What is necessary for sorting TObjArray's
   Bool_t IsSortable() const { return kTRUE; }
@@ -66,7 +67,7 @@ class AliMUONDetElement : public TObject
   AliMUONHitMapA1 *fHitMap[2]; ///< map of digits
   TObjArray *fDigits[2]; ///< container of digits from this det. elem.
   TObjArray *fRawClus; ///< raw clusters
-  TClonesArray *fHitsForRec; ///< HitForRec's
+  TObjArray *fHitsForRec; ///< HitForRec's
   AliMUONClusterFinderAZ *fRecModel; ///< cluster finder
 
   // Functions