//-----------------------------------------------------//
// //
-// Source File : PMDClustering.cxx, Version 00 //
+// Source File : PMDClusteringV1.cxx, Version 00 //
// //
// Date : September 26 2002 //
// //
It is assumed that for data we have 0 <= id <= 48
and 0 <= jd <=96
*/
+
int i, i1, i2, j, nmx1, incr, id, jd;
Int_t celldataX[15], celldataY[15];
Float_t clusdata[7];
const float ktwobysqrt3 = 1.1547; // 2./sqrt(3.)
- for (id = 0; id < kNDIMXr; id++)
+ // ndimXr and ndimYr are different because of different module size
+
+ Int_t ndimXr =0;
+ Int_t ndimYr =0;
+
+ if (ismn < 12)
+ {
+ ndimXr = 96;
+ ndimYr = 48;
+ }
+ else if (ismn >= 12 && ismn <= 23)
+ {
+ ndimXr = 48;
+ ndimYr = 96;
+ }
+
+ for (Int_t i =0; i < kNDIMX; i++)
{
- for (jd = 0; jd < kNDIMYr; jd++)
+ for (Int_t j =0; j < kNDIMY; j++)
+ {
+ fEdepCell[i][j] = 0;
+ fCellTrNo[i][j] = -1;
+ }
+ }
+
+ for (id = 0; id < ndimXr; id++)
+ {
+ for (jd = 0; jd < ndimYr; jd++)
{
j=jd;
- i=id+(kNDIMYr/2-1)-(jd/2);
- fEdepCell[i][j] = celladc[id][jd];
+ i=id+(ndimYr/2-1)-(jd/2);
+
+ if (ismn < 12)
+ {
+ fEdepCell[i][j] = celladc[jd][id];
+ fCellTrNo[i][j] = jd*10000+id; /* for association */
+ }
+ else if (ismn >= 12 && ismn <= 23)
+ {
+ fEdepCell[i][j] = celladc[id][jd];
+ fCellTrNo[i][j] = id*10000+jd; /* for association */
+ }
+
}
}
Order(); // order the data
if (fEdepCell[i1][i2] > 0.) {ave = ave + fEdepCell[i1][i2];}
if (fEdepCell[i1][i2] > cutoff ) nmx1 = nmx1 + 1;
}
- // nmx1 --- number of cells having ener dep >= cutoff
AliDebug(1,Form("Number of cells having energy >= %f are %d",cutoff,nmx1));
- // if (nmx1 == 0 | nmx1 == -1) return;
-
if (nmx1 == 0) nmx1 = 1;
ave=ave/nmx1;
RefClust(incr);
AliDebug(1,Form("Detector Plane = %d Serial Module No = %d Number of clusters = %d",idet, ismn, fClno));
+
for(i1=0; i1<=fClno; i1++)
{
Float_t cluRAD = (Float_t) fClusters[4][i1];
Float_t cluY0 = ktwobysqrt3*cluYC;
Float_t cluX0 = cluXC - cluY0/2.;
+
+
//
// Cluster X centroid is back transformed
//
- clusdata[0] = cluX0 - (48-1) + cluY0/2.;
- clusdata[1] = cluY0;
- clusdata[2] = cluADC;
- clusdata[3] = cluCELLS;
- clusdata[4] = cluRAD;
- clusdata[5] = 0.;
+ if (ismn < 12)
+ {
+ clusdata[0] = cluX0 - (24-1) + cluY0/2.;
+ }
+ else if (ismn >= 12 && ismn <= 23)
+ {
+ clusdata[0] = cluX0 - (48-1) + cluY0/2.;
+ }
+
+ clusdata[1] = cluY0;
+ clusdata[2] = cluADC;
+ clusdata[3] = cluCELLS;
+ clusdata[4] = cluRAD;
+ clusdata[5] = 0.;
//
// Cells associated with a cluster
//
for (Int_t ihit = 0; ihit < 15; ihit++)
{
- celldataX[ihit] = 2;
- celldataY[ihit] = 5;
+
+ if (ismn < 12)
+ {
+ celldataX[ihit] = fClTr[ihit][i1]%10000;
+ celldataY[ihit] = fClTr[ihit][i1]/10000;
+ }
+ else if (ismn >= 12 && ismn <= 23)
+ {
+ celldataX[ihit] = fClTr[ihit][i1]/10000;
+ celldataY[ihit] = fClTr[ihit][i1]%10000;
+ }
}
pmdcl = new AliPMDcluster(idet, ismn, clusdata, celldataX, celldataY);
int ig, nsupcl, lev1[20], lev2[20];
double x[4500], y[4500], z[4500], x1, y1, z1, x2, y2, z2, dist;
double xc[4500], yc[4500], zc[4500], cells[4500], sum, rc[4500], rr;
+
+
+ //asso
+ Int_t t[4500],cellCount[4500];
+ for(i=0; i<4500; i++)
+ {
+ t[i]=-1;
+ cellCount[i]=0;
+ }
+
+
// fClno counts the final clusters
// nsupcl = # of superclusters; ncl[i]= # of cells in supercluster i
// x, y and z store (x,y) coordinates of and energy deposited in a cell
fClusters[3][fClno] = 1.;
fClusters[4][fClno] = 0.5;
- //ofl1 << icl << " " << fCoord[0][i1][i2] << " " << fCoord[1][i1][i2] <<
+
+ //asso
+
+ fClTr[0][fClno]=fCellTrNo[i1][i2];
+ for(Int_t icltr=1;icltr<14;icltr++)
+ {
+ fClTr[icltr][fClno]=-1;
+ }
+
+ //ofl1 << icl << " " << fCoord[0][i1][i2] << " " << fCoord[1][i1][i2] <<
//" " << fEdepCell[i1][i2] << " " << fClusters[3][fClno] <<endl;
}else if(ncl[i] == 1){
// two cell super-cluster --> single cluster
y1 = fCoord[1][i1][i2];
z1 = fEdepCell[i1][i2];
+ //asso
+ fClTr[0][fClno]=fCellTrNo[i1][i2];
+ //
id = id+1;
i1 = fInfcl[1][id];
y2 = fCoord[1][i1][i2];
z2 = fEdepCell[i1][i2];
+ //asso
+
+ fClTr[1][fClno]=fCellTrNo[i1][i2];
+ for(Int_t icltr=2;icltr<14;icltr++)
+ {
+ fClTr[icltr][fClno] = -1;
+ }
+ //
+
fClusters[0][fClno] = (x1*z1+x2*z2)/(z1+z2);
fClusters[1][fClno] = (y1*z1+y2*z2)/(z1+z2);
fClusters[2][fClno] = z1+z2;
// << " " << fClusters[2][fClno] << " " <<fClusters[3][fClno] <<endl;
}
else{
+
+ //asso
+ for(Int_t icg=0;icg<4500;icg++)
+ {
+ cellCount[icg]=0;
+ }
+ //
+
id = id + 1;
iord[0] = 0;
// super-cluster of more than two cells - broken up into smaller
x[0] = fCoord[0][i1][i2];
y[0] = fCoord[1][i1][i2];
z[0] = fEdepCell[i1][i2];
+
+ //asso
+ t[0]=fCellTrNo[i1][i2];
+ //
+
iord[0] = 0;
for(j=1;j<=ncl[i];j++){
x[j] = fCoord[0][i1][i2];
y[j] = fCoord[1][i1][i2];
z[j] = fEdepCell[i1][i2];
+
+ //asso
+ t[j]=fCellTrNo[i1][i2];
+ //
+
+
}
// arranging cells within supercluster in decreasing order
for(j=1;j<=ncl[i];j++){
for(k=0; k<=ig; k++){
dist=Distance(x[j], y[j], xc[k], yc[k]);
if(dist < sqrt(3.) ){
+
+ //asso
+ fClTr[cellCount[k]][fClno+k+1]=t[j];
+ cellCount[k]++;
+ //
+
lev1[0]++;
i1=lev1[0];
lev1[i1]=k;
}
}
}
+
+ // zero rest of the cell array
+ //asso
+ for(k=0; k<=ig; k++)
+ {
+ for(Int_t icltr=cellCount[k];icltr<14;icltr++)
+ {
+ fClTr[icltr][fClno]=-1;
+ }
+ }
+ //
+
+
+
for(j=0; j<=ig; j++){
if (fClno >= 5000) {
AliWarning("RefClust: Too many clusters! more than 5000");
* See cxx source for full Copyright notice */
//-----------------------------------------------------//
// //
-// Header File : PMDClustering.h, Version 00 //
+// Header File : PMDClusteringV1.h, Version 00 //
// //
// Date : September 26 2002 //
// //
protected:
static const Double_t fgkSqroot3by2; // fgkSqroot3by2 = sqrt(3.)/2.
- /*enum {
- kNMX = 4608,
- kNDIMX = 48,
- kNDIMY = 96
- };*/
- /*
- Proposed changes inNMX, kNDIMX and kNDIMY by S. C. Phatak to account
- for rectangular ( vs rhomboid ) geometry.
- To keep the clustering functional, we define a rhomboid which
- superscribes the rectangle. So we need to pad up dummy cells in x
- direction. The number of these cells is 96/2-1=47 in each row ( value
- of x ). For first two rows, all dummy cells are to the left. For
- every two rows add one cell to right and subtract one from left.
- So previous (i,j) values go over to ( i',j) i'=i+(96-j)/2-1
- Note we use C++ convention so i and j run from 0 to 47 or 95.
- */
-
+
enum {
- kNMX = 9120,
- kNDIMX = 95,
- kNDIMY = 96,
- kNDIMXr = 48,
- kNDIMYr = 96
+ kNMX = 11424,
+ kNDIMX = 119,
+ kNDIMY = 96
};
+
/*
kNMX : # of cells in a supermodule
kNDIMX : maximum number of cells along x direction (origin at one corner)
kNDIMY : maximum number of cells along axis at 60 degrees with x axis
*/
- Double_t fEdepCell[kNDIMX][kNDIMY]; //energy(ADC) in each cell of the supermodule
- Double_t fClusters[5][5000]; // Cluster informations
- Int_t fClno; // number of clusters in a supermodule
+ Double_t fEdepCell[kNDIMX][kNDIMY]; //energy(ADC) in each cell
+ Double_t fClusters[5][5000]; // Cluster informations
+ Int_t fClno; // number of clusters in a supermodule
/*
clusters[0][i] --- x position of the cluster center
clusters[4][i] --- cluster radius
*/
- Int_t fIord[2][kNMX]; // ordered list of i and j according to decreasing energy dep.
+ //Variables for association
+ Int_t fCellTrNo[kNDIMX][kNDIMY]; //id x-y value of cells
+ Int_t fClTr[15][5000]; // 1d x-y cell info of attached cells
+
+ Int_t fIord[2][kNMX]; // ordered list of i and j according
+ // to decreasing energy dep.
Int_t fInfocl[2][kNDIMX][kNDIMY]; // cellwise information on the cluster to which the cell
- Int_t fInfcl[3][kNMX]; // cluster information [0][i] -- cluster number
+ Int_t fInfcl[3][kNMX]; // cluster information [0][i]
+ // -- cluster number
Double_t fCoord[2][kNDIMX][kNDIMY];
/*
Float_t fCutoff; // Energy(ADC) cutoff per cell before clustering
- ClassDef(AliPMDClusteringV1,1) // Does clustering for PMD
+ ClassDef(AliPMDClusteringV1,2) // Does clustering for PMD
};
#endif
branch1->SetAddress(&fRecpoints);
/**********************************************************************
* det : Detector, 0: PRE & 1:CPV *
- * smn : Serial Module Number from which Super Module Number *
- * and Unit Module Numbers are extracted *
+ * smn : Serial Module Number from 0 to 23 for both detector *
* xpos : x-position of the cluster *
* ypos : y-position of the cluster *
* THESE xpos & ypos are not the true xpos and ypos *
* adc : ADC contained in the cluster *
* ncell : Number of cells contained in the cluster *
* rad : radius of the cluster (1d fit) *
- * ism : Supermodule number extracted from smn *
- * ium : Unit module number extracted from smn *
* xpad : TRUE x-position of the cluster *
* ypad : TRUE y-position of the cluster *
**********************************************************************/
- Int_t ism, ium;
+
Int_t det,smn;
Float_t xpos,ypos, xpad, ypad;
Float_t adc, ncell, sigx, sigy;
ncell = pmdrecpoint->GetClusCells();
sigx = pmdrecpoint->GetClusSigmaX();
sigy = pmdrecpoint->GetClusSigmaY();
- //
- // Now change the xpos and ypos to its original values
- // for the unit modules which are earlier changed.
- // xpad and ypad are the real positions.
- //
- if(det == 0 || det == 1)
- {
- if(smn < 12)
- {
- ism = smn/6;
- ium = smn - ism*6;
- xpad = ypos;
- ypad = xpos;
- }
- else if( smn >= 12 && smn < 24)
- {
- ism = smn/6;
- ium = smn - ism*6;
- xpad = xpos;
- ypad = ypos;
- }
- }
+
//
// User has to plug in his analysis code here
//
- fprintf(fpw,"%d %d %d %d %f %f %f %f %f %f\n",
- det,smn,ism,ium,xpad,ypad,adc,ncell,sigx,sigy);
+ fprintf(fpw,"%d %d %d %d\n",
+ det,smn,xpos,ypos);
//
// Plot the cluster centroid to see the PMD geometry
// using the PMD Utility class
//
- if (det == 1)
+ if (det == 0)
{
// Draw only for PRE plane
- cc->RectGeomCellPos(ism,ium,xpad,ypad,xx,yy);
+ //cc->RectGeomCellPos(ism,xpad,ypad,xx,yy);
+ cc->RectGeomCellPos(smn,xpos,ypos,xx,yy);
h2->Fill(xx,yy);
}
// Default destructor
}
-void AliPMDUtility::RectGeomCellPos(Int_t ism, Int_t ium, Int_t xpad, Int_t ypad, Float_t &xpos, Float_t &ypos)
+void AliPMDUtility::RectGeomCellPos(Int_t ism, Int_t xpad, Int_t ypad, Float_t &xpos, Float_t &ypos)
{
// This routine finds the cell eta,phi for the new PMD rectangular
// geometry in ALICE
// Authors : Bedanga Mohanty and Dipak Mishra - 29.4.2003
- // modified by B. K. Nnadi for change of coordinate sys
+ // modified by B. K. Nandi for change of coordinate sys
//
// SMA ---> Supermodule Type A ( SM - 0)
// SMAR ---> Supermodule Type A ROTATED ( SM - 1)
// SMB ---> Supermodule Type B ( SM - 2)
// SMBR ---> Supermodule Type B ROTATED ( SM - 3)
//
- // ism : number of supermodules in one plane = 4
- // ium : number of unitmodules in one SM = 6
- // gbum : (global) unit module numbering in a supermodule
- //
-
- Int_t gbum = ism*6 + ium;
- Int_t irow = xpad;
- Int_t icol = ypad;
+ // ism : Serial module number from 0 to 23 for each plane
+
// Corner positions (x,y) of the 24 unit moudles in ALICE PMD
-
-
double xcorner[24] =
{
74.8833, 53.0045, 31.1255, //Type-A
//
Float_t cellRadius = 0.25;
Float_t shift = 0.0;
- if(irow%2 == 0)
+ if(xpad%2 == 0)
{
shift = -cellRadius/2.0;
}
shift = 0.0;
}
- if(ism == 0)
+
+ if(ism < 6)
{
- ypos = ycorner[gbum] - irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] - icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] - (Float_t) xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] - (Float_t) ypad*kSqroot3*kCellRadius;
}
- else if(ism == 1)
+ else if(ism >=6 && ism < 12)
{
- ypos = ycorner[gbum] + irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] + icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] + (Float_t) xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] + (Float_t) ypad*kSqroot3*kCellRadius;
}
- else if(ism == 2)
+ else if(ism >= 12 && ism < 18)
{
- ypos = ycorner[gbum] - irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] - icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] - (Float_t) xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] - (Float_t) ypad*kSqroot3*kCellRadius;
}
- else if(ism == 3)
+ else if(ism >= 18 && ism < 24)
{
- ypos = ycorner[gbum] + irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] + icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] + (Float_t) xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] + (Float_t) ypad*kSqroot3*kCellRadius;
}
}
-void AliPMDUtility::RectGeomCellPos(Int_t ism, Int_t ium, Float_t xpad, Float_t ypad, Float_t &xpos, Float_t &ypos)
+void AliPMDUtility::RectGeomCellPos(Int_t ism, Float_t xpad, Float_t ypad, Float_t &xpos, Float_t &ypos)
{
// If the xpad and ypad inputs are float, then 0.5 is added to it
// to find the layer which is shifted.
// SMB ---> Supermodule Type B ( SM - 2)
// SMBR ---> Supermodule Type B ROTATED ( SM - 3)
//
- // ism : number of supermodules in one plane = 4
- // ium : number of unitmodules in one SM = 6
- // gbum : (global) unit module numbering in a supermodule
- //
-
- Int_t gbum = ism*6 + ium;
- Float_t irow = xpad;
- Float_t icol = ypad;
+ // ism : Serial Module number from 0 to 23 for each plane
// Corner positions (x,y) of the 24 unit moudles in ALICE PMD
-
double xcorner[24] =
{
74.8833, 53.0045, 31.1255, //Type-A
//
Float_t cellRadius = 0.25;
Float_t shift = 0.0;
- Int_t iirow = (Int_t) (irow+0.5);
+ Int_t iirow = (Int_t) (xpad+0.5);
if(iirow%2 == 0)
{
shift = -cellRadius/2.0;
shift = 0.0;
}
-
- if(ism == 0)
+ if(ism < 6)
{
- ypos = ycorner[gbum] - irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] - icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] - xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] - ypad*kSqroot3*kCellRadius;
}
- else if(ism == 1)
+ else if(ism >=6 && ism < 12)
{
- ypos = ycorner[gbum] + irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] + icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] + xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] + ypad*kSqroot3*kCellRadius;
}
- else if(ism == 2)
+ else if(ism >= 12 && ism < 18)
{
- ypos = ycorner[gbum] - irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] - icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] - xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] - ypad*kSqroot3*kCellRadius;
}
- else if(ism == 3)
+ else if(ism >= 18 && ism < 24)
{
- ypos = ycorner[gbum] + irow*kCellRadius*2.0 + shift;
- xpos = xcorner[gbum] + icol*kSqroot3*kCellRadius;
+ ypos = ycorner[ism] + xpad*kCellRadius*2.0 + shift;
+ xpos = xcorner[ism] + ypad*kSqroot3*kCellRadius;
}
-
-
}
void AliPMDUtility::SetPxPyPz(Float_t px, Float_t py, Float_t pz)
AliPMDUtility(Float_t px, Float_t py, Float_t pz);
virtual ~AliPMDUtility();
- void RectGeomCellPos(Int_t ism, Int_t ium,
- Int_t xpad, Int_t ypad,
+ void RectGeomCellPos(Int_t ism, Int_t xpad, Int_t ypad,
Float_t & xpos, Float_t & ypos);
- void RectGeomCellPos(Int_t ism, Int_t ium,
- Float_t xpad, Float_t ypad,
+ void RectGeomCellPos(Int_t ism, Float_t xpad, Float_t ypad,
Float_t & xpos, Float_t & ypos);
void SetPxPyPz(Float_t px, Float_t py, Float_t pz);
void SetXYZ(Float_t xpos, Float_t ypos, Float_t zpos);
Float_t fEta; // Pseudo-rapidity
Float_t fPhi; // Azimuthal angle in radian
- ClassDef(AliPMDUtility,2) // Utility class for the detector set:PMD
+ ClassDef(AliPMDUtility,3) // Utility class for the detector set:PMD
};
#endif
pmddiscriminator->Discrimination(fPMDcontin,fPMDcontout);
const Float_t kzpos = 361.5; // middle of the PMD
- Int_t ism =0, ium=0;
+
Int_t det,smn;
Float_t xpos,ypos;
- Float_t xpad = 0, ypad = 0;
Float_t adc, ncell, rad;
- Float_t xglobal, yglobal, zglobal = 0;
+ Float_t xglobal = 0., yglobal = 0., zglobal = 0;
Float_t pid;
rad = fPMDclout->GetClusRadius();
pid = fPMDclout->GetClusPID();
- //
- // Now change the xpos and ypos to its original values
- // for the unit modules which are earlier changed.
- // xpad and ypad are the real positions.
//
/**********************************************************************
* det : Detector, 0: PRE & 1:CPV *
- * smn : Serial Module Number from which Super Module Number *
- * and Unit Module Numbers are extracted *
+ * smn : Serial Module Number 0 to 23 for each plane *
* xpos : x-position of the cluster *
* ypos : y-position of the cluster *
* THESE xpos & ypos are not the true xpos and ypos *
* adc : ADC contained in the cluster *
* ncell : Number of cells contained in the cluster *
* rad : radius of the cluster (1d fit) *
- * ism : Supermodule number extracted from smn *
- * ium : Unit module number extracted from smn *
- * xpad : TRUE x-position of the cluster *
- * ypad : TRUE y-position of the cluster *
**********************************************************************/
//
- if(det == 0 || det == 1)
- {
- if(smn < 12)
- {
- ism = smn/6;
- ium = smn - ism*6;
- xpad = ypos;
- ypad = xpos;
- }
- else if( smn >= 12 && smn < 24)
- {
- ism = smn/6;
- ium = smn - ism*6;
- xpad = xpos;
- ypad = ypos;
- }
- }
-
- fPMDutil->RectGeomCellPos(ism,ium,xpad,ypad,xglobal,yglobal);
+
+ fPMDutil->RectGeomCellPos(smn,xpos,ypos,xglobal,yglobal);
if (det == 0)
{
zglobal = kzpos - 1.7; // CPV plane
}
-
// Fill ESD
AliESDPmdTrack *esdpmdtr = new AliESDPmdTrack();