Int_t i, j, nmx1, incr, id, jd;
Int_t celldataX[kNmaxCell], celldataY[kNmaxCell];
Int_t celldataTr[kNmaxCell], celldataPid[kNmaxCell];
+ Float_t celldataAdc[kNmaxCell];
Float_t clusdata[6];
Double_t cutoff, ave;
Double_t edepcell[kNMX];
- Double_t *cellenergy = new Double_t [11424];
+ Double_t cellenergy[11424];
// call the isolated cell search method
- CalculateIsoCell(idet, ismn, celladc, pmdisocell);
+ FindIsoCell(idet, ismn, celladc, pmdisocell);
// ndimXr and ndimYr are different because of different module size
}
for (i = 0; i < kNMX; i++)
- {
- edepcell[i] = cellenergy[i];
- }
-
- delete [] cellenergy;
+ {
+ edepcell[i] = cellenergy[i];
+ }
Int_t iord1[kNMX];
TMath::Sort((Int_t)kNMX,edepcell,iord1);// order the data
clusdata[3] = cluCELLS;
clusdata[4] = cluSIGX;
clusdata[5] = cluSIGY;
-
+
//
// Cells associated with a cluster
//
Int_t irow = celldataX[ihit];
Int_t icol = celldataY[ihit];
- if ((irow >= 0 && irow < 48) && (icol >= 0 && icol < 96))
+ if (ismn < 12)
{
- celldataTr[ihit] = celltrack[irow][icol];
- celldataPid[ihit] = cellpid[irow][icol];
+ if ((irow >= 0 && irow < 96) && (icol >= 0 && icol < 48))
+ {
+ celldataTr[ihit] = celltrack[icol][irow];
+ celldataPid[ihit] = cellpid[icol][irow];
+ celldataAdc[ihit] = (Float_t) celladc[icol][irow];
+ }
+ else
+ {
+ celldataTr[ihit] = -1;
+ celldataPid[ihit] = -1;
+ celldataAdc[ihit] = -1;
+ }
}
- else
+ else if (ismn >= 12 && ismn < 24)
{
- celldataTr[ihit] = -1;
- celldataPid[ihit] = -1;
+ if ((irow >= 0 && irow < 48) && (icol >= 0 && icol < 96))
+ {
+ celldataTr[ihit] = celltrack[irow][icol];
+ celldataPid[ihit] = cellpid[irow][icol];
+ celldataAdc[ihit] = (Float_t) celladc[irow][icol];
+
+ }
+ else
+ {
+ celldataTr[ihit] = -1;
+ celldataPid[ihit] = -1;
+ celldataAdc[ihit] = -1;
+ }
}
+
}
-
+
pmdcl = new AliPMDcluster(idet, ismn, clusdata, celldataX, celldataY,
- celldataTr, celldataPid);
+ celldataTr, celldataPid, celldataAdc);
pmdcont->Add(pmdcl);
}
{
clxy[icltr] = -1;
}
+
pmdcludata = new AliPMDcludata(clusdata,clxy);
fPMDclucont->Add(pmdcludata);
}
if(rr==0)
{
- ncell[kcl] = 1.;
+ ncell[kcl] = 1.;
totaladc[kcl] = z1;
- totaladc2[kcl] = pow(z1,2);
- ax[kcl] = x1 * z1;
- ay[kcl] = y1 * z1;
- ax2[kcl]= 0.;
- ay2[kcl]= 0.;
- status[j] = 1;
+ totaladc2[kcl] = z1*z1;
+ ax[kcl] = x1 * z1;
+ ay[kcl] = y1 * z1;
+ ax2[kcl] = 0.;
+ ay2[kcl] = 0.;
+ status[j] = 1;
}
}
}
cellCount[maxweight]++;
- for(Int_t kcl = 0; kcl <= ig; kcl++)
- {
- x2 = xc[kcl];
- y2 = yc[kcl];
- rr = Distance(x1,y1,x2,y2);
- t2 = tc[kcl];
- // For calculating weighted mean:
- // Weighted_mean = (\sum w_i x_i) / (\sum w_i)
-
- if(sumweight>0.)
- {
- totaladc[kcl] = totaladc[kcl] + z1*(weight[kcl]/sumweight);
- ncell[kcl] = ncell[kcl] + (weight[kcl]/sumweight);
- ax[kcl] = ax[kcl] + (x1 * z1 *weight[kcl]/sumweight);
- ay[kcl] = ay[kcl] + (y1 * z1 *weight[kcl]/sumweight);
-
- // For calculating weighted sigma:
- // Wieghted_sigma= ((\sum w_i)/((\sum w_i)^2 - \sum (w_i^2))) \sum w_i (x_i - \Hat\mu)^2
- // I assume here x1,y1 are cluster centers, otherwise life becomes too difficult
- totaladc2[kcl] = totaladc2[kcl] + pow((z1*(weight[kcl]/sumweight)),2);
- ax2[kcl] = ax2[kcl] + (z1 *weight[kcl]/sumweight) * pow((x1-x2),2);
- ay2[kcl] = ay2[kcl] + (z1 *weight[kcl]/sumweight) * pow((y1-y2),2);
- }
- }
+ x2 = xc[maxweight];
+ y2 = yc[maxweight];
+ totaladc[maxweight] += z1;
+ ax[maxweight] += x1*z1;
+ ay[maxweight] += y1*z1;
+ totaladc2[maxweight] += z1*z1;
+ ax2[maxweight] += z1*(x1-x2)*(x1-x2);
+ ay2[maxweight] += z1*(y1-y2)*(y1-y2);
+ ncell[maxweight]++;
+
}
}
for(Int_t kcl = 0; kcl <= ig; kcl++)
{
-
- if(totaladc[kcl]>0){
- if(totaladc[kcl]>0.)xclust[kcl] = (ax[kcl])/ totaladc[kcl];
- if(totaladc[kcl]>0.)yclust[kcl] = (ay[kcl])/ totaladc[kcl];
-
- if(totaladc[kcl]>0.)sigxclust[kcl] = (totaladc[kcl]/(pow(totaladc[kcl],2)-totaladc2[kcl]))*ax2[kcl];
- if(totaladc[kcl]>0.)sigyclust[kcl] = (totaladc[kcl]/(pow(totaladc[kcl],2)-totaladc2[kcl]))*ay2[kcl];
- }
-
+
+ if(totaladc[kcl] > 0.)
+ {
+ xclust[kcl] = (ax[kcl])/ totaladc[kcl];
+ yclust[kcl] = (ay[kcl])/ totaladc[kcl];
+
+ //natasha
+ Float_t sqtotadc = totaladc[kcl]*totaladc[kcl];
+ if(totaladc2[kcl] >= sqtotadc)
+ {
+ sigxclust[kcl] = 0.25;
+ sigyclust[kcl] = 0.25;
+ }
+ else
+ {
+ sigxclust[kcl] = (totaladc[kcl]/(sqtotadc-totaladc2[kcl]))*ax2[kcl];
+ sigyclust[kcl] = (totaladc[kcl]/(sqtotadc-totaladc2[kcl]))*ay2[kcl];
+ }
+ }
+
for(j = 0; j < cellCount[kcl]; j++) clno++;
if (clno >= 4608)
clusdata[1] = yclust[kcl];
clusdata[2] = totaladc[kcl];
clusdata[3] = ncell[kcl];
+
+
if(sigxclust[kcl] > sigyclust[kcl])
{
- clusdata[4] = pow(sigxclust[kcl],0.5);
- clusdata[5] = pow(sigyclust[kcl],0.5);
+ clusdata[4] = TMath::Sqrt(sigxclust[kcl]);
+ clusdata[5] = TMath::Sqrt(sigyclust[kcl]);
}
else
{
- clusdata[4] = pow(sigyclust[kcl],0.5);
- clusdata[5] = pow(sigxclust[kcl],0.5);
+ clusdata[4] = TMath::Sqrt(sigyclust[kcl]);
+ clusdata[5] = TMath::Sqrt(sigxclust[kcl]);
}
-
+
clxy[0] = tc[kcl];
-
+
Int_t Ncell=1;
for (Int_t ii = 0; ii < cellCount[kcl]; ii++)
{
Ncell++;
}
}
-
+
pmdcludata = new AliPMDcludata(clusdata,clxy);
fPMDclucont->Add(pmdcludata);
}
return TMath::Sqrt((x1-x2)*(x1-x2) + (y1-y2)*(y1-y2));
}
// ------------------------------------------------------------------------ //
-void AliPMDClusteringV1::CalculateIsoCell(Int_t idet, Int_t ismn, Double_t celladc[][96], TObjArray *pmdisocell)
+void AliPMDClusteringV1::FindIsoCell(Int_t idet, Int_t ismn, Double_t celladc[][96], TObjArray *pmdisocell)
{
// Does isolated cell search for offline calibration
{
id1 = irow + neibx[ii];
jd1 = icol + neiby[ii];
+ if (id1 < 0) id1 = 0;
+ if (id1 > kMaxRow-1) id1 = kMaxRow - 1;
+ if (jd1 < 0) jd1 = 0;
+ if (jd1 > kMaxCol-1) jd1 = kMaxCol - 1;
Float_t adc = (Float_t) celladc[id1][jd1];
- if(adc == 0.)
+ if(adc < 1.)
{
isocount++;
if(isocount == kCellNeighbour)