Fixes for #86059: Install data when ALICE_ROOT does not point to source (Christian)
[u/mrichter/AliRoot.git] / MUON / AliMUONClusterFinderMLEM.cxx
index d863182..21d26f5 100644 (file)
@@ -51,8 +51,6 @@
 ClassImp(AliMUONClusterFinderMLEM)
 /// \endcond
  
-const Double_t AliMUONClusterFinderMLEM::fgkZeroSuppression = 6; // average zero suppression value
-//const Double_t AliMUONClusterFinderMLEM::fgkDistancePrecision = 1e-6; // (cm) used to check overlaps and so on
 const Double_t AliMUONClusterFinderMLEM::fgkDistancePrecision = 1e-3; // (cm) used to check overlaps and so on
 const TVector2 AliMUONClusterFinderMLEM::fgkIncreaseSize(-AliMUONClusterFinderMLEM::fgkDistancePrecision,-AliMUONClusterFinderMLEM::fgkDistancePrecision);
 const TVector2 AliMUONClusterFinderMLEM::fgkDecreaseSize(AliMUONClusterFinderMLEM::fgkDistancePrecision,AliMUONClusterFinderMLEM::fgkDistancePrecision);
@@ -81,10 +79,13 @@ fDebug(0),
 fPlot(plot),
 fSplitter(0x0),
 fNClusters(0),
-fNAddVirtualPads(0)
+fNAddVirtualPads(0),
+fLowestPixelCharge(0),
+fLowestPadCharge(0),
+fLowestClusterCharge(0)
 {
   /// Constructor
+
   fkSegmentation[1] = fkSegmentation[0] = 0x0; 
 
   if (fPlot) fDebug = 1;
@@ -105,7 +106,7 @@ AliMUONClusterFinderMLEM::~AliMUONClusterFinderMLEM()
 //_____________________________________________________________________________
 Bool_t 
 AliMUONClusterFinderMLEM::Prepare(Int_t detElemId,
-                                  TClonesArray* pads[2],
+                                  TObjArray* pads[2],
                                   const AliMpArea& area,
                                   const AliMpVSegmentation* seg[2])
 {
@@ -121,7 +122,11 @@ AliMUONClusterFinderMLEM::Prepare(Int_t detElemId,
   fDetElemId = detElemId;
   
   delete fSplitter;
-  fSplitter = new AliMUONClusterSplitterMLEM(fDetElemId,fPixArray);
+  fSplitter = new AliMUONClusterSplitterMLEM(fDetElemId,
+                                             fPixArray,
+                                             fLowestPixelCharge,
+                                             fLowestPadCharge,
+                                             fLowestClusterCharge);
   fSplitter->SetDebug(fDebug);
     
   // find out current event number, and reset the cluster number
@@ -129,7 +134,8 @@ AliMUONClusterFinderMLEM::Prepare(Int_t detElemId,
   fEventNumber = runLoader ? runLoader->GetEventNumber() : 0;
   fClusterNumber = -1;
   fClusterList.Delete();
-  
+  fPixArray->Delete();
+
   AliDebug(3,Form("EVT %d DE %d",fEventNumber,fDetElemId));
   
   if ( fPreClusterFinder->NeedSegmentation() )
@@ -150,7 +156,14 @@ AliMUONClusterFinderMLEM::NextCluster()
 //  AliCodeTimerAuto("",0)
   
   // if the list of clusters is not void, pick one from there
-  TObject* o = fClusterList.At(++fClusterNumber);
+  TObject* o(0x0);
+  
+  // do we have clusters in our list ?
+  if ( fClusterNumber < fClusterList.GetLast() )
+  {
+    o = fClusterList.At(++fClusterNumber);
+  }
+  
   if ( o != 0x0 ) return static_cast<AliMUONCluster*>(o);
   
   //FIXME : at this point, must check whether we've used all the digits
@@ -162,6 +175,7 @@ AliMUONClusterFinderMLEM::NextCluster()
 
   fPreCluster = fPreClusterFinder->NextCluster();
 
+  fPixArray->Delete();
   fClusterList.Delete(); // reset the list of clusters for this pre-cluster
   fClusterNumber = -1; //AZ
     
@@ -297,7 +311,7 @@ AliMUONClusterFinderMLEM::CheckPrecluster(const AliMUONCluster& origCluster)
 
   // Disregard small clusters (leftovers from splitting or noise)
   if ((origCluster.Multiplicity()==1 || origCluster.Multiplicity()==2) &&
-      origCluster.Charge(0)+origCluster.Charge(1) < 10) 
+      origCluster.Charge(0)+origCluster.Charge(1) < fLowestClusterCharge )
   { 
     return 0x0;
   }
@@ -367,8 +381,7 @@ AliMUONClusterFinderMLEM::CheckPreclusterTwoCathodes(AliMUONCluster* cluster)
       fkSegmentation[cath1]->PadByPosition(pad->Position().X(),
                                            pad->Position().Y(),kFALSE);
       if (!mpPad.IsValid()) continue;
-      //if (nFlags == 1 && pad->Charge() < fgkZeroSuppression * 3) continue;
-      if (nFlags == 1 && pad->Charge() < 20) continue;
+      if (nFlags == 1 && pad->Charge() < fLowestPadCharge) continue; 
       AliDebug(2,Form("Releasing the following pad : de,cath,ix,iy %d,%d,%d,%d charge %e",
                       fDetElemId,pad->Cathode(),pad->Ix(),pad->Iy(),pad->Charge()));
       toBeRemoved.AddLast(pad);
@@ -739,20 +752,20 @@ void AliMUONClusterFinderMLEM::PadOverHist(Int_t idir, Int_t ix0, Int_t iy0, Ali
 
 //_____________________________________________________________________________
 void
-AliMUONClusterFinderMLEM::Plot(const char* basename)
+AliMUONClusterFinderMLEM::Plot(const char* /*basename*/)
 {
   /// Make a plot and save it as png
   
   return; //AZ
-  if (!fPlot) return;
-  
-  TCanvas* c = new TCanvas("MLEM","MLEM",800,600);
-  c->Draw();
-  Draw();
-  c->Modified();
-  c->Update();
-  c->Print(Form("%s.EVT%d.DE%d.CLU%d.png",basename,fEventNumber,
-                fDetElemId,fClusterNumber));
+//  if (!fPlot) return;
+//  
+//  TCanvas* c = new TCanvas("MLEM","MLEM",800,600);
+//  c->Draw();
+//  Draw();
+//  c->Modified();
+//  c->Update();
+//  c->Print(Form("%s.EVT%d.DE%d.CLU%d.png",basename,fEventNumber,
+//                fDetElemId,fClusterNumber));
 }
 
 //_____________________________________________________________________________
@@ -841,10 +854,11 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
   {
     ++lc;
     delete fHistMlem;
+    fHistMlem = 0x0;
     
     AliDebug(2,Form("lc %d nPix %d(%d) npadTot %d npadOK %d",lc,nPix,fPixArray->GetLast()+1,npadTot,npadOK));
     AliDebug(2,Form("EVT%d PixArray=",fEventNumber));
-    //StdoutToAliDebug(2,fPixArray->Print("","full"));
+    //StdoutToAliDebug(2,fPixArray->Print("full"));
         
     coef = new Double_t [npadTot*nPix];
     probi = new Double_t [nPix];
@@ -894,6 +908,10 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
                     xylim[0],-xylim[1],xylim[2],-xylim[3]
                     ));
     
+    AliDebug(2,Form("LowestPadCharge=%e",fLowestPadCharge));
+    
+    delete fHistMlem;
+    
     fHistMlem = new TH2D("mlem","mlem",nx,xylim[0],-xylim[1],ny,xylim[2],-xylim[3]);
 
     for (Int_t ipix = 0; ipix < nPix; ++ipix) 
@@ -909,7 +927,7 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
       qTot += Pixel(i)->Charge();
     }
     
-    if ( qTot < 1.e-4 || ( npadOK < 3 && qTot < 7 ) ) 
+    if ( qTot < 1.e-4 || ( npadOK < 3 && qTot < fLowestClusterCharge ) )
     {
       AliDebug(1,Form("Deleting the above cluster (charge %e too low, npadOK=%d)",qTot,npadOK));
       delete [] coef; 
@@ -946,7 +964,7 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
     fPixArray->Sort();
     MaskPeaks(0); // unmask local maxima
     Double_t pixMin = 0.01*Pixel(0)->Charge();
-    pixMin = TMath::Min(pixMin,50.);
+    pixMin = TMath::Min(pixMin,100*fLowestPixelCharge);
 
     // Decrease pixel size and shift pixels to make them centered at 
     // the maximum one
@@ -972,7 +990,7 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
       {
         if (!i) 
         {
-          pixPtr2->SetCharge(10);
+          pixPtr2->SetCharge(fLowestPadCharge);
           pixPtr2->SetSize(indx, pixPtr2->Size(indx)/2);
           width = -pixPtr2->Size(indx);
           pixPtr2->Shift(indx, width);
@@ -1052,8 +1070,8 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
 
   // remove pixels with low signal or low visibility
   // Cuts are empirical !!!
-  Double_t thresh = TMath::Max (fHistMlem->GetMaximum()/100.,1.);
-  thresh = TMath::Min (thresh,50.);
+  Double_t thresh = TMath::Max (fHistMlem->GetMaximum()/100.,2.0*fLowestPixelCharge);
+  thresh = TMath::Min (thresh,100.0*fLowestPixelCharge);
   Double_t charge = 0;
 
   // Mark pixels which should be removed
@@ -1097,7 +1115,7 @@ Bool_t AliMUONClusterFinderMLEM::MainLoop(AliMUONCluster& cluster, Int_t iSimple
 
   // Try to split into clusters
   Bool_t ok = kTRUE;
-  if (fHistMlem->GetSum() < 1) 
+  if (fHistMlem->GetSum() < 2.0*fLowestPixelCharge) 
   {
     ok = kFALSE;
   }
@@ -1296,7 +1314,7 @@ Int_t AliMUONClusterFinderMLEM::FindNearest(const AliMUONPad *pixPtr0)
 
   for (Int_t i = 0; i < nPix; ++i) {
     pixPtr = (AliMUONPad*) fPixArray->UncheckedAt(i);
-    if (pixPtr == pixPtr0 || pixPtr->Charge() < 0.5) continue;
+    if (pixPtr == pixPtr0 || pixPtr->Charge() < fLowestPixelCharge) continue;
     dx = (xc - pixPtr->Coord(0)) / pixPtr->Size(0);
     dy = (yc - pixPtr->Coord(1)) / pixPtr->Size(1);
     r = dx *dx + dy * dy;
@@ -1376,6 +1394,9 @@ Int_t AliMUONClusterFinderMLEM::FindLocalMaxima(TObjArray *pixArray, Int_t *loca
   Double_t xylim[4] = {999, 999, 999, 999};
 
   Int_t nPix = pixArray->GetEntriesFast();
+  
+  if ( nPix <= 0 ) return 0;
+  
   AliMUONPad *pixPtr = 0;
   for (Int_t ipix = 0; ipix < nPix; ++ipix) {
     pixPtr = (AliMUONPad*) pixArray->UncheckedAt(ipix);
@@ -1401,7 +1422,7 @@ Int_t AliMUONClusterFinderMLEM::FindLocalMaxima(TObjArray *pixArray, Int_t *loca
   for (Int_t i = 1; i <= ny; ++i) {
     indx = (i-1) * nx;
     for (Int_t j = 1; j <= nx; ++j) {
-      if (fHistAnode->GetCellContent(j,i) < 0.5) continue;
+      if (fHistAnode->GetCellContent(j,i) < fLowestPixelCharge) continue;
       //if (isLocalMax[indx+j-1] < 0) continue;
       if (isLocalMax[indx+j-1] != 0) continue;
       FlagLocalMax(fHistAnode, i, j, isLocalMax);
@@ -1551,7 +1572,7 @@ AliMUONClusterFinderMLEM::AddBinSimple(TH2D *hist, Int_t ic, Int_t jc)
     for (Int_t j = TMath::Max(jc-1,1); j <= je; ++j) {
       cont1 = hist->GetCellContent(j,i);
       if (cont1 > cont) continue;
-      if (cont1 < 0.5) continue;
+      if (cont1 < fLowestPixelCharge) continue;
       pixPtr = new AliMUONPad (hist->GetXaxis()->GetBinCenter(j), 
                               hist->GetYaxis()->GetBinCenter(i), 0, 0, cont1);
       fPixArray->Add(pixPtr);
@@ -1671,10 +1692,10 @@ void AliMUONClusterFinderMLEM::AddVirtualPad(AliMUONCluster& cluster)
       AliMUONPad muonPad(fDetElemId, nonb[inb], mppad.GetIx(), mppad.GetIy(), 
                         mppad.GetPositionX(), mppad.GetPositionY(), 
                         mppad.GetDimensionX(), mppad.GetDimensionY(), 0);
-      if (inb == 0) muonPad.SetCharge(TMath::Min (amax[j]/100, 5.));
+      if (inb == 0) muonPad.SetCharge(TMath::Min (amax[j]/100, fLowestPadCharge));
       //else muonPad.SetCharge(TMath::Min (amax[j]/15, fgkZeroSuppression));
-      else muonPad.SetCharge(TMath::Min (amax[j]/15, 6.));
-      if (muonPad.Charge() < 1.) muonPad.SetCharge(1.);
+      else muonPad.SetCharge(TMath::Min (amax[j]/15, fLowestPadCharge));
+      if (muonPad.Charge() < 2.0*fLowestPixelCharge) muonPad.SetCharge(2.0*fLowestPixelCharge);
       muonPad.SetReal(kFALSE);
       if (fDebug) printf(" ***** Add virtual pad in %d direction ***** %f %f %f %3d %3d %f %f \n",
                         inb, muonPad.Charge(), muonPad.X(), muonPad.Y(), muonPad.Ix(), 
@@ -1768,4 +1789,14 @@ AliMUONClusterFinderMLEM::Print(Option_t* what) const
   }
 }
 
+//_____________________________________________________________________________
+void 
+AliMUONClusterFinderMLEM::SetChargeHints(Double_t lowestPadCharge, Double_t lowestClusterCharge)
+{
+  /// Set some thresholds we use later on in the algorithm
+  fLowestPadCharge=lowestPadCharge;
+  fLowestClusterCharge=lowestClusterCharge;
+  fLowestPixelCharge=fLowestPadCharge/12.0; 
+}
+