// * provided "as is" without express or implied warranty. *
// **************************************************************************
+#include "AliHMPIDCluster.h" //class header
#include <TVirtualFitter.h> //Solve()
#include <TMinuit.h> //Solve()
#include <TClonesArray.h> //Solve()
#include <TMarker.h> //Draw()
-#include "AliLog.h" //FitFunc()
-#include "AliHMPIDCluster.h" //class header
+#include "AliLog.h" //FindCusterSize()
Bool_t AliHMPIDCluster::fgDoCorrSin=kTRUE;
//get L->T cs matrix for a given chamber
const TGeoHMatrix *t2l= AliGeomManager::GetTracking2LocalMatrix(volId);
- if(fParam->GetInstType()) //if there is no geometry we cannot retrieve the volId (only for monitoring)
+ if(!fParam->GetInstType()) //if there is no geometry we cannot retrieve the volId (only for monitoring)
{
new(this) AliCluster3D(); return;
}
// Returns: none
Int_t minPadX=999,minPadY=999,maxPadX=-1,maxPadY=-1; //for box finding
if(fDigs==0) return; //no digits in this cluster
- fXX=fYY=fQRaw=0; //init summable parameters
+ fXX=fYY=fQRaw=0; //init summable parameters
+ fCh = -1; //init chamber
Int_t maxQpad=-1,maxQ=-1; //to calculate the pad with the highest charge
AliHMPIDDigit *pDig=0x0;
for(Int_t iDig=0;iDig<fDigs->GetEntriesFast();iDig++){ //digits loop
pDig=(AliHMPIDDigit*)fDigs->At(iDig); //get pointer to next digit
-
+ if(!pDig) continue; //protection
if(pDig->PadPcX() > maxPadX) maxPadX = pDig->PadPcX(); // find the minimum box that contain the cluster MaxX
if(pDig->PadPcY() > maxPadY) maxPadY = pDig->PadPcY(); // MaxY
if(pDig->PadPcX() < minPadX) minPadX = pDig->PadPcX(); // MinX
fXX += pDig->LorsX()*q;fYY +=pDig->LorsY()*q; //add digit center weighted by QDC
fQRaw+=q; //increment total charge
if(q>maxQ) {maxQpad = pDig->Pad();maxQ=(Int_t)q;} // to find pad with highest charge
+ fCh=pDig->Ch(); //initialize chamber number
}//digits loop
fBox=(maxPadX-minPadX+1)*100+maxPadY-minPadY+1; // dimension of the box: format Xdim*100+Ydim
if(fDigs->GetEntriesFast()>1&&fgDoCorrSin)CorrSin(); //correct it by sinoid
fQ = fQRaw; // Before starting fit procedure, Q and QRaw must be equal
- fCh=pDig->Ch(); //initialize chamber number
fMaxQpad = maxQpad; fMaxQ=maxQ; //store max charge pad to the field
fChi2=0; // no Chi2 to find
fNlocMax=0; // proper status from this method
for(Int_t i=0;i<iNpars;i++) delete [] derivPart[i]; delete [] derivPart;
}
//---gradient calculations ended
+
+// fit ended. Final calculations
+
}//FitFunction()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
case kSi1 : status="size 1 (cog)" ;break;
case kNoLoc: status="no LocMax(fit)" ;break;
case kAbn : status="Abnormal fit " ;break;
+ case kBig : status="Big Clu(>100) " ;break;
default: status="??????" ;break;
}
Double_t ratio=0;
if(Q()>0&&QRaw()>0) ratio = Q()/QRaw()*100;
- Printf("%sCLU: ch=%i (%7.3f,%7.3f) Q=%8.3f Qraw=%8.3f(%3.0f%%) Size=%2i DimBox=%i LocMax=%i Chi2=%7.3f %s",
+ Printf("%sCLU: ch=%i (%7.3f,%7.3f) Q=%8.3f Qraw=%8.3f(%3.0f%%) Size=%2i DimBox=%i LocMax=%i Chi2=%7.3f %s",
opt,Ch(),X(),Y(),Q(),QRaw(),ratio,Size(),fBox,fNlocMax,fChi2,status);
if(fDigs) fDigs->Print();
}//Print()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
+Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Int_t *pSigmaCut, Bool_t isTryUnfold)
{
//This methode is invoked when the cluster is formed to solve it. Solve the cluster means to try to unfold the cluster
//into the local maxima number of clusters. This methode is invoked by AliHMPIDRconstructor::Dig2Clu() on cluster by cluster basis.
const Int_t kMaxLocMax=6; //max allowed number of loc max for fitting
//
CoG(); //First calculate CoG for the given cluster
+
Int_t iCluCnt=pCluLst->GetEntriesFast(); //get current number of clusters already stored in the list by previous operations
- if(isTryUnfold==kFALSE || Size()==1) { //if cluster contains single pad there is no way to improve the knowledge
- fSt = (isTryUnfold)? kSi1: kNot;
- if(fParam->GetInstType()) SetClusterParams(fXX,fYY,fCh);
+
+ Int_t rawSize = Size(); //get current raw cluster size
+
+ if(rawSize>100) {
+ fSt = kBig;
+ } else if(isTryUnfold==kFALSE) {
+ fSt = kNot;
+ } else if(rawSize==1) {
+ fSt = kSi1;
+ }
+
+ if(rawSize>100 || isTryUnfold==kFALSE || rawSize==1) { //No deconv if: 1 - big cluster (also avoid no zero suppression!)
+ // 2 - flag is set to FALSE
+ if(fParam->GetInstType()) SetClusterParams(fXX,fYY,fCh); // 3 - size = 1
new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this); //add this raw cluster
return 1;
+
}
//Phase 0. Initialise Fitter
Int_t ierflg = 0;
TVirtualFitter *fitter = TVirtualFitter::Fitter(this,3*6); //initialize Fitter
+ delete fitter; //temporary solution to avoid the inteference with previous instances
+ fitter = TVirtualFitter::Fitter(this,3*6); //initialize Fitter
+
arglist[0] = -1;
ierflg = fitter->ExecuteCommand("SET PRI", arglist, 1); // no printout
ierflg = fitter->ExecuteCommand("SET NOW", arglist, 0); //no warning messages
//Phase 1. Find number of local maxima. Strategy is to check if the current pad has QDC more then all neigbours. Also find the box contaning the cluster
fNlocMax=0;
- for(Int_t iDig1=0;iDig1<Size();iDig1++) { //first digits loop
+ for(Int_t iDig1=0;iDig1<rawSize;iDig1++) { //first digits loop
AliHMPIDDigit *pDig1 = Dig(iDig1); //take next digit
Int_t iCnt = 0; //counts how many neighbouring pads has QDC more then current one
- for(Int_t iDig2=0;iDig2<Size();iDig2++) { //loop on all digits again
+ for(Int_t iDig2=0;iDig2<rawSize;iDig2++) { //loop on all digits again
if(iDig1==iDig2) continue; //the same digit, no need to compare
AliHMPIDDigit *pDig2 = Dig(iDig2); //take second digit to compare with the first one
Double_t dummy; char sName[80]; //vars to get results from Minuit
Double_t edm, errdef;
Int_t nvpar, nparx;
-
+
for(Int_t i=0;i<fNlocMax;i++){ //store the local maxima parameters
- fitter->GetParameter(3*i ,sName, fXX, fErrX , dummy, dummy); // X
- fitter->GetParameter(3*i+1 ,sName, fYY, fErrY , dummy, dummy); // Y
+ fitter->GetParameter(3*i ,sName, fXX, fErrX , dummy, dummy); // X
+ fitter->GetParameter(3*i+1 ,sName, fYY, fErrY , dummy, dummy); // Y
fitter->GetParameter(3*i+2 ,sName, fQ, fErrQ , dummy, dummy); // Q
fitter->GetStats(fChi2, edm, errdef, nvpar, nparx); //get fit infos
- if(fSt!=kAbn) {
- if(fNlocMax!=1)fSt=kUnf; // if unfolded
- if(fNlocMax==1&&fSt!=kNoLoc) fSt=kLo1; // if only 1 loc max
- if ( !IsInPc()) fSt = kEdg; // if Out of Pc
- if(fSt==kNoLoc) fNlocMax=0; // if with no loc max (pads with same charge..)
- }
- if(fParam->GetInstType()) SetClusterParams(fXX,fYY,fCh); //need to fill the AliCluster3D part
- new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this); //add new unfolded cluster
-
+
+ if(fNlocMax>1)FindClusterSize(i,pSigmaCut); //find clustersize for deconvoluted clusters
+ //after this call, fSi temporarly is the calculated size. Later is set again
+ //to its original value
+ if(fSt!=kAbn) {
+ if(fNlocMax!=1)fSt=kUnf; // if unfolded
+ if(fNlocMax==1&&fSt!=kNoLoc) fSt=kLo1; // if only 1 loc max
+ if ( !IsInPc()) fSt = kEdg; // if Out of Pc
+ if(fSt==kNoLoc) fNlocMax=0; // if with no loc max (pads with same charge..)
+ }
+ if(fParam->GetInstType()) SetClusterParams(fXX,fYY,fCh); //need to fill the AliCluster3D part
+ new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this); //add new unfolded cluster
+ if(fNlocMax>1)SetSize(rawSize); //Original raw size is set again to its proper value
}
}
}//Solve()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDCluster::FindClusterSize(Int_t i,Int_t *pSigmaCut)
+{
+
+//Estimate of the clustersize for a deconvoluted cluster
+ Int_t size = 0;
+ for(Int_t iDig=0;iDig<Size();iDig++) { //digits loop
+ AliHMPIDDigit *pDig = Dig(iDig); //take digit
+ Int_t iCh = pDig->Ch();
+ Double_t qPad = Q()*pDig->IntMathieson(X(),Y()); //pad charge
+ AliDebug(1,Form("Chamber %i X %i Y %i SigmaCut %i pad %i qpadMath %8.2f qPadRaw %8.2f Qtotal %8.2f cluster n.%i",iCh,pDig->PadChX(),pDig->PadChY(),
+ pSigmaCut[iCh],iDig,qPad,pDig->Q(),QRaw(),i));
+ if(qPad>pSigmaCut[iCh]) size++;
+ }
+ AliDebug(1,Form(" Calculated size %i",size));
+ if(size>0) SetSize(size); //in case of size == 0, original raw clustersize used
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff