* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-/*
-$Log$
-Revision 1.18 2001/01/26 21:37:53 morsch
-Use access functions to AliMUONDigit member data.
-
-Revision 1.17 2001/01/23 18:58:19 hristov
-Initialisation of some pointers
-
-Revision 1.16 2000/12/21 23:27:30 morsch
-Error in argument list of AddRawCluster corrected.
-
-Revision 1.15 2000/12/21 22:14:38 morsch
-Clean-up of coding rule violations.
-
-Revision 1.14 2000/10/23 16:03:45 morsch
-Correct z-position of all clusters created "on the flight".
-
-Revision 1.13 2000/10/23 13:38:23 morsch
-Set correct z-coordinate when cluster is split.
-
-Revision 1.12 2000/10/18 11:42:06 morsch
-- AliMUONRawCluster contains z-position.
-- Some clean-up of useless print statements during initialisations.
-
-Revision 1.11 2000/10/06 09:04:05 morsch
-- Dummy z-arguments in GetPadI, SetHit, FirstPad replaced by real z-coordinate
- to make code work with slat chambers (AM)
-- Replace GetPadI calls with unchecked x,y coordinates by pad iterator calls wherever possible.
-
-Revision 1.10 2000/10/03 13:51:57 egangler
-Removal of useless dependencies via forward declarations
-
-Revision 1.9 2000/10/02 16:58:29 egangler
-Cleaning of the code :
--> coding conventions
--> void Streamers
--> some useless includes removed or replaced by "class" statement
-
-Revision 1.8 2000/07/03 11:54:57 morsch
-AliMUONSegmentation and AliMUONHitMap have been replaced by AliSegmentation and AliHitMap in STEER
-The methods GetPadIxy and GetPadXxy of AliMUONSegmentation have changed name to GetPadI and GetPadC.
-
-Revision 1.7 2000/06/28 15:16:35 morsch
-(1) Client code adapted to new method signatures in AliMUONSegmentation (see comments there)
-to allow development of slat-muon chamber simulation and reconstruction code in the MUON
-framework. The changes should have no side effects (mostly dummy arguments).
-(2) Hit disintegration uses 3-dim hit coordinates to allow simulation
-of chambers with overlapping modules (MakePadHits, Disintegration).
-
-Revision 1.6 2000/06/28 12:19:18 morsch
-More consequent seperation of global input data services (AliMUONClusterInput singleton) and the
-cluster and hit reconstruction algorithms in AliMUONClusterFinderVS.
-AliMUONClusterFinderVS becomes the base class for clustering and hit reconstruction.
-It requires two cathode planes. Small modifications in the code will make it usable for
-one cathode plane and, hence, more general (for test beam data).
-AliMUONClusterFinder is now obsolete.
-
-Revision 1.5 2000/06/28 08:06:10 morsch
-Avoid global variables in AliMUONClusterFinderVS by seperating the input data for the fit from the
-algorithmic part of the class. Input data resides inside the AliMUONClusterInput singleton.
-It also naturally takes care of the TMinuit instance.
-
-Revision 1.4 2000/06/27 16:18:47 gosset
-Finally correct implementation of xm, ym, ixm, iym sizes
-when at least three local maxima on cathode 1 or on cathode 2
-
-Revision 1.3 2000/06/22 14:02:45 morsch
-Parameterised size of xm[], ym[], ixm[], iym[] correctly implemented (PH)
-Some HP scope problems corrected (PH)
-
-Revision 1.2 2000/06/15 07:58:48 morsch
-Code from MUON-dev joined
-
-Revision 1.1.2.3 2000/06/09 21:58:33 morsch
-Most coding rule violations corrected.
-
-Revision 1.1.2.2 2000/02/15 08:33:52 morsch
-Error in calculation of contribution map for double clusters (Split method) corrected (A.M.)
-Error in determination of track list for double cluster (FillCluster method) corrected (A.M.)
-Revised and extended SplitByLocalMaxima method (Isabelle Chevrot):
- - For clusters with more than 2 maxima on one of the cathode planes all valid
- combinations of maxima on the two cathodes are preserved. The position of the maxima is
- taken as the hit position.
- - New FillCluster method with 2 arguments to find tracks associated to the clusters
- defined above added. (Method destinction by argument list not very elegant in this case,
- should be revides (A.M.)
- - Bug in if-statement to handle maximum 1 maximum per plane corrected
- - Two cluster per cathode but only 1 combination valid is handled.
- - More rigerous treatment of 1-2 and 2-1 combinations of maxima.
-*/
+/* $Id$ */
#include "AliMUONClusterFinderVS.h"
#include "AliMUONDigit.h"
#include <TF1.h>
#include <stdio.h>
-#include <iostream.h>
+#include <Riostream.h>
//_____________________________________________________________________
// This function is minimized in the double-Mathieson fit
}
}
-AliMUONClusterFinderVS::AliMUONClusterFinderVS(
- const AliMUONClusterFinderVS & clusterFinder)
+AliMUONClusterFinderVS::AliMUONClusterFinderVS(const AliMUONClusterFinderVS & clusterFinder):TObject(clusterFinder)
{
// Dummy copy Constructor
;
} else if (fNLocal[0]==2 && fNLocal[1]==2) {
//
// Let's look for ghosts first
-//
+
Float_t xm[4][2], ym[4][2];
Float_t dpx, dpy, dx, dy;
Int_t ixm[4][2], iym[4][2];
// Analyse the combinations and keep those that are possible !
// For each combination check consistency in x and y
- Int_t iacc;
- Bool_t accepted[4];
+ Int_t iacc;
+ Bool_t accepted[4];
+ Float_t dr[4] = {1.e4, 1.e4, 1.e4, 1.e4};
iacc=0;
-
+
+// In case of staggering maxima are displaced by exactly half the pad-size in y.
+// We have to take into account the numerical precision in the consistency check;
+ Float_t eps = 1.e-5;
+//
for (ico=0; ico<4; ico++) {
accepted[ico]=kFALSE;
// cathode one: x-coordinate
isec=fSeg[1]->Sector(ixm[ico][1], iym[ico][1]);
dpy=fSeg[1]->Dpy(isec)/2.;
dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
- if (fDebugLevel)
- printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
- if ((dx <= dpx) && (dy <= dpy)) {
+ if (fDebugLevel>1)
+ printf("\n %i %f %f %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy, dx, dpx );
+ if ((dx <= dpx) && (dy <= dpy+eps)) {
// consistent
accepted[ico]=kTRUE;
+ dr[ico] = TMath::Sqrt(dx*dx+dy*dy);
iacc++;
} else {
// reject
accepted[ico]=kFALSE;
}
}
+ printf("\n iacc= %d:\n", iacc);
+ if (iacc == 3) {
+ if (accepted[0] && accepted[1]) {
+ if (dr[0] >= dr[1]) {
+ accepted[0]=kFALSE;
+ } else {
+ accepted[1]=kFALSE;
+ }
+ }
+ if (accepted[2] && accepted[3]) {
+ if (dr[2] >= dr[3]) {
+ accepted[2]=kFALSE;
+ } else {
+ accepted[3]=kFALSE;
+ }
+ }
+/*
+// eliminate one candidate
+ Float_t drmax = 0;
+ Int_t icobad = -1;
+
+ for (ico=0; ico<4; ico++) {
+ if (accepted[ico] && dr[ico] > drmax) {
+ icobad = ico;
+ drmax = dr[ico];
+ }
+ }
+
+ accepted[icobad] = kFALSE;
+*/
+ iacc = 2;
+ }
+
+
+ printf("\n iacc= %d:\n", iacc);
if (fDebugLevel) {
if (iacc==2) {
fprintf(stderr,"\n iacc=2: No problem ! \n");
if ((accepted[0]&&accepted[1]) || (accepted[2]&&accepted[3])) {
fprintf(stderr,"\n Maximum taken twice !!!\n");
-// Have a try !! with that
+// Have a try !! with that
if (accepted[0]&&accepted[3]) {
fXInit[0]=xm[0][1];
fYInit[0]=ym[0][0];
Int_t iacc;
Bool_t accepted[4];
iacc=0;
+ // In case of staggering maxima are displaced by exactly half the pad-size in y.
+ // We have to take into account the numerical precision in the consistency check;
+ Float_t eps = 1.e-5;
+
for (ico=0; ico<2; ico++) {
accepted[ico]=kFALSE;
isec=fSeg[0]->Sector(ixm[ico][0], iym[ico][0]);
isec=fSeg[1]->Sector(ixm[ico][1], iym[ico][1]);
dpy=fSeg[1]->Dpy(isec)/2.;
dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
- if (fDebugLevel)
+ if (fDebugLevel>1)
printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
- if ((dx <= dpx) && (dy <= dpy)) {
+ if ((dx <= dpx) && (dy <= dpy+eps)) {
// consistent
accepted[ico]=kTRUE;
iacc++;
Float_t chi21 = 100;
Float_t chi22 = 100;
+ Float_t chi23 = 100;
+
+ // Initial value for charge ratios
+ fQrInit[0]=Float_t(fQ[fIndLocal[0][0]][0])/
+ Float_t(fQ[fIndLocal[0][0]][0]+fQ[fIndLocal[1][0]][0]);
+ fQrInit[1]=fQrInit[0];
- if (accepted[0]) {
+ if (accepted[0] && accepted[1]) {
+
+ fXInit[0]=0.5*(xm[0][1]+xm[0][0]);
+ fYInit[0]=ym[0][0];
+ fXInit[1]=0.5*(xm[0][1]+xm[1][0]);
+ fYInit[1]=ym[1][0];
+ fQrInit[0]=0.5;
+ fQrInit[1]=0.5;
+ chi23=CombiDoubleMathiesonFit(c);
+ if (chi23<10) {
+ Split(c);
+ Float_t yst;
+ yst = fYFit[0];
+ fYFit[0] = fYFit[1];
+ fYFit[1] = yst;
+ Split(c);
+ }
+ } else if (accepted[0]) {
fXInit[0]=xm[0][1];
fYInit[0]=ym[0][0];
fXInit[1]=xm[1][0];
if (chi22<10) Split(c);
}
- if (chi21 > 10 && chi22 > 10) {
+ if (chi21 > 10 && chi22 > 10 && chi23 > 10) {
// We keep only the combination found (X->cathode 2, Y->cathode 1)
for (Int_t ico=0; ico<2; ico++) {
if (accepted[ico]) {
Int_t iacc;
Bool_t accepted[4];
iacc=0;
+ // In case of staggering maxima are displaced by exactly half the pad-size in y.
+ // We have to take into account the numerical precision in the consistency check;
+ Float_t eps = 1.e-5;
+
for (ico=0; ico<2; ico++) {
accepted[ico]=kFALSE;
isec=fSeg[1]->Sector(ixm[ico][1], iym[ico][1]);
dpy=fSeg[1]->Dpy(isec)/2.;
dy=TMath::Abs(ym[ico][0]-ym[ico][1]);
- if (fDebugLevel)
+ if (fDebugLevel>0)
printf("\n %i %f %f %f %f \n", ico, ym[ico][0], ym[ico][1], dy, dpy );
- if ((dx <= dpx) && (dy <= dpy)) {
+ if ((dx <= dpx) && (dy <= dpy+eps)) {
// consistent
accepted[ico]=kTRUE;
fprintf(stderr,"ico %d\n",ico);
Float_t chi21 = 100;
Float_t chi22 = 100;
+ Float_t chi23 = 100;
- if (accepted[0]) {
+ fQrInit[1]=Float_t(fQ[fIndLocal[0][1]][1])/
+ Float_t(fQ[fIndLocal[0][1]][1]+fQ[fIndLocal[1][1]][1]);
+
+ fQrInit[0]=fQrInit[1];
+
+
+ if (accepted[0] && accepted[1]) {
+ fXInit[0]=xm[0][1];
+ fYInit[0]=0.5*(ym[0][0]+ym[0][1]);
+ fXInit[1]=xm[1][1];
+ fYInit[1]=0.5*(ym[0][0]+ym[1][1]);
+ fQrInit[0]=0.5;
+ fQrInit[1]=0.5;
+ chi23=CombiDoubleMathiesonFit(c);
+ if (chi23<10) {
+ Split(c);
+ Float_t yst;
+ yst = fYFit[0];
+ fYFit[0] = fYFit[1];
+ fYFit[1] = yst;
+ Split(c);
+ }
+ } else if (accepted[0]) {
fXInit[0]=xm[0][0];
fYInit[0]=ym[0][1];
fXInit[1]=xm[1][1];
if (chi22<10) Split(c);
}
- if (chi21 > 10 && chi22 > 10) {
+ if (chi21 > 10 && chi22 > 10 && chi23 > 10) {
//We keep only the combination found (X->cathode 2, Y->cathode 1)
for (Int_t ico=0; ico<2; ico++) {
if (accepted[ico]) {
}
}
-void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* c)
+void AliMUONClusterFinderVS::FindLocalMaxima(AliMUONRawCluster* /*c*/)
{
// Find all local maxima of a cluster
if (fDebugLevel)
// Two local maxima on cathode 1 and one maximum on cathode 2
// Look for local maxima considering left and right neighbours on the 2nd cathode only
cath=1;
- Int_t cath1=0;
+ Int_t cath1 = 0;
+ Float_t eps = 1.e-5;
+
//
// Loop over cluster digits
for (i=0; i<fMul[cath]; i++) {
dpy=fSeg[cath]->Dpy(isec);
if (isLocal[i][cath]) continue;
// Pad position should be consistent with position of local maxima on the opposite cathode
- if ((TMath::Abs(fY[i][cath]-fY[fIndLocal[0][cath1]][cath1]) > dpy/2.) &&
- (TMath::Abs(fY[i][cath]-fY[fIndLocal[1][cath1]][cath1]) > dpy/2.))
+ if ((TMath::Abs(fY[i][cath]-fY[fIndLocal[0][cath1]][cath1]) > dpy/2.+eps) &&
+ (TMath::Abs(fY[i][cath]-fY[fIndLocal[1][cath1]][cath1]) > dpy/2.+eps))
continue;
+
//
// get neighbours for that digit and assume that it is local maximum
isLocal[i][cath]=kTRUE;
// iNN counts the number of neighbours with signal, it should be 1 or 2
Int_t iNN=0;
for (fSeg[cath]
- ->FirstPad(fX[i][cath], fY[i][cath], fZPlane, 0., dpx);
+ ->FirstPad(fX[i][cath], fY[i][cath], fZPlane, dpx, 0.);
fSeg[cath]
->MorePads();
fSeg[cath]
->NextPad())
{
+
ix = fSeg[cath]->Ix();
iy = fSeg[cath]->Iy();
-
+
// skip the current pad
if (ix == fIx[i][cath]) continue;
} else c->fPhysicsMap[i]=1;
//
//
- if (fDebugLevel)
+ if (fDebugLevel>1)
fprintf(stderr,"q %d c->fPeakSignal[cath] %d\n",q,c->fPeakSignal[cath]);
// peak signal and track list
if (q>c->fPeakSignal[cath]) {
iy=yList[in];
if (fHitMap[cath]->TestHit(ix,iy)==kUnused) {
- if (fDebugLevel)
+ if (fDebugLevel>1)
printf("\n Neighbours %d %d %d", cath, ix, iy);
FindCluster(ix, iy, cath, c);
}
{
ix = fSeg[iop]->Ix(); iy = fSeg[iop]->Iy();
- if (fDebugLevel)
+ if (fDebugLevel > 1)
printf("\n ix, iy: %f %f %f %d %d %d", x,y,z,ix, iy, fSector);
if (fHitMap[iop]->TestHit(ix,iy)==kUnused){
iXopp[nOpp]=ix;
iYopp[nOpp++]=iy;
- if (fDebugLevel)
+ if (fDebugLevel > 1)
printf("\n Opposite %d %d %d", iop, ix, iy);
}
fSector= fSeg[cath]->Sector(i,j)/100;
if (fDebugLevel)
printf("\n New Seed %d %d ", i,j);
-
+
+
FindCluster(i,j,cath,c);
// ^^^^^^^^^^^^^^^^^^^^^^^^
// center of gravity
- c.fX[0] /= c.fQ[0];
+ if (c.fX[0]!=0.) c.fX[0] /= c.fQ[0];
// Force on anod
c.fX[0]=fSeg[0]->GetAnod(c.fX[0]);
- c.fY[0] /= c.fQ[0];
- c.fX[1] /= c.fQ[1];
-// Force on anod
+ if (c.fY[0]!=0.) c.fY[0] /= c.fQ[0];
+
+ if(c.fQ[1]!=0.) c.fX[1] /= c.fQ[1];
+
+ // Force on anod
c.fX[1]=fSeg[0]->GetAnod(c.fX[1]);
- c.fY[1] /= c.fQ[1];
+ if(c.fQ[1]!=0.) c.fY[1] /= c.fQ[1];
c.fZ[0] = fZPlane;
c.fZ[1] = fZPlane;
return fmin;
}
-Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster *c)
+Float_t AliMUONClusterFinderVS::CombiSingleMathiesonFit(AliMUONRawCluster * /*c*/)
{
// Perform combined Mathieson fit on both cathode planes
//
return fmin;
}
-Bool_t AliMUONClusterFinderVS::DoubleMathiesonFit(AliMUONRawCluster *c, Int_t cath)
+Bool_t AliMUONClusterFinderVS::DoubleMathiesonFit(AliMUONRawCluster * /*c*/, Int_t cath)
{
// Performs a double Mathieson fit on one cathode
//
return kTRUE;
}
-Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster *c)
+Float_t AliMUONClusterFinderVS::CombiDoubleMathiesonFit(AliMUONRawCluster * /*c*/)
{
//
// Perform combined double Mathieson fit on both cathode planes
fSeg[1]->MorePads(); fSeg[1]->NextPad())
{
ix=fSeg[1]->Ix(); iy=fSeg[1]->Iy();
+// if (fHitMap[1]->TestHit(ix, iy) == kEmpty) continue;
fSeg[1]->GetPadC(ix,iy,upper[0],ydum,zdum);
if (icount ==0) lower[0]=upper[0];
icount++;
}
if (lower[0]>upper[0]) {xdum=lower[0]; lower[0]=upper[0]; upper[0]=xdum;}
+// vstart[0] = 0.5*(lower[0]+upper[0]);
+
+
icount=0;
for (fSeg[0]->FirstPad(fXInit[0], fYInit[0], fZPlane, 0., dpy);
fSeg[0]->MorePads(); fSeg[0]->NextPad())
{
ix=fSeg[0]->Ix(); iy=fSeg[0]->Iy();
+// if (fHitMap[0]->TestHit(ix, iy) == kEmpty) continue;
fSeg[0]->GetPadC(ix,iy,xdum,upper[1],zdum);
if (icount ==0) lower[1]=upper[1];
icount++;
}
+
if (lower[1]>upper[1]) {xdum=lower[1]; lower[1]=upper[1]; upper[1]=xdum;}
+// vstart[1] = 0.5*(lower[1]+upper[1]);
+
fSeg[1]->GetPadI(fXInit[1], fYInit[1], fZPlane, ix, iy);
isec=fSeg[1]->Sector(ix, iy);
fSeg[1]->MorePads(); fSeg[1]->NextPad())
{
ix=fSeg[1]->Ix(); iy=fSeg[1]->Iy();
+// if (fHitMap[1]->TestHit(ix, iy) == kEmpty) continue;
fSeg[1]->GetPadC(ix,iy,upper[2],ydum,zdum);
if (icount ==0) lower[2]=upper[2];
icount++;
}
if (lower[2]>upper[2]) {xdum=lower[2]; lower[2]=upper[2]; upper[2]=xdum;}
+ // vstart[2] = 0.5*(lower[2]+upper[2]);
icount=0;
fSeg[0]-> MorePads(); fSeg[0]->NextPad())
{
ix=fSeg[0]->Ix(); iy=fSeg[0]->Iy();
+// if (fHitMap[0]->TestHit(ix, iy) != kEmpty) continue;
+
fSeg[0]->GetPadC(ix,iy,xdum,upper[3],zdum);
if (icount ==0) lower[3]=upper[3];
icount++;
+
}
if (lower[3]>upper[3]) {xdum=lower[3]; lower[3]=upper[3]; upper[3]=xdum;}
-
+
+// vstart[3] = 0.5*(lower[3]+upper[3]);
+
lower[4]=0.;
upper[4]=1.;
lower[5]=0.;
//
// Minimisation functions
// Single Mathieson
-void fcnS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
+void fcnS1(Int_t & /*npar*/, Double_t * /*gin*/, Double_t &f, Double_t *par, Int_t /*iflag*/)
{
AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
Int_t i;
f=chisq;
}
-void fcnCombiS1(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
+void fcnCombiS1(Int_t & /*npar*/, Double_t * /*gin*/, Double_t &f, Double_t *par, Int_t /*iflag*/)
{
AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
Int_t i, cath;
}
// Double Mathieson
-void fcnS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
+void fcnS2(Int_t & /*npar*/, Double_t * /*gin*/, Double_t &f, Double_t *par, Int_t /*iflag*/)
{
AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
Int_t i;
}
// Double Mathieson
-void fcnCombiS2(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
+void fcnCombiS2(Int_t & /*npar*/, Double_t * /*gin*/, Double_t &f, Double_t *par, Int_t /*iflag*/)
{
AliMUONClusterInput& clusterInput = *(AliMUONClusterInput::Instance());
Int_t i, cath;
// Add a raw cluster copy to the list
//
AliMUON *pMUON=(AliMUON*)gAlice->GetModule("MUON");
- pMUON->AddRawCluster(fInput->Chamber(),c);
+ pMUON->GetMUONData()->AddRawCluster(fInput->Chamber(),c);
fNRawClusters++;
// if (fDebugLevel)
fprintf(stderr,"\nfNRawClusters %d\n",fNRawClusters);
}
AliMUONClusterFinderVS& AliMUONClusterFinderVS
-::operator = (const AliMUONClusterFinderVS& rhs)
+::operator = (const AliMUONClusterFinderVS& /*rhs*/)
{
// Dummy assignment operator
return *this;
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff