1 // ---- ROOT system ---
2 #include "TDirectory.h"
8 // --- AliRoot header files ---
9 #include "AliPHOSEvalRecPoint.h"
11 #include "AliPHOSGetter.h"
12 #include "AliPHOSRecCpvManager.h"
13 #include "AliPHOSRecEmcManager.h"
14 #include "AliPHOSDigitizer.h"
16 // --- Standard library ---
19 ClassImp(AliPHOSEvalRecPoint)
21 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(): AliPHOSCpvRecPoint()
29 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Bool_t cpv, AliPHOSEvalRecPoint* parent) : AliPHOSCpvRecPoint()
36 TObjArray* wPool = (TObjArray*)GetWorkingPool();
38 Error("AliPHOSEvalRecPoint", "Couldn't find working pool. Exit.") ;
42 fParent = wPool->IndexOf((TObject*)parent);
43 fChi2Dof = parent->Chi2Dof();
54 // Add itself to working pool
55 AddToWorkingPool((TObject*)this);
59 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Int_t i, Bool_t cpv) : AliPHOSCpvRecPoint()
65 AliPHOSEmcRecPoint* rp=0;
67 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
70 rp = (AliPHOSCpvRecPoint*)fGetter->CpvRecPoints()->At(i);
75 rp = (AliPHOSEmcRecPoint*)fGetter->EmcRecPoints()->At(i);
80 Int_t* Digits = rp->GetDigitsList();
81 Float_t* Energies = rp->GetEnergiesList();
82 Int_t nDigits = rp->GetMultiplicity();
84 for(Int_t iDigit=0; iDigit<nDigits; iDigit++) {
85 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
86 Float_t eDigit = Energies[iDigit];
87 this->AddDigit(*digit,eDigit);
91 rp->GetLocalPosition(locpos);
94 // Add itself to working pool
95 AddToWorkingPool((TObject*)this);
99 AliPHOSClusterizer* AliPHOSEvalRecPoint::GetClusterizer()
101 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
102 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
103 AliPHOSClusterizer* clu = (AliPHOSClusterizer*)wPoolF->FindObject("PHOS:clu-v1");
105 Error("GetClusterizer", "Couldn't find Clusterizer. Exit.") ;
112 Bool_t AliPHOSEvalRecPoint::TooClose(AliPHOSRecPoint* pt) const
116 pt->GetLocalPosition(her_pos);
117 this->GetLocalPosition(my_pos);
118 Float_t dx = her_pos.X() - my_pos.X();
119 Float_t dz = her_pos.Z() - my_pos.Z();
120 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
122 if(dr<GetReconstructionManager()->MergeGammasMinDistanceCut())
129 Bool_t AliPHOSEvalRecPoint::NeedToSplit() const
134 void AliPHOSEvalRecPoint::DeleteParent()
139 void AliPHOSEvalRecPoint::UpdateWorkingPool()
142 Int_t i; //loop variable
144 for(i=0; i<this->InWorkingPool(); i++) {
145 AliPHOSEvalRecPoint* parent = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
147 Int_t nChild = parent->HasChild(children);
148 for(Int_t iChild=0; iChild<nChild; iChild++)
150 ((AliPHOSEvalRecPoint*)children.At(iChild))->DeleteParent();
154 RemoveFromWorkingPool(parent);
160 for(i=0; i<this->InWorkingPool(); i++) {
161 AliPHOSEvalRecPoint* weak = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
162 if (weak->KillWeakPoint()) delete weak;
165 for(i=0; i<this->InWorkingPool(); i++) {
166 AliPHOSEvalRecPoint* close = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
167 close->MergeClosePoint();
170 for(i=0; i<this->InWorkingPool(); i++)
171 ((AliPHOSEvalRecPoint*)AliPHOSEvalRecPoint::GetFromWorkingPool(i))->SetIndexInList(i);
175 Int_t AliPHOSEvalRecPoint::HasChild(TObjArray& children)
177 for(Int_t iChild=0; iChild<InWorkingPool(); iChild++)
179 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)GetFromWorkingPool(iChild);
180 TObject* parent = (TObject*)child->Parent();
181 TObject* me = (TObject*)this;
182 if(parent==me) children.Add(child);
185 return children.GetEntries();
188 void AliPHOSEvalRecPoint::Init()
191 AliPHOSClusterizer* clusterizer = GetClusterizer();
193 Error("Init", "Cannot get clusterizer. Exit.") ;
197 TClonesArray* digits = AliPHOSGetter::GetInstance()->Digits();
201 LogWeight = clusterizer->GetEmcLogWeight();
204 LogWeight = clusterizer->GetCpvLogWeight();
207 EvalLocalPosition(LogWeight,digits); // evaluate initial position
211 void AliPHOSEvalRecPoint::MakeJob()
213 // Reconstruction algoritm implementation.
215 AliPHOSRecManager* recMng = GetReconstructionManager();
222 Int_t nChild = HasChild(children);
226 if(Chi2Dof()>recMng->OneGamChisqCut()) SplitMergedShowers();
229 for(Int_t iChild=0; iChild<nChild; iChild++) {
230 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)children.At(iChild);
231 child->EvaluatePosition();
233 if(child->Chi2Dof()>recMng->OneGamChisqCut())
234 child->SplitMergedShowers();
239 void AliPHOSEvalRecPoint::InitTwoGam(Float_t* gamma1, Float_t* gamma2)
241 //Compute start values for two gamma fit algorithm.
242 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed "gammas".
245 TVector3 lpos; // start point choosing.
246 GetLocalPosition(lpos);
247 Float_t xx = lpos.Z();
248 Float_t yy = lpos.X();
249 Float_t E = GetEnergy();
251 Info("InitTwoGam", "(x,z,E)[old] = (%f, %f, %f)", yy, xx, E) ;
258 AliPHOSDigit * digit ;
259 Int_t nDigits = GetMultiplicity();
260 Int_t * Digits = GetDigitsList();
261 Float_t * Energies = GetEnergiesList();
273 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
274 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
276 Int_t iDigit; //loop variable
278 for(iDigit = 0 ; iDigit < nDigits ; iDigit ++)
280 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
281 eDigit = Energies[iDigit];
282 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
283 fGeom->RelPosInModule(relid, iy, ix);
285 Float_t dx = ix - xx;
286 Float_t dy = iy - yy;
293 sqr = TMath::Sqrt(4.*exy*exy + (exx-eyy)*(exx-eyy));
294 Float_t euu = (exx+eyy+sqr)/2.;
295 if(sqr>1.e-10) cos2fi = (exx - eyy)/sqr;
296 Float_t cosfi = TMath::Sqrt((1.+cos2fi)/2.);
297 Float_t sinfi = TMath::Sqrt((1.-cos2fi)/2.);
298 if(exy<0) sinfi = -sinfi;
301 for(iDigit = 0 ; iDigit < nDigits ; iDigit ++)
303 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
304 eDigit = Energies[iDigit];
305 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
306 fGeom->RelPosInModule(relid, iy, ix);
308 Float_t dx = ix - xx;
309 Float_t dy = iy - yy;
310 Float_t du = dx*cosfi + dy*sinfi;
311 eu3 += eDigit*du*du*du;
314 Float_t c = E*eu3*eu3/(euu*euu*euu)/2.;
316 if(eu3<0) sign = -1.;
317 Float_t r = 1.+ c + sign*TMath::Sqrt(2.*c + c*c);
319 if(TMath::Abs(r-1.)>0.1) u = eu3/euu*(r+1.)/(r-1.);
320 if(TMath::Abs(r-1.)<0.1) u = TMath::Sqrt(sqr/E/r)*(1.+r);
322 Float_t e2c = E/(1.+r);
324 Float_t u1 = -u/(1.+r);
327 Float_t x1c = xx+u1*cosfi;
328 Float_t y1c = yy+u1*sinfi;
329 Float_t x2c = xx+u2*cosfi;
330 Float_t y2c = yy+u2*sinfi;
332 // printf("e1c -> %f\n",e1c);
333 // printf("x1c -> %f\n",x1c);
334 // printf("y1c -> %f\n",y1c);
335 // printf("e2c -> %f\n",e2c);
336 // printf("x2c -> %f\n",x2c);
337 // printf("y2c -> %f\n",y2c);
349 void AliPHOSEvalRecPoint::TwoGam(Float_t* gamma1, Float_t* gamma2)
351 //Fitting algorithm for the two very closed gammas
352 //that merged into the cluster with one maximum.
353 // Starting points gamma1 and gamma2 must be provided before ( see InitTwoGam)
354 //Set chisquare of the fit.
357 Float_t st0 = GetReconstructionManager()->TwoGamInitialStep();
358 Float_t chmin = GetReconstructionManager()->TwoGamChisqMin();
359 Float_t emin = GetReconstructionManager()->TwoGamEmin();
360 Float_t stpmin = GetReconstructionManager()->TwoGamStepMin();
361 Int_t Niter = GetReconstructionManager()->TwoGamNumOfIterations(); // Number of iterations.
363 Float_t chisq = 100.; //Initial chisquare.
365 Int_t nadc = GetMultiplicity();
369 Int_t dof = nadc - 5;
371 Float_t chstop = chmin*dof;
381 Float_t ee1=0,xx1=0,yy1=0,ee2=0,xx2=0,yy2=0,cosi;
383 Float_t e1c = gamma1[0];
384 Float_t y1c = gamma1[1];
385 Float_t x1c = gamma1[2];
387 Float_t e2c = gamma2[0];
388 Float_t y2c = gamma2[1];
389 Float_t x2c = gamma2[2];
391 Float_t E = GetEnergy();
394 AliPHOSDigit * digit ;
395 Int_t nDigits = GetMultiplicity();
396 Int_t * Digits = GetDigitsList();
397 Float_t * Energies = GetEnergiesList();
403 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
404 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
406 for(Int_t Iter=0; Iter<Niter; Iter++)
415 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
417 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
418 eDigit = Energies[iDigit];
419 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
420 fGeom->RelPosInModule(relid, iy, ix);
425 Float_t dx1 = x1c - ix;
426 Float_t dy1 = y1c - iy;
427 // Info("TwoGam", "Mult %d dx1 %f dy1 %f", nDigits, dx1, dy1) ;
428 // AG(e1c,dx1,dy1,a1,gx1,gy1);
429 GetReconstructionManager()->AG(e1c,dx1,dy1,a1,gx1,gy1);
431 Float_t dx2 = x2c - ix;
432 Float_t dy2 = y2c - iy;
433 // Info("TwoGam", " dx2 %f dy2 %f", dx2, dy2) ;
434 // AG(e2c,dx2,dy2,a2,gx2,gy2);
435 GetReconstructionManager()->AG(e2c,dx2,dy2,a2,gx2,gy2);
438 // Float_t D = Const*A*(1. - A/E);
441 // // D = 9.; // ????
443 // Float_t da = A - eDigit;
444 // chisqc += da*da/D;
445 // Float_t dd = da/D;
446 // dd = dd*(2.-dd*Const*(1.-2.*A/E));
449 chisqc += GetReconstructionManager()->TwoGamChi2(A,eDigit,E,dd);
455 grec += (a1/e1c - a2/e2c)*E*dd;
459 Float_t grc = TMath::Sqrt(grx1c*grx1c + gry1c*gry1c + grx2c*grx2c + gry2c*gry2c);
460 if(grc<1.e-10) grc=1.e-10;
462 Float_t sc = 1. + chisqc/ch;
467 cosi = (grx1*grx1c + gry1*gry1c + grx2*grx2c + gry2*gry2c + gre*grec)/gr/grc;
468 st = st*sc/(1.4 - cosi);
490 Float_t step = st*gr;
492 Info("TwoGam", "Iteration %d dof %d chisq/dof %f chstop/dof %f step %d stpmin %d",
493 Iter, dof, ch/dof, chstop/dof, step, stpmin) ;
499 Float_t de = st*gre*E;
503 if( (e1c>emin) && (e2c>emin) )
519 // if(ch>chisq*(nadc-2)-delch)
531 Float_t x1_new = yy1;
532 Float_t z1_new = xx1;
533 Float_t e1_new = ee1;
534 Float_t x2_new = yy2;
535 Float_t z2_new = xx2;
536 Float_t e2_new = ee2;
539 message = " (x,z,E)[1 fit] = (%f, %f, %f)\n" ;
540 message = " (x,z,E)[2 fit] = (%f, %f, %f)\n" ;
541 Info("TwoGam", message.Data(),
542 x1_new, z1_new, e1_new,
543 x2_new, z2_new, e2_new) ;
549 void AliPHOSEvalRecPoint::UnfoldTwoMergedPoints(Float_t* gamma1, Float_t* gamma2)
551 //Unfold TWO merged rec. points in the case when cluster has only one maximum,
552 //but it's fitting to the one gamma shower is too bad.
553 // Use TwoGam() to estimate the positions and energies of merged points.
559 Int_t* Digits = GetDigitsList();
560 Int_t nDigits = GetMultiplicity();
561 Float_t* Energies = GetEnergiesList();
562 Float_t* eFit = new Float_t[nDigits];
564 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
565 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
567 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
569 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
571 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
573 fGeom->RelPosInModule(relid, x, z);
576 for(Int_t iMax=0; iMax<Nmax; iMax++)
583 Float_t eMax = gamma[0];
584 Float_t xMax = gamma[1];
585 Float_t zMax = gamma[2];
587 Float_t dx = xMax - x;
588 Float_t dz = zMax - z;
590 GetReconstructionManager()->AG(eMax,dz,dx,Amp,gx,gy);
598 for( Int_t iMax=0; iMax<Nmax; iMax++)
605 Float_t eMax = gamma[0];
606 Float_t xMax = gamma[1];
607 Float_t zMax = gamma[2];
609 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
610 TVector3 newpos(xMax,0,zMax);
611 newRP->SetLocalPosition(newpos);
613 for( Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
615 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
616 Float_t eDigit = Energies[iDigit];
618 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
620 fGeom->RelPosInModule(relid, ix, iz);
622 Float_t dx = xMax - ix;
623 Float_t dz = zMax - iz;
624 Float_t single_shower_gain,gxMax,gyMax;
625 GetReconstructionManager()->AG(eMax,dz,dx,single_shower_gain,gxMax,gyMax);
626 Float_t total_gain = eFit[iDigit];
627 Float_t ratio = single_shower_gain/total_gain;
628 eDigit = eDigit*ratio;
629 newRP->AddDigit(*digit,eDigit);
631 Info("UnfoldTwoMergedPoints", "======= Split: daughter rec point %d =================", iMax) ;
642 void AliPHOSEvalRecPoint::EvaluatePosition()
644 // One gamma fit algorithm.
645 // Set chisq/dof of the fit.
646 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed gamma, respectively.
648 Int_t nDigits = GetMultiplicity();
652 Int_t Niter = GetReconstructionManager()->OneGamNumOfIterations(); // number of iterations
653 Float_t St0 = GetReconstructionManager()->OneGamInitialStep(); // initial step
654 // const Float_t Stpmin = 0.005;
655 Float_t Stpmin = GetReconstructionManager()->OneGamStepMin();
656 Float_t chmin = GetReconstructionManager()->OneGamChisqMin();
664 GetLocalPosition(locpos);
665 Float_t E = GetEnergy();
666 Float_t xc = locpos.Z();
667 Float_t yc = locpos.X();
668 Float_t dx,dy,gx,gy,grxc,gryc;
670 Float_t chisq = 1.e+20;
674 Int_t dof = GetMultiplicity() - 2;
677 Float_t chstop = chmin*dof;
678 Float_t cosi,x1=0,y1=0;
681 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
682 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
684 for(Int_t Iter=0; Iter<Niter; Iter++)
690 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
693 Float_t* Energies = GetEnergiesList();
694 Int_t* Digits = GetDigitsList();
695 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
696 eDigit = Energies[iDigit];
697 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
698 fGeom->RelPosInModule(relid, iy, ix);
707 GetReconstructionManager()->AG(E,dx,dy,A,gx,gy);
709 Info("EvaluatePosition", " (ix iy xc yc dx dy) %f %f %f %f %f %f", ix, iy, xc, yc, dx, dy) ;
710 Float_t chi2dg = GetReconstructionManager()->OneGamChi2(A,eDigit,E,dd);
712 // Exclude digit with too large chisquare.
713 if(chi2dg > 10) { continue; }
721 Float_t grc = TMath::Sqrt(grxc*grxc + gryc*gryc);
722 if(grc<1.e-10) grc=1.e-10;
723 Float_t sc = 1. + chisqc/chisq;
724 Info("EvaluatePosition", " chisq: %f", chisq) ;
728 cosi = (grx*grxc + gry*gryc)/gr/grc;
729 st = st*sc/(1.4 - cosi);
741 Float_t step = st*gr;
743 Info("EvaluatePosition", " Iteration %d dof %d chisq/dof %f chstop/dof %f step %d Stpmin %d",
744 Iter, dof, chisq/dof, chisq/dof, chstop/dof, step, Stpmin) ;
769 TVector3 newpos(gamma1[1],0,gamma1[2]);
770 //SetLocalPosition(newpos);
776 // RP->GetLocalPosition(pos);
785 Bool_t AliPHOSEvalRecPoint::KillWeakPoint()
787 // Remove this point from further procession
788 // if it's energy is too small.
790 Float_t Thr0 = GetReconstructionManager()->KillGamMinEnergy();
792 if(GetEnergy()<Thr0) {
793 Info("KillWeakPoint", "+++++++ Killing this rec point ++++++++++") ;
794 RemoveFromWorkingPool(this);
802 void AliPHOSEvalRecPoint::SplitMergedShowers()
804 // Split two very close rec. points.
806 if(GetMultiplicity()<2)
812 InitTwoGam(gamma1,gamma2); // start points for fit
813 TwoGam(gamma1,gamma2); // evaluate the positions and energies
814 UnfoldTwoMergedPoints(gamma1,gamma2); // make the job
819 void AliPHOSEvalRecPoint::MergeClosePoint()
822 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
823 // AliPHOSDigitizer* digitizer = fGetter->Digitizer();
824 // Float_t fPedestal = digitizer->GetPedestal(); // YVK 30.09.2001
825 // Float_t fSlope = digitizer->GetSlope();
827 for (Int_t i=0;i<InWorkingPool(); i++)
829 TObject* obj = GetFromWorkingPool(i);
830 if(!((TObject*)this)->IsEqual(obj))
832 AliPHOSRecPoint* rp = (AliPHOSRecPoint*)obj;
833 if(GetPHOSMod() == rp->GetPHOSMod())
837 Info("MergeClosePoint", "+++++++ Merging point 1: ++++++") ;
839 Info("MergeClosePoint", "+++++++ and point 2: ++++++++++") ;
840 ((AliPHOSEvalRecPoint*)rp)->Print("");
842 //merge two rec. points
845 this->GetLocalPosition(lpos1);
846 rp->GetLocalPosition(lpos2);
848 Int_t* Digits = rp->GetDigitsList();
849 Float_t dE = rp->GetEnergy()/(rp->GetEnergy()+this->GetEnergy());
850 Float_t new_x = lpos1.X()*dE + lpos2.X()*(1.-dE);
851 Float_t new_z = lpos1.Z()*dE + lpos2.Z()*(1.-dE);
852 Float_t new_E = rp->GetEnergy()+this->GetEnergy();
853 Float_t* Energies = ((AliPHOSEmcRecPoint*)rp)->GetEnergiesList();
855 for(Int_t iDigit=0; iDigit<rp->GetDigitsMultiplicity(); iDigit++)
857 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At(Digits[iDigit]);
858 Float_t eDigit = Energies[iDigit];
859 this->AddDigit(*digit,eDigit);
862 TVector3 newpos(new_x,0,new_z);
865 RemoveFromWorkingPool(rp);
868 Info("MergeClosePoint", "++++++ Resulting point: ++++++++") ;
879 Int_t AliPHOSEvalRecPoint::UnfoldLocalMaxima()
881 // Make unfolding in the reconstruction point with several local maxima.
882 // Return the number of local maxima.
884 // if multiplicity less then 2 - nothing to unfold
885 if(GetMultiplicity()<2) return 1;
887 AliPHOSDigit * maxAt[1000];
888 Float_t maxAtEnergy[1000];
889 Float_t LocMaxCut, LogWeight;
894 // AliPHOSClusterizer* clusterizer = fGetter->Clusterizer("PHOSclu-v1");
895 AliPHOSClusterizer* clusterizer = GetClusterizer();
897 Error("UnfoldLocalMaxima", "Cannot get clusterizer. Exit.") ;
902 LocMaxCut = clusterizer->GetEmcLocalMaxCut();
903 LogWeight = clusterizer->GetEmcLogWeight();
906 LocMaxCut = clusterizer->GetCpvLocalMaxCut();
907 LogWeight = clusterizer->GetCpvLogWeight();
910 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
911 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
912 TClonesArray* digits = fGetter->Digits();
914 // if number of local maxima less then 2 - nothing to unfold
915 Int_t Nmax = GetNumberOfLocalMax(maxAt,maxAtEnergy,LocMaxCut,digits);
916 if(Nmax<2) return Nmax;
918 Int_t* Digits = GetDigitsList();
919 Int_t nDigits = GetMultiplicity();
920 Float_t* Energies = GetEnergiesList();
921 Float_t* eFit = new Float_t[nDigits];
923 Int_t iDigit; //loop variable
925 for(iDigit=0; iDigit<nDigits; iDigit++)
930 for(iDigit=0; iDigit<nDigits; iDigit++)
933 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
934 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
936 fGeom->RelPosInModule(relid, x, z);
938 for(Int_t iMax=0; iMax<Nmax; iMax++)
940 AliPHOSDigit* digitMax = maxAt[iMax];
941 Float_t eMax = maxAtEnergy[iMax];
942 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
943 fGeom->RelPosInModule(relid, xMax, zMax);
944 Float_t dx = xMax - x;
945 Float_t dz = zMax - z;
947 GetReconstructionManager()->AG(eMax,dz,dx,Amp,gx,gy);
954 for(Int_t iMax=0; iMax<Nmax; iMax++)
956 AliPHOSDigit* digitMax = maxAt[iMax];
957 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
958 fGeom->RelPosInModule(relid, xMax, zMax);
959 Float_t eMax = maxAtEnergy[iMax];
961 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
962 newRP->AddDigit(*digitMax,maxAtEnergy[iMax]);
964 //Neighbous ( matrix 3x3 around the local maximum)
965 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
967 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
968 Float_t eDigit = Energies[iDigit];
969 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
971 fGeom->RelPosInModule(relid, ix, iz);
973 if(AreNeighbours(digitMax,digit))
975 Float_t dx = xMax - ix;
976 Float_t dz = zMax - iz;
977 Float_t single_shower_gain,gxMax,gyMax;
978 GetReconstructionManager()->AG(eMax,dz,dx,single_shower_gain,gxMax,gyMax);
979 Float_t total_gain = eFit[iDigit];
980 Float_t ratio = single_shower_gain/total_gain;
981 Info("UnfoldLocalMaxima", " ratio -> %f", ratio) ;
982 eDigit = eDigit*ratio;
983 newRP->AddDigit(*digit,eDigit);
987 newRP->EvalLocalPosition(LogWeight,digits);
988 Info("UnfoldLocalMaxima", "======= Unfold: daughter rec point %d =================", iMax) ;
992 // RemoveFromWorkingPool(this);
1000 void AliPHOSEvalRecPoint::PrintPoint(Option_t* opt)
1002 AliPHOSCpvRecPoint::Print(opt);
1005 GetLocalPosition(lpos);
1008 message = " Chi2/dof = %f" ;
1009 message += " Local (x,z) = (%f, %f) in module %d" ;
1010 Info("Print", message.Data(), Chi2Dof(), lpos.X(), lpos.Z(), GetPHOSMod()) ;
1014 AliPHOSRecManager* AliPHOSEvalRecPoint::GetReconstructionManager() const
1017 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1018 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1019 AliPHOSRecManager* recMng = (AliPHOSRecManager*)wPoolF->FindObject("AliPHOSRecManager");
1021 Error("GetReconstructionManager", "Couldn't find Reconstruction Manager. Exit.") ;
1029 AliPHOSRecPoint* AliPHOSEvalRecPoint::Parent()
1031 if(fParent<0) return NULL;
1033 return (AliPHOSRecPoint*)GetFromWorkingPool(fParent);
1037 Float_t AliPHOSEvalRecPoint::Chi2Dof() const
1042 const TObject* AliPHOSEvalRecPoint::GetWorkingPool()
1045 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1046 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1047 TObject* wPool = wPoolF->FindObject("SmartPoints");
1049 Error("GetWorkingPool", "Couldn't find Working Pool. Exit.") ;
1056 void AliPHOSEvalRecPoint::AddToWorkingPool(TObject* obj)
1058 ((TObjArray*)GetWorkingPool())->Add(obj);
1061 TObject* AliPHOSEvalRecPoint::GetFromWorkingPool(Int_t i)
1063 // return fWorkingPool->At(i);
1064 return ((TObjArray*)GetWorkingPool())->At(i);
1067 Int_t AliPHOSEvalRecPoint::InWorkingPool()
1069 return ((TObjArray*)GetWorkingPool())->GetEntries();
1072 void AliPHOSEvalRecPoint::RemoveFromWorkingPool(TObject* obj)
1074 ((TObjArray*)GetWorkingPool())->Remove(obj);
1075 ((TObjArray*)GetWorkingPool())->Compress();
1078 void AliPHOSEvalRecPoint::PrintWorkingPool()
1080 ((TObjArray*)GetWorkingPool())->Print("");