/* $Id$ */
-// Clusterizer class developed by A. Zinchenko (Dubna), based on the
-// Expectation-Maximization algorithm
+// -------------------------------
+// Class AliMUONClusterFinderAZ
+// -------------------------------
+// Clusterizer class based on the Expectation-Maximization algorithm
+// Author: Alexander Zinchenko, JINR Dubna
#include <stdlib.h>
#include <Riostream.h>
#include <TH2.h>
#include <TMinuit.h>
#include <TMatrixD.h>
+#include <TRandom.h>
+#include <TROOT.h>
#include "AliMUONClusterFinderAZ.h"
#include "AliMUONClusterDrawAZ.h"
#include "AliMUONVGeometryDESegmentation.h"
#include "AliMUONGeometryModuleTransformer.h"
-#include "AliRun.h"
-#include "AliMUON.h"
+#include "AliMUONGeometrySegmentation.h"
#include "AliMUONDigit.h"
#include "AliMUONRawCluster.h"
#include "AliMUONClusterInput.h"
#include "AliMUONMathieson.h"
#include "AliLog.h"
+/// \cond CLASSIMP
ClassImp(AliMUONClusterFinderAZ)
+/// \endcond
const Double_t AliMUONClusterFinderAZ::fgkCouplMin = 1.e-3; // threshold on coupling
const Double_t AliMUONClusterFinderAZ::fgkZeroSuppression = 6; // average zero suppression value
//_____________________________________________________________________________
AliMUONClusterFinderAZ::AliMUONClusterFinderAZ(Bool_t draw)
- : AliMUONClusterFinderVS()
+ : AliMUONClusterFinderVS(),
+ fZpad(0),
+ fNpar(0),
+ fQtot(0),
+ fReco(1),
+ fCathBeg(0),
+ fDraw(0x0),
+ fPixArray(0x0),
+ fnCoupled(0),
+ fDebug(0)
{
-// Constructor
+/// Constructor
fnPads[0]=fnPads[1]=0;
for (Int_t i=0; i<7; i++)
fSegmentation[1] = fSegmentation[0] = 0x0;
- fZpad = 0;
- fQtot = 0;
- fPadBeg[0] = fPadBeg[1] = fCathBeg = fNpar = fnCoupled = 0;
+ fPadBeg[0] = fPadBeg[1] = 0;
if (!fgMinuit) fgMinuit = new TMinuit(8);
if (!fgClusterFinder) fgClusterFinder = this;
fPixArray = new TObjArray(20);
- fDebug = 0; //0;
- fReco = 1;
- fDraw = 0x0;
if (draw) {
fDebug = 1;
fReco = 0;
cout << " *** Running AZ cluster finder *** " << endl;
}
-//_____________________________________________________________________________
-AliMUONClusterFinderAZ::AliMUONClusterFinderAZ(const AliMUONClusterFinderAZ& rhs)
- : AliMUONClusterFinderVS(rhs)
-{
-// Protected copy constructor
-
- AliFatal("Not implemented.");
-}
-
//_____________________________________________________________________________
AliMUONClusterFinderAZ::~AliMUONClusterFinderAZ()
{
- // Destructor
+/// Destructor
delete fgMinuit; fgMinuit = 0; delete fPixArray; fPixArray = 0;
delete fDraw;
}
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::FindRawClusters()
{
-// To provide the same interface as in AliMUONClusterFinderVS
+/// To provide the same interface as in AliMUONClusterFinderVS
ResetRawClusters();
- EventLoop (gAlice->GetEvNumber(), fInput->Chamber());
+ EventLoop (fEvtNumber, fInput->Chamber());
}
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::EventLoop(Int_t nev, Int_t ch)
{
-// Loop over digits
+/// Loop over digits
if (fDraw && !fDraw->FindEvCh(nev, ch)) return;
*/
for (Int_t i=0; i<nMax; i++) {
if (nMax > 1) FindCluster(localMax, maxPos[i]);
- if (!MainLoop(iSimple)) AliWarning(Form(" MainLoop failed "));
+ MainLoop(iSimple);
if (i < nMax-1) {
for (Int_t j=0; j<fnPads[0]+fnPads[1]; j++) {
if (fPadIJ[1][j] == 0) continue; // pad charge was not modified
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AddPad(Int_t cath, Int_t digit)
{
- // Add pad to the cluster
+/// Add pad to the cluster
+
AliMUONDigit *mdig = fInput->Digit(cath,digit);
- Int_t charge = mdig->Signal();
+ Float_t charge = mdig->Signal();
// get the center of the pad
Float_t xpad, ypad, zpad0;
//if (!fSegmentation[cath]->GetPadC(fInput->DetElemId(),mdig->PadX(),mdig->PadY(),xpad,ypad,zpad0)) { // Handle "non-existing" pads
fPadIJ[2][nPads] = mdig->PadX();
fPadIJ[3][nPads] = mdig->PadY();
fUsed[cath][digit] = kTRUE;
- if (fDebug) printf(" bbb %d %d %f %f %f %f %f %4d %3d %3d \n", nPads, cath, xpad, ypad, zpad0, fXyq[3][nPads]*2, fXyq[4][nPads]*2, charge, mdig->PadX(), mdig->PadY());
+ if (fDebug) printf(" bbb %d %d %f %f %f %f %f %f %3d %3d \n", nPads, cath,
+ xpad, ypad, zpad0, fXyq[3][nPads]*2, fXyq[4][nPads]*2,
+ charge, mdig->PadX(), mdig->PadY());
fnPads[cath]++;
// Check neighbours
//_____________________________________________________________________________
Bool_t AliMUONClusterFinderAZ::Overlap(Int_t cath, AliMUONDigit *mdig)
{
- // Check if the pad from one cathode overlaps with a pad
- // in the precluster on the other cathode
+/// Check if the pad from one cathode overlaps with a pad
+/// in the precluster on the other cathode
Float_t xpad, ypad, zpad;
fSegmentation[cath]->GetPadC(mdig->PadX(), mdig->PadY(), xpad, ypad, zpad);
//_____________________________________________________________________________
Bool_t AliMUONClusterFinderAZ::Overlap(Float_t *xy1, Int_t iPad, Float_t *xy12, Int_t iSkip)
{
- // Check if the pads xy1 and iPad overlap and return overlap area
+/// Check if the pads xy1 and iPad overlap and return overlap area
Float_t xy2[4];
xy2[0] = fXyq[0][iPad] - fXyq[3][iPad];
//_____________________________________________________________________________
Bool_t AliMUONClusterFinderAZ::CheckPrecluster(Int_t *nShown)
{
- // Check precluster in order to attempt to simplify it (mostly for
- // two-cathode preclusters)
+/// Check precluster in order to attempt to simplify it (mostly for
+/// two-cathode preclusters)
Int_t i1, i2, cath=0, digit=0;
Float_t xy1[4], xy12[4];
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::BuildPixArray()
{
- // Build pixel array for MLEM method
+/// Build pixel array for MLEM method
Int_t nPix=0, i1, i2;
Float_t xy1[4], xy12[4];
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AdjustPixel(Float_t width, Int_t ixy)
{
- // Check if some pixels have small size (adjust if necessary)
+/// Check if some pixels have small size (adjust if necessary)
AliMUONPixel *pixPtr, *pixPtr1 = 0;
Int_t ixy1 = TMath::Even(ixy);
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AdjustPixel(Float_t wxmin, Float_t wymin)
{
- // Check if some pixels have large size (adjust if necessary)
+/// Check if some pixels have large size (adjust if necessary)
Int_t n1[2], n2[2], iOK = 1, nPix = fPixArray->GetEntriesFast();
AliMUONPixel *pixPtr, pix;
- Double_t xy0[2] = {9999, 9999}, wxy[2], dist[2];
+ Double_t xy0[2] = {9999, 9999}, wxy[2], dist[2] = {0};
// Check if large pixel size
for (Int_t i = 0; i < nPix; i++) {
//_____________________________________________________________________________
Bool_t AliMUONClusterFinderAZ::MainLoop(Int_t iSimple)
{
- // Repeat MLEM algorithm until pixel size becomes sufficiently small
+/// Repeat MLEM algorithm until pixel size becomes sufficiently small
TH2D *mlem;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Mlem(Double_t *coef, Double_t *probi, Int_t nIter)
{
- // Use MLEM to find pixel charges
+/// Use MLEM to find pixel charges
Int_t nPix = fPixArray->GetEntriesFast();
Int_t npad = fnPads[0] + fnPads[1];
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::FindCOG(TH2D *mlem, Double_t *xyc)
{
- // Calculate position of the center-of-gravity around the maximum pixel
+/// Calculate position of the center-of-gravity around the maximum pixel
Int_t ixmax, iymax, ix, nsumx=0, nsumy=0, nsum=0;
Int_t i1 = -9, j1 = -9;
//_____________________________________________________________________________
Int_t AliMUONClusterFinderAZ::FindNearest(AliMUONPixel *pixPtr0)
{
- // Find the pixel nearest to the given one
- // (algorithm may be not very efficient)
+/// Find the pixel nearest to the given one
+/// (algorithm may be not very efficient)
Int_t nPix = fPixArray->GetEntriesFast(), imin = 0;
Double_t rmin = 99999, dx = 0, dy = 0, r = 0;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Split(TH2D *mlem, Double_t *coef)
{
- // The main steering function to work with clusters of pixels in anode
- // plane (find clusters, decouple them from each other, merge them (if
- // necessary), pick up coupled pads, call the fitting function)
+/// The main steering function to work with clusters of pixels in anode
+/// plane (find clusters, decouple them from each other, merge them (if
+/// necessary), pick up coupled pads, call the fitting function)
Int_t nx = mlem->GetNbinsX();
Int_t ny = mlem->GetNbinsY();
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AddBin(TH2D *mlem, Int_t ic, Int_t jc, Int_t mode, Bool_t *used, TObjArray *pix)
{
- // Add a bin to the cluster
+/// Add a bin to the cluster
Int_t nx = mlem->GetNbinsX();
Int_t ny = mlem->GetNbinsY();
//_____________________________________________________________________________
TObject* AliMUONClusterFinderAZ::BinToPix(TH2D *mlem, Int_t jc, Int_t ic)
{
- // Translate histogram bin to pixel
+/// Translate histogram bin to pixel
Double_t yc = mlem->GetYaxis()->GetBinCenter(ic);
Double_t xc = mlem->GetXaxis()->GetBinCenter(jc);
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AddCluster(Int_t ic, Int_t nclust, TMatrixD *aijcluclu, Bool_t *used, Int_t *clustNumb, Int_t &nCoupled)
{
- // Add a cluster to the group of coupled clusters
+/// Add a cluster to the group of coupled clusters
for (Int_t i=0; i<nclust; i++) {
if (used[i]) continue;
//_____________________________________________________________________________
Double_t AliMUONClusterFinderAZ::MinGroupCoupl(Int_t nCoupled, Int_t *clustNumb, TMatrixD *aijcluclu, Int_t *minGroup)
{
- // Find group of clusters with minimum coupling to all the others
+/// Find group of clusters with minimum coupling to all the others
Int_t i123max = TMath::Min(3,nCoupled/2);
Int_t indx, indx1, indx2, indx3, nTot = 0;
//_____________________________________________________________________________
Int_t AliMUONClusterFinderAZ::SelectPad(Int_t nCoupled, Int_t nForFit, Int_t *clustNumb, Int_t *clustFit, TMatrixD *aijclupad)
{
- // Select pads for fit. If too many coupled clusters, find pads giving
- // the strongest coupling with the rest of clusters and exclude them from the fit.
+/// Select pads for fit. If too many coupled clusters, find pads giving
+/// the strongest coupling with the rest of clusters and exclude them from the fit.
Int_t npad = fnPads[0] + fnPads[1];
Double_t *padpix = 0;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Merge(Int_t nForFit, Int_t nCoupled, Int_t *clustNumb, Int_t *clustFit, TObjArray **clusters, TMatrixD *aijcluclu, TMatrixD *aijclupad)
{
- // Merge the group of clusters with the one having the strongest coupling with them
+/// Merge the group of clusters with the one having the strongest coupling with them
Int_t indx, indx1, npxclu, npxclu1, imax=0;
TObjArray *pix, *pix1;
//_____________________________________________________________________________
Int_t AliMUONClusterFinderAZ::Fit(Int_t iSimple, Int_t nfit, Int_t *clustFit, TObjArray **clusters, Double_t *parOk)
{
- // Find selected clusters to selected pad charges
+/// Find selected clusters to selected pad charges
TH2D *mlem = (TH2D*) gROOT->FindObject("mlem");
Double_t xmin = mlem->GetXaxis()->GetXmin() - mlem->GetXaxis()->GetBinWidth(1);
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Fcn1(Int_t & /*npar*/, Double_t * /*gin*/, Double_t &f, Double_t *par, Int_t /*iflag*/)
{
- // Fit for one track
- //AZ for Muinuit AliMUONClusterFinderAZ& c = *(AliMUONClusterFinderAZ::fgClusterFinder);
+/// Fit for one track
+/// AZ for Muinuit AliMUONClusterFinderAZ& c = *(AliMUONClusterFinderAZ::fgClusterFinder);
+
AliMUONClusterFinderAZ& c = *this; //AZ
Int_t cath, ix, iy, indx, npads=0;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::UpdatePads(Int_t /*nfit*/, Double_t *par)
{
- // Subtract the fitted charges from pads with strong coupling
+/// Subtract the fitted charges from pads with strong coupling
Int_t cath, ix, iy, indx;
Double_t charge, coef=0;
}
//_____________________________________________________________________________
-Bool_t AliMUONClusterFinderAZ::TestTrack(Int_t /*t*/) const {
-// Test if track was user selected
+Bool_t AliMUONClusterFinderAZ::TestTrack(Int_t /*t*/) const
+{
+/// Test if track was user selected
+
return kTRUE;
/*
if (fTrack[0]==-1 || fTrack[1]==-1) {
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AddRawCluster(Double_t x, Double_t y, Double_t qTot, Double_t fmin, Int_t nfit, Int_t *tracks, Double_t /*sigx*/, Double_t /*sigy*/, Double_t /*dist*/)
{
- //
- // Add a raw cluster copy to the list
- //
+/// Add a raw cluster copy to the list
+
if (qTot <= 0.501) return;
AliMUONRawCluster cnew;
if (fPadIJ[1][j] == -9) nover[cath]++;
if (fPadIJ[1][j] != 1 && fPadIJ[1][j] != -9) continue;
cnew.SetMultiplicity(cath,cnew.GetMultiplicity(cath)+1);
- if (fXyq[2][j] > cnew.GetPeakSignal(cath)) cnew.SetPeakSignal(cath,TMath::Nint (fXyq[2][j]));
+ if (fXyq[2][j] > cnew.GetPeakSignal(cath)) cnew.SetPeakSignal(cath,fXyq[2][j]);
//cnew.SetCharge(cath,cnew.GetCharge(cath) + TMath::Nint (fXyq[2][j]));
cnew.SetContrib(npads[cath],cath,fXyq[2][j]);
cnew.SetIndex(npads[cath],cath,TMath::Nint (fXyq[5][j]));
//_____________________________________________________________________________
Int_t AliMUONClusterFinderAZ::FindLocalMaxima(TObjArray *pixArray, Int_t *localMax, Double_t *maxVal)
{
- // Find local maxima in pixel space for large preclusters in order to
- // try to split them into smaller pieces (to speed up the MLEM procedure)
- // or to find additional fitting seeds if clusters were not completely resolved
+/// Find local maxima in pixel space for large preclusters in order to
+/// try to split them into smaller pieces (to speed up the MLEM procedure)
+/// or to find additional fitting seeds if clusters were not completely resolved
TH2D *hist = NULL;
//if (pixArray == fPixArray) hist = (TH2D*) gROOT->FindObject("anode");
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::FlagLocalMax(TH2D *hist, Int_t i, Int_t j, Int_t *isLocalMax)
{
- // Flag pixels (whether or not local maxima)
+/// Flag pixels (whether or not local maxima)
Int_t nx = hist->GetNbinsX();
Int_t ny = hist->GetNbinsY();
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::FindCluster(Int_t *localMax, Int_t iMax)
{
- // Find pixel cluster around local maximum #iMax and pick up pads
- // overlapping with it
+/// Find pixel cluster around local maximum \a iMax and pick up pads
+/// overlapping with it
TH2D *hist = (TH2D*) gROOT->FindObject("anode");
Int_t nx = hist->GetNbinsX();
delete [] used; used = 0;
}
-//_____________________________________________________________________________
-AliMUONClusterFinderAZ&
-AliMUONClusterFinderAZ::operator=(const AliMUONClusterFinderAZ& rhs)
-{
-// Protected assignement operator
-
- if (this == &rhs) return *this;
-
- AliFatal("Not implemented.");
-
- return *this;
-}
-
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::AddVirtualPad()
{
- // Add virtual pad (with small charge) to improve fit for some
- // clusters (when pad with max charge is at the extreme of the cluster)
+/// Add virtual pad (with small charge) to improve fit for some
+/// clusters (when pad with max charge is at the extreme of the cluster)
// Get number of pads in X and Y-directions
Int_t nInX = -1, nInY;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::PadsInXandY(Int_t &nInX, Int_t &nInY)
{
- // Find number of pads in X and Y-directions (excluding virtual ones and
- // overflows)
+/// Find number of pads in X and Y-directions (excluding virtual ones and
+/// overflows)
static Int_t nXsaved = 0, nYsaved = 0;
nXsaved = nYsaved = 0;
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Simple()
{
- // Process simple cluster (small number of pads) without EM-procedure
+/// Process simple cluster (small number of pads) without EM-procedure
Int_t nForFit = 1, clustFit[1] = {0}, nfit;
Double_t parOk[3] = {0.};
//_____________________________________________________________________________
void AliMUONClusterFinderAZ::Errors(AliMUONRawCluster *clus)
{
- // Correct reconstructed coordinates for some clusters and evaluate errors
+/// Correct reconstructed coordinates for some clusters and evaluate errors
Double_t qTot = clus->GetCharge(0), fmin = clus->GetChi2(0);
Double_t xreco = clus->GetX(0), yreco = clus->GetY(0), zreco = clus->GetZ(0);
Double_t dyc, Double_t /*dxc*/, Double_t qtot,
Double_t &yrec, Double_t &xrec, Double_t &erry, Double_t &errx)
{
- // Correct reconstructed coordinates for some clusters and evaluate errors
+/// Correct reconstructed coordinates for some clusters and evaluate errors
erry = 0.01;
errx = 0.144;
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff