{
/// Constructor
- fSegmentation[1] = fSegmentation[0] = 0x0;
+ fkSegmentation[1] = fkSegmentation[0] = 0x0;
if (fPlot) fDebug = 1;
}
const AliMpVSegmentation* seg[2])
{
/// Prepare for clustering
-// AliCodeTimerAuto("")
+// AliCodeTimerAuto("",)
for ( Int_t i = 0; i < 2; ++i )
{
- fSegmentation[i] = seg[i];
+ fkSegmentation[i] = seg[i];
}
// Find out the DetElemId
fSplitter->SetDebug(fDebug);
// find out current event number, and reset the cluster number
- AliRunLoader *runLoader = AliRunLoader::GetRunLoader();
+ AliRunLoader *runLoader = AliRunLoader::Instance();
fEventNumber = runLoader ? runLoader->GetEventNumber() : 0;
fClusterNumber = -1;
fClusterList.Delete();
AliMUONClusterFinderMLEM::NextCluster()
{
/// Return next cluster
-// AliCodeTimerAuto("")
+// AliCodeTimerAuto("",)
// if the list of clusters is not void, pick one from there
TObject* o = fClusterList.At(++fClusterNumber);
/// Starting from a precluster, builds a pixel array, and then
/// extract clusters from this array
- // AliCodeTimerAuto("")
+ // AliCodeTimerAuto("",)
if (fDebug) {
cout << " *** Event # " << fEventNumber
/// Check precluster in order to attempt to simplify it (mostly for
/// two-cathode preclusters)
- AliCodeTimerAuto("")
+ AliCodeTimerAuto("",)
// Disregard small clusters (leftovers from splitting or noise)
if ((origCluster.Multiplicity()==1 || origCluster.Multiplicity()==2) &&
Int_t cath = pad->Cathode();
Int_t cath1 = TMath::Even(cath);
// Check for edge effect (missing pads on the _other_ cathode)
- AliMpPad mpPad = fSegmentation[cath1]->PadByPosition(pad->Position(),kFALSE);
+ AliMpPad mpPad =
+ 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;
cout << " Surface = " << pixelI->Size(0)*pixelI->Size(1)*4 << endl;
pixelJ->Print();
cout << " Surface = " << pixelJ->Size(0)*pixelJ->Size(1)*4 << endl;
- cout << " Area surface = " << area.Dimensions().X()*area.Dimensions().Y()*4 << endl;
+ cout << " Area surface = " << area.GetDimensionX()*area.GetDimensionY()*4 << endl;
cout << "-------" << endl;
);
*/
// AliDebug(2,Form("cluster.Multiplicity=%d",cluster.Multiplicity()));
TVector2 dim = cluster.MinPadDimensions (-1, kFALSE);
- Double_t width[2] = {dim.X(), dim.Y()}, xy0[2];
- Int_t found[2] = {0}, mult = cluster.Multiplicity();
+ Double_t width[2] = {dim.X(), dim.Y()}, xy0[2]={99999,99999};
+ Int_t found[2] = {0,0}, mult = cluster.Multiplicity();
for ( Int_t i = 0; i < mult; ++i) {
AliMUONPad* pad = cluster.Pad(i);
if (cluster.Multiplicity(0) == 0) cath0 = 1;
else if (cluster.Multiplicity(1) == 0) cath1 = 0;
- TVector2 leftDown = cluster.Area(cath0).LeftDownCorner();
- TVector2 rightUp = cluster.Area(cath0).RightUpCorner();
- min[0] = leftDown.X();
- min[1] = leftDown.Y();
- max[0] = rightUp.X();
- max[1] = rightUp.Y();
+
+ Double_t leftDownX, leftDownY;
+ cluster.Area(cath0).LeftDownCorner(leftDownX, leftDownY);
+ Double_t rightUpX, rightUpY;
+ cluster.Area(cath0).RightUpCorner(rightUpX, rightUpY);
+ min[0] = leftDownX;
+ min[1] = leftDownY;
+ max[0] = rightUpX;
+ max[1] = rightUpY;;
if (cath1 != cath0) {
- leftDown = cluster.Area(cath1).LeftDownCorner();
- rightUp = cluster.Area(cath1).RightUpCorner();
- min[0] = TMath::Max (min[0], leftDown.X());
- min[1] = TMath::Max (min[1], leftDown.Y());
- max[0] = TMath::Min (max[0], rightUp.X());
- max[1] = TMath::Min (max[1], rightUp.Y());
+ cluster.Area(cath1).LeftDownCorner(leftDownX, leftDownY);
+ cluster.Area(cath1).RightUpCorner(rightUpX, rightUpY);
+ min[0] = TMath::Max (min[0], leftDownX);
+ min[1] = TMath::Max (min[1], leftDownY);
+ max[0] = TMath::Min (max[0], rightUpX);
+ max[1] = TMath::Min (max[1], rightUpY);
}
// Adjust limits
//width[0] /= 2; width[1] /= 2; // just for check
- Int_t nbins[2];
+ Int_t nbins[2]={0,0};
for (Int_t i = 0; i < 2; ++i) {
Double_t dist = (min[i] - xy0[i]) / width[i] / 2;
if (TMath::Abs(dist) < 1.e-6) dist = -1.e-6;
{
/// Compute coefficients needed for MLEM algorithm
+ Int_t npadTot = cluster.Multiplicity();
Int_t nPix = fPixArray->GetLast()+1;
//memset(probi,0,nPix*sizeof(Double_t));
+ for (Int_t j = 0; j < npadTot*nPix; ++j) coef[j] = 0.;
for (Int_t j = 0; j < nPix; ++j) probi[j] = 0.;
Int_t mult = cluster.Multiplicity();
{
/// Repeat MLEM algorithm until pixel size becomes sufficiently small
- // AliCodeTimerAuto("")
+ // AliCodeTimerAuto("",)
Int_t nPix = fPixArray->GetLast()+1;
for (Int_t ipix = 0; ipix < nPix; ++ipix)
{
- AliMUONPad* pixPtr = Pixel(ipix);
- fHistMlem->Fill(pixPtr->Coord(0),pixPtr->Coord(1),pixPtr->Charge());
+ AliMUONPad* pixPtr2 = Pixel(ipix);
+ fHistMlem->Fill(pixPtr2->Coord(0),pixPtr2->Coord(1),pixPtr2->Charge());
}
// Check if the total charge of pixels is too low
nPix = 0;
for (Int_t ipix = 0; ipix < nPix1; ++ipix)
{
- AliMUONPad* pixPtr = Pixel(ipix);
+ AliMUONPad* pixPtr2 = Pixel(ipix);
if ( nPix >= npadOK // too many pixels already
||
- pixPtr->Charge() < pixMin && pixPtr->Status() != fgkMustKeep // too low charge
+ ((pixPtr2->Charge() < pixMin) && (pixPtr2->Status() != fgkMustKeep)) // too low charge
)
{
RemovePixel(ipix);
{
if (!i)
{
- pixPtr->SetCharge(10);
- pixPtr->SetSize(indx, pixPtr->Size(indx)/2);
- width = -pixPtr->Size(indx);
- pixPtr->Shift(indx, width);
+ pixPtr2->SetCharge(10);
+ pixPtr2->SetSize(indx, pixPtr2->Size(indx)/2);
+ width = -pixPtr2->Size(indx);
+ pixPtr2->Shift(indx, width);
// Shift pixel position
if (ix)
{
ix = 0;
for (Int_t j = 0; j < 2; ++j)
{
- shift[j] = pixPtr->Coord(j) - xyCOG[j];
- shift[j] -= ((Int_t)(shift[j]/pixPtr->Size(j)/2))*pixPtr->Size(j)*2;
+ shift[j] = pixPtr2->Coord(j) - xyCOG[j];
+ shift[j] -= ((Int_t)(shift[j]/pixPtr2->Size(j)/2))*pixPtr2->Size(j)*2;
}
} // if (ix)
- pixPtr->Shift(0, -shift[0]);
- pixPtr->Shift(1, -shift[1]);
+ pixPtr2->Shift(0, -shift[0]);
+ pixPtr2->Shift(1, -shift[1]);
++nPix;
}
else if (nPix < npadOK)
{
- pixPtr = new AliMUONPad(*pixPtr);
- pixPtr->Shift(indx, -2*width);
- pixPtr->SetStatus(fgkZero);
- fPixArray->Add(pixPtr);
+ pixPtr2 = new AliMUONPad(*pixPtr2);
+ pixPtr2->Shift(indx, -2*width);
+ pixPtr2->SetStatus(fgkZero);
+ fPixArray->Add(pixPtr2);
++nPix;
}
else continue; // skip adjustment of histo limits
for (Int_t j = 0; j < 4; ++j)
{
- xylim[j] = TMath::Min (xylim[j], (j%2 ? -1 : 1)*pixPtr->Coord(j/2));
+ xylim[j] = TMath::Min (xylim[j], (j%2 ? -1 : 1)*pixPtr2->Coord(j/2));
}
} // for (Int_t i=0; i<2;
} // for (Int_t ipix=0;
if (nPix < npadOK)
{
- AliMUONPad* pixPtr = Pixel(0);
+ AliMUONPad* pixPtr2 = Pixel(0);
// add pixels if the maximum is at the limit of pixel area:
// start from Y-direction
Int_t j = 0;
for (Int_t i = 3; i > -1; --i)
{
if (nPix < npadOK &&
- TMath::Abs((i%2 ? -1 : 1)*xylim[i]-xyCOG[i/2]) < pixPtr->Size(i/2))
+ TMath::Abs((i%2 ? -1 : 1)*xylim[i]-xyCOG[i/2]) < pixPtr2->Size(i/2))
{
//AliMUONPad* p = static_cast<AliMUONPad*>(pixPtr->Clone());
- AliMUONPad* p = new AliMUONPad(*pixPtr);
- p->SetCoord(i/2, xyCOG[i/2]+(i%2 ? 2:-2)*pixPtr->Size(i/2));
+ AliMUONPad* p = new AliMUONPad(*pixPtr2);
+ p->SetCoord(i/2, xyCOG[i/2]+(i%2 ? 2:-2)*pixPtr2->Size(i/2));
xylim[i] = p->Coord(i/2) * (i%2 ? -1 : 1); // update histo limits
j = TMath::Even (i/2);
p->SetCoord(j, xyCOG[j]);
// Mark pixels which should be removed
for (Int_t i = 0; i < nPix; ++i)
{
- AliMUONPad* pixPtr = Pixel(i);
- charge = pixPtr->Charge();
+ AliMUONPad* pixPtr2 = Pixel(i);
+ charge = pixPtr2->Charge();
if (charge < thresh)
{
- pixPtr->SetCharge(-charge);
+ pixPtr2->SetCharge(-charge);
}
}
Int_t near = 0;
for (Int_t i = 0; i < nPix; ++i)
{
- AliMUONPad* pixPtr = Pixel(i);
- charge = pixPtr->Charge();
+ AliMUONPad* pixPtr2 = Pixel(i);
+ charge = pixPtr2->Charge();
if (charge > 0) continue;
- near = FindNearest(pixPtr);
- pixPtr->SetCharge(0);
+ near = FindNearest(pixPtr2);
+ pixPtr2->SetCharge(0);
probi[i] = 0; // make it "invisible"
AliMUONPad* pnear = Pixel(near);
pnear->SetCharge(pnear->Charge() + (-charge));
// Update histogram
for (Int_t i = 0; i < nPix; ++i)
{
- AliMUONPad* pixPtr = Pixel(i);
- Int_t ix = fHistMlem->GetXaxis()->FindBin(pixPtr->Coord(0));
- Int_t iy = fHistMlem->GetYaxis()->FindBin(pixPtr->Coord(1));
- fHistMlem->SetBinContent(ix, iy, pixPtr->Charge());
+ AliMUONPad* pixPtr2 = Pixel(i);
+ Int_t ix = fHistMlem->GetXaxis()->FindBin(pixPtr2->Coord(0));
+ Int_t iy = fHistMlem->GetYaxis()->FindBin(pixPtr2->Coord(1));
+ fHistMlem->SetBinContent(ix, iy, pixPtr2->Charge());
}
// Try to split into clusters
//_____________________________________________________________________________
void AliMUONClusterFinderMLEM::Mlem(AliMUONCluster& cluster,
- Double_t* coef, Double_t* probi,
+ const Double_t* coef, Double_t* probi,
Int_t nIter)
{
/// Use MLEM to find pixel charges
{
Pixel(ipix)->SetChargeBackup(0);
// Correct each pixel
+ probi1[ipix] = 0;
if (probi[ipix] < 0.01) continue; // skip "invisible" pixel
Double_t sum = 0;
probi1[ipix] = probMax;
}
//_____________________________________________________________________________
-Int_t AliMUONClusterFinderMLEM::FindNearest(AliMUONPad *pixPtr0)
+Int_t AliMUONClusterFinderMLEM::FindNearest(const AliMUONPad *pixPtr0)
{
/// Find the pixel nearest to the given one
/// (algorithm may be not very efficient)
//_____________________________________________________________________________
void AliMUONClusterFinderMLEM::FindCluster(AliMUONCluster& cluster,
- Int_t *localMax, Int_t iMax)
+ const Int_t *localMax, Int_t iMax)
{
/// Find pixel cluster around local maximum \a iMax and pick up pads
/// overlapping with it
Bool_t same = kFALSE;
if (TMath::Abs(dim0.Y()-dim1.Y()) < fgkDistancePrecision) same = kTRUE; // the same pad size on both planes
- AliMpIntPair cn;
+ Long_t cn;
Bool_t check[2] = {kFALSE, kFALSE};
Int_t nxy[2];
nxy[0] = nxy[1] = 0;
for (Int_t inb = 0; inb < 2; ++inb) {
cn = cluster.NofPads(nonb[inb], 0, kTRUE);
- if (inb == 0 && cn.GetFirst() == 2) check[inb] = kTRUE; // check non-bending plane
- else if (inb == 1 && cn.GetSecond() == 2) check[inb] = kTRUE; // check bending plane
+ if (inb == 0 && AliMp::PairFirst(cn) == 2) check[inb] = kTRUE; // check non-bending plane
+ else if (inb == 1 && AliMp::PairSecond(cn) == 2) check[inb] = kTRUE; // check bending plane
if (same) {
- nxy[0] = TMath::Max (nxy[0], cn.GetFirst());
- nxy[1] = TMath::Max (nxy[1], cn.GetSecond());
+ nxy[0] = TMath::Max (nxy[0], AliMp::PairFirst(cn));
+ nxy[1] = TMath::Max (nxy[1], AliMp::PairSecond(cn));
if (inb == 0 && nxy[0] < 2) nonb[inb] = !nonb[inb];
- else if (inb == 1 && cn.GetSecond() < 2) nonb[inb] = !nonb[inb];
+ else if (inb == 1 && AliMp::PairSecond(cn) < 2) nonb[inb] = !nonb[inb];
}
}
if (same) {
TVector2 pos;
if (inb == 0) pos.Set (pad->X() + idir * (pad->DX()+fgkDistancePrecision), pad->Y());
else pos.Set (pad->X(), pad->Y() + idir * (pad->DY()+fgkDistancePrecision));
- //AliMpPad mppad = fSegmentation[nonb[inb]]->PadByPosition(pos,kTRUE);
- AliMpPad mppad = fSegmentation[nonb[inb]]->PadByPosition(pos,kFALSE);
+ //AliMpPad mppad = fkSegmentation[nonb[inb]]->PadByPosition(pos,kTRUE);
+ AliMpPad mppad = fkSegmentation[nonb[inb]]->PadByPosition(pos.X(), pos.Y(),kFALSE);
if (!mppad.IsValid()) continue; // non-existing pad
- cn = mppad.GetIndices();
- AliMUONPad muonPad(fDetElemId, nonb[inb], cn.GetFirst(), cn.GetSecond(),
- mppad.Position().X(), mppad.Position().Y(),
- mppad.Dimensions().X(), mppad.Dimensions().Y(), 0);
+ 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.));
//else muonPad.SetCharge(TMath::Min (amax[j]/15, fgkZeroSuppression));
else muonPad.SetCharge(TMath::Min (amax[j]/15, 6.));
Bool_t mustMatch(kTRUE);
- AliMpIntPair cn = cluster.NofPads(statusToTest,mustMatch);
+ Long_t cn = cluster.NofPads(statusToTest,mustMatch);
- nInX = cn.GetFirst();
- nInY = cn.GetSecond();
+ nInX = AliMp::PairFirst(cn);
+ nInY = AliMp::PairSecond(cn);
}
//_____________________________________________________________________________
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff