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 ---
20 ClassImp(AliPHOSEvalRecPoint)
22 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(): AliPHOSCpvRecPoint()
30 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Bool_t cpv, AliPHOSEvalRecPoint* parent) : AliPHOSCpvRecPoint()
37 TObjArray* wPool = (TObjArray*)GetWorkingPool();
39 cout<<" Couldn't find working pool. Exit."<<endl;
43 fParent = wPool->IndexOf((TObject*)parent);
44 fChi2Dof = parent->Chi2Dof();
55 // Add itself to working pool
56 AddToWorkingPool((TObject*)this);
60 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Int_t i, Bool_t cpv) : AliPHOSCpvRecPoint()
66 AliPHOSEmcRecPoint* rp=0;
68 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
71 rp = (AliPHOSCpvRecPoint*)fGetter->CpvRecPoints()->At(i);
76 rp = (AliPHOSEmcRecPoint*)fGetter->EmcRecPoints()->At(i);
81 Int_t* Digits = rp->GetDigitsList();
82 Float_t* Energies = rp->GetEnergiesList();
83 Int_t nDigits = rp->GetMultiplicity();
85 for(Int_t iDigit=0; iDigit<nDigits; iDigit++) {
86 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
87 Float_t eDigit = Energies[iDigit];
88 this->AddDigit(*digit,eDigit);
92 rp->GetLocalPosition(locpos);
95 // Add itself to working pool
96 AddToWorkingPool((TObject*)this);
100 AliPHOSClusterizer* AliPHOSEvalRecPoint::GetClusterizer()
102 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
103 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
104 AliPHOSClusterizer* clu = (AliPHOSClusterizer*)wPoolF->FindObject("PHOS:clu-v1");
106 cout<<" Couldn't find Clusterizer. Exit."<<endl;
113 Bool_t AliPHOSEvalRecPoint::TooClose(AliPHOSRecPoint* pt) const
117 pt->GetLocalPosition(her_pos);
118 this->GetLocalPosition(my_pos);
119 Float_t dx = her_pos.X() - my_pos.X();
120 Float_t dz = her_pos.Z() - my_pos.Z();
121 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
123 if(dr<GetReconstructionManager()->MergeGammasMinDistanceCut())
130 Bool_t AliPHOSEvalRecPoint::NeedToSplit() const
135 void AliPHOSEvalRecPoint::DeleteParent()
140 void AliPHOSEvalRecPoint::UpdateWorkingPool()
143 for(Int_t i=0; i<this->InWorkingPool(); i++) {
144 AliPHOSEvalRecPoint* parent = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
146 Int_t nChild = parent->HasChild(children);
147 for(Int_t iChild=0; iChild<nChild; iChild++)
149 ((AliPHOSEvalRecPoint*)children.At(iChild))->DeleteParent();
153 RemoveFromWorkingPool(parent);
159 for(Int_t i=0; i<this->InWorkingPool(); i++) {
160 AliPHOSEvalRecPoint* weak = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
161 if (weak->KillWeakPoint()) delete weak;
164 for(Int_t i=0; i<this->InWorkingPool(); i++) {
165 AliPHOSEvalRecPoint* close = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
166 close->MergeClosePoint();
169 for(Int_t i=0; i<this->InWorkingPool(); i++)
170 ((AliPHOSEvalRecPoint*)AliPHOSEvalRecPoint::GetFromWorkingPool(i))->SetIndexInList(i);
174 Int_t AliPHOSEvalRecPoint::HasChild(TObjArray& children)
176 for(Int_t iChild=0; iChild<InWorkingPool(); iChild++)
178 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)GetFromWorkingPool(iChild);
179 TObject* parent = (TObject*)child->Parent();
180 TObject* me = (TObject*)this;
181 if(parent==me) children.Add(child);
184 return children.GetEntries();
187 void AliPHOSEvalRecPoint::Init()
190 AliPHOSClusterizer* clusterizer = GetClusterizer();
192 cout<<" Cannot get clusterizer. Exit."<<endl;
196 TClonesArray* digits = AliPHOSGetter::GetInstance()->Digits();
200 LogWeight = clusterizer->GetEmcLogWeight();
203 LogWeight = clusterizer->GetCpvLogWeight();
206 EvalLocalPosition(LogWeight,digits); // evaluate initial position
210 void AliPHOSEvalRecPoint::MakeJob()
212 // Reconstruction algoritm implementation.
214 AliPHOSRecManager* recMng = GetReconstructionManager();
221 Int_t nChild = HasChild(children);
225 if(Chi2Dof()>recMng->OneGamChisqCut()) SplitMergedShowers();
228 for(Int_t iChild=0; iChild<nChild; iChild++) {
229 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)children.At(iChild);
230 child->EvaluatePosition();
232 if(child->Chi2Dof()>recMng->OneGamChisqCut())
233 child->SplitMergedShowers();
238 void AliPHOSEvalRecPoint::InitTwoGam(Float_t* gamma1, Float_t* gamma2)
240 //Compute start values for two gamma fit algorithm.
241 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed "gammas".
244 TVector3 lpos; // start point choosing.
245 GetLocalPosition(lpos);
246 Float_t xx = lpos.Z();
247 Float_t yy = lpos.X();
248 Float_t E = GetEnergy();
250 cout<<" (x,z,E)[old] = ("<<yy<<","<<xx<<","<<E<<")"<<endl;
257 AliPHOSDigit * digit ;
258 Int_t nDigits = GetMultiplicity();
259 Int_t * Digits = GetDigitsList();
260 Float_t * Energies = GetEnergiesList();
272 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
273 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
275 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
277 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
278 eDigit = Energies[iDigit];
279 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
280 fGeom->RelPosInModule(relid, iy, ix);
282 Float_t dx = ix - xx;
283 Float_t dy = iy - yy;
290 sqr = TMath::Sqrt(4.*exy*exy + (exx-eyy)*(exx-eyy));
291 Float_t euu = (exx+eyy+sqr)/2.;
292 if(sqr>1.e-10) cos2fi = (exx - eyy)/sqr;
293 Float_t cosfi = TMath::Sqrt((1.+cos2fi)/2.);
294 Float_t sinfi = TMath::Sqrt((1.-cos2fi)/2.);
295 if(exy<0) sinfi = -sinfi;
298 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
300 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
301 eDigit = Energies[iDigit];
302 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
303 fGeom->RelPosInModule(relid, iy, ix);
305 Float_t dx = ix - xx;
306 Float_t dy = iy - yy;
307 Float_t du = dx*cosfi + dy*sinfi;
308 eu3 += eDigit*du*du*du;
311 Float_t c = E*eu3*eu3/(euu*euu*euu)/2.;
313 if(eu3<0) sign = -1.;
314 Float_t r = 1.+ c + sign*TMath::Sqrt(2.*c + c*c);
316 if(TMath::Abs(r-1.)>0.1) u = eu3/euu*(r+1.)/(r-1.);
317 if(TMath::Abs(r-1.)<0.1) u = TMath::Sqrt(sqr/E/r)*(1.+r);
319 Float_t e2c = E/(1.+r);
321 Float_t u1 = -u/(1.+r);
324 Float_t x1c = xx+u1*cosfi;
325 Float_t y1c = yy+u1*sinfi;
326 Float_t x2c = xx+u2*cosfi;
327 Float_t y2c = yy+u2*sinfi;
329 // printf("e1c -> %f\n",e1c);
330 // printf("x1c -> %f\n",x1c);
331 // printf("y1c -> %f\n",y1c);
332 // printf("e2c -> %f\n",e2c);
333 // printf("x2c -> %f\n",x2c);
334 // printf("y2c -> %f\n",y2c);
346 void AliPHOSEvalRecPoint::TwoGam(Float_t* gamma1, Float_t* gamma2)
348 //Fitting algorithm for the two very closed gammas
349 //that merged into the cluster with one maximum.
350 // Starting points gamma1 and gamma2 must be provided before ( see InitTwoGam)
351 //Set chisquare of the fit.
354 Float_t st0 = GetReconstructionManager()->TwoGamInitialStep();
355 Float_t chmin = GetReconstructionManager()->TwoGamChisqMin();
356 Float_t emin = GetReconstructionManager()->TwoGamEmin();
357 Float_t stpmin = GetReconstructionManager()->TwoGamStepMin();
358 Int_t Niter = GetReconstructionManager()->TwoGamNumOfIterations(); // Number of iterations.
360 Float_t chisq = 100.; //Initial chisquare.
362 Int_t nadc = GetMultiplicity();
366 Int_t dof = nadc - 5;
368 Float_t chstop = chmin*dof;
378 Float_t ee1=0,xx1=0,yy1=0,ee2=0,xx2=0,yy2=0,cosi;
380 Float_t e1c = gamma1[0];
381 Float_t y1c = gamma1[1];
382 Float_t x1c = gamma1[2];
384 Float_t e2c = gamma2[0];
385 Float_t y2c = gamma2[1];
386 Float_t x2c = gamma2[2];
388 Float_t E = GetEnergy();
391 AliPHOSDigit * digit ;
392 Int_t nDigits = GetMultiplicity();
393 Int_t * Digits = GetDigitsList();
394 Float_t * Energies = GetEnergiesList();
400 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
401 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
403 for(Int_t Iter=0; Iter<Niter; Iter++)
412 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
414 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
415 eDigit = Energies[iDigit];
416 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
417 fGeom->RelPosInModule(relid, iy, ix);
422 Float_t dx1 = x1c - ix;
423 Float_t dy1 = y1c - iy;
424 // cout<<" Mult "<<nDigits<<" dx1 "<<dx1<<" dy1 "<<dy1<<endl;
425 // AG(e1c,dx1,dy1,a1,gx1,gy1);
426 GetReconstructionManager()->AG(e1c,dx1,dy1,a1,gx1,gy1);
428 Float_t dx2 = x2c - ix;
429 Float_t dy2 = y2c - iy;
430 // cout<<" "<<" dx2 "<<dx2<<" dy2 "<<dy2<<endl;
431 // AG(e2c,dx2,dy2,a2,gx2,gy2);
432 GetReconstructionManager()->AG(e2c,dx2,dy2,a2,gx2,gy2);
435 // Float_t D = Const*A*(1. - A/E);
438 // // D = 9.; // ????
440 // Float_t da = A - eDigit;
441 // chisqc += da*da/D;
442 // Float_t dd = da/D;
443 // dd = dd*(2.-dd*Const*(1.-2.*A/E));
446 chisqc += GetReconstructionManager()->TwoGamChi2(A,eDigit,E,dd);
452 grec += (a1/e1c - a2/e2c)*E*dd;
456 Float_t grc = TMath::Sqrt(grx1c*grx1c + gry1c*gry1c + grx2c*grx2c + gry2c*gry2c);
457 if(grc<1.e-10) grc=1.e-10;
459 Float_t sc = 1. + chisqc/ch;
464 cosi = (grx1*grx1c + gry1*gry1c + grx2*grx2c + gry2*gry2c + gre*grec)/gr/grc;
465 st = st*sc/(1.4 - cosi);
487 Float_t step = st*gr;
490 cout<<" Iteration "<<Iter<<" dof "<<dof<<" chisq/dof "<<ch/dof<<" chstop/dof "<<chstop/dof<<" step " <<step<<" stpmin "<<stpmin<<endl;
496 Float_t de = st*gre*E;
500 if( (e1c>emin) && (e2c>emin) )
516 // if(ch>chisq*(nadc-2)-delch)
528 Float_t x1_new = yy1;
529 Float_t z1_new = xx1;
530 Float_t e1_new = ee1;
531 Float_t x2_new = yy2;
532 Float_t z2_new = xx2;
533 Float_t e2_new = ee2;
535 cout<<" (x,z,E)[1 fit] = ("<<x1_new<<","<<z1_new<<","<<e1_new<<")"<<endl;
536 cout<<" (x,z,E)[2 fit] = ("<<x2_new<<","<<z2_new<<","<<e2_new<<")"<<endl;
542 void AliPHOSEvalRecPoint::UnfoldTwoMergedPoints(Float_t* gamma1, Float_t* gamma2)
544 //Unfold TWO merged rec. points in the case when cluster has only one maximum,
545 //but it's fitting to the one gamma shower is too bad.
546 // Use TwoGam() to estimate the positions and energies of merged points.
552 Int_t* Digits = GetDigitsList();
553 Int_t nDigits = GetMultiplicity();
554 Float_t* Energies = GetEnergiesList();
555 Float_t* eFit = new Float_t[nDigits];
557 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
558 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
560 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
562 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
564 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
566 fGeom->RelPosInModule(relid, x, z);
569 for(Int_t iMax=0; iMax<Nmax; iMax++)
576 Float_t eMax = gamma[0];
577 Float_t xMax = gamma[1];
578 Float_t zMax = gamma[2];
580 Float_t dx = xMax - x;
581 Float_t dz = zMax - z;
583 GetReconstructionManager()->AG(eMax,dz,dx,Amp,gx,gy);
591 for( Int_t iMax=0; iMax<Nmax; iMax++)
598 Float_t eMax = gamma[0];
599 Float_t xMax = gamma[1];
600 Float_t zMax = gamma[2];
602 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
603 TVector3 newpos(xMax,0,zMax);
604 newRP->SetLocalPosition(newpos);
606 for( Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
608 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
609 Float_t eDigit = Energies[iDigit];
611 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
613 fGeom->RelPosInModule(relid, ix, iz);
615 Float_t dx = xMax - ix;
616 Float_t dz = zMax - iz;
617 Float_t single_shower_gain,gxMax,gyMax;
618 GetReconstructionManager()->AG(eMax,dz,dx,single_shower_gain,gxMax,gyMax);
619 Float_t total_gain = eFit[iDigit];
620 Float_t ratio = single_shower_gain/total_gain;
621 // cout<<" ratio -> "<<ratio<<endl;
622 eDigit = eDigit*ratio;
623 newRP->AddDigit(*digit,eDigit);
626 cout<<"======= Split: daughter rec point "<<iMax<<" ================="<<endl;
638 void AliPHOSEvalRecPoint::EvaluatePosition()
640 // One gamma fit algorithm.
641 // Set chisq/dof of the fit.
642 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed gamma, respectively.
644 Int_t nDigits = GetMultiplicity();
648 Int_t Niter = GetReconstructionManager()->OneGamNumOfIterations(); // number of iterations
649 Float_t St0 = GetReconstructionManager()->OneGamInitialStep(); // initial step
650 // const Float_t Stpmin = 0.005;
651 Float_t Stpmin = GetReconstructionManager()->OneGamStepMin();
652 Float_t chmin = GetReconstructionManager()->OneGamChisqMin();
660 GetLocalPosition(locpos);
661 Float_t E = GetEnergy();
662 Float_t xc = locpos.Z();
663 Float_t yc = locpos.X();
664 Float_t dx,dy,gx,gy,grxc,gryc;
666 Float_t chisq = 1.e+20;
670 Int_t dof = GetMultiplicity() - 2;
673 Float_t chstop = chmin*dof;
674 Float_t cosi,x1=0,y1=0;
677 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
678 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
680 for(Int_t Iter=0; Iter<Niter; Iter++)
686 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
689 Float_t* Energies = GetEnergiesList();
690 Int_t* Digits = GetDigitsList();
691 digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
692 eDigit = Energies[iDigit];
693 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
694 fGeom->RelPosInModule(relid, iy, ix);
703 GetReconstructionManager()->AG(E,dx,dy,A,gx,gy);
705 cout<<" (ix iy xc yc dx dy) "<<ix<<" "<<iy<<" "<<xc<<" "<<yc<<" "<<dx<<" "<<dy<<endl;
706 Float_t chi2dg = GetReconstructionManager()->OneGamChi2(A,eDigit,E,dd);
708 // Exclude digit with too large chisquare.
709 if(chi2dg > 10) { continue; }
717 Float_t grc = TMath::Sqrt(grxc*grxc + gryc*gryc);
718 if(grc<1.e-10) grc=1.e-10;
719 Float_t sc = 1. + chisqc/chisq;
720 cout<<" chisq: "<<chisq<<endl;
724 cosi = (grx*grxc + gry*gryc)/gr/grc;
725 st = st*sc/(1.4 - cosi);
737 Float_t step = st*gr;
740 cout<<" Iteration "<<Iter<<" dof "<<dof<<" chisq/dof "<<chisq/dof<<" chstop/dof "<<chstop/dof<<" step " <<step<<" Stpmin "<<Stpmin<<endl;
765 TVector3 newpos(gamma1[1],0,gamma1[2]);
766 //SetLocalPosition(newpos);
772 // RP->GetLocalPosition(pos);
773 // cout<<" (x,z)[old] = ("<<x_old<<","<<z_old<<")"<<endl;
774 // cout<<" (x,z)[new] = ("<<pos.X()<<","<<pos.Z()<<")"<<endl;
777 // cout<<" gamma[1]: "<<gamma1[1]<<" gamma1[2]: "<<gamma1[2]<<endl;
784 Bool_t AliPHOSEvalRecPoint::KillWeakPoint()
786 // Remove this point from further procession
787 // if it's energy is too small.
789 Float_t Thr0 = GetReconstructionManager()->KillGamMinEnergy();
791 if(GetEnergy()<Thr0) {
792 cout<<"+++++++ Killing this rec point ++++++++++"<<endl;
793 RemoveFromWorkingPool(this);
801 void AliPHOSEvalRecPoint::SplitMergedShowers()
803 // Split two very close rec. points.
805 if(GetMultiplicity()<2)
811 InitTwoGam(gamma1,gamma2); // start points for fit
812 TwoGam(gamma1,gamma2); // evaluate the positions and energies
813 UnfoldTwoMergedPoints(gamma1,gamma2); // make the job
818 void AliPHOSEvalRecPoint::MergeClosePoint()
821 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
822 // AliPHOSDigitizer* digitizer = fGetter->Digitizer();
823 // Float_t fPedestal = digitizer->GetPedestal(); // YVK 30.09.2001
824 // Float_t fSlope = digitizer->GetSlope();
826 for (Int_t i=0;i<InWorkingPool(); i++)
828 TObject* obj = GetFromWorkingPool(i);
829 if(!((TObject*)this)->IsEqual(obj))
831 AliPHOSRecPoint* rp = (AliPHOSRecPoint*)obj;
832 if(GetPHOSMod() == rp->GetPHOSMod())
836 cout<<"+++++++ Merging point 1: ++++++"<<endl;
838 cout<<"+++++++ and point 2: ++++++++++"<<endl;
839 ((AliPHOSEvalRecPoint*)rp)->Print("");
841 //merge two rec. points
844 this->GetLocalPosition(lpos1);
845 rp->GetLocalPosition(lpos2);
847 Int_t* Digits = rp->GetDigitsList();
848 Float_t dE = rp->GetEnergy()/(rp->GetEnergy()+this->GetEnergy());
849 Float_t new_x = lpos1.X()*dE + lpos2.X()*(1.-dE);
850 Float_t new_z = lpos1.Z()*dE + lpos2.Z()*(1.-dE);
851 Float_t new_E = rp->GetEnergy()+this->GetEnergy();
852 Float_t* Energies = ((AliPHOSEmcRecPoint*)rp)->GetEnergiesList();
854 for(Int_t iDigit=0; iDigit<rp->GetDigitsMultiplicity(); iDigit++)
856 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At(Digits[iDigit]);
857 Float_t eDigit = Energies[iDigit];
858 this->AddDigit(*digit,eDigit);
861 TVector3 newpos(new_x,0,new_z);
864 RemoveFromWorkingPool(rp);
867 cout<<"++++++ Resulting point: ++++++++"<<endl;
878 Int_t AliPHOSEvalRecPoint::UnfoldLocalMaxima()
880 // Make unfolding in the reconstruction point with several local maxima.
881 // Return the number of local maxima.
883 // if multiplicity less then 2 - nothing to unfold
884 if(GetMultiplicity()<2) return 1;
887 Float_t maxAtEnergy[1000];
888 Float_t LocMaxCut, LogWeight;
893 // AliPHOSClusterizer* clusterizer = fGetter->Clusterizer("PHOSclu-v1");
894 AliPHOSClusterizer* clusterizer = GetClusterizer();
896 cout<<" Cannot get clusterizer. Exit."<<endl;
901 LocMaxCut = clusterizer->GetEmcLocalMaxCut();
902 LogWeight = clusterizer->GetEmcLogWeight();
905 LocMaxCut = clusterizer->GetCpvLocalMaxCut();
906 LogWeight = clusterizer->GetCpvLogWeight();
909 AliPHOSGetter* fGetter = AliPHOSGetter::GetInstance();
910 const AliPHOSGeometry* fGeom = fGetter->PHOSGeometry();
911 TClonesArray* digits = fGetter->Digits();
913 // if number of local maxima less then 2 - nothing to unfold
914 Int_t Nmax = GetNumberOfLocalMax(maxAt,maxAtEnergy,LocMaxCut,digits);
915 if(Nmax<2) return Nmax;
917 Int_t* Digits = GetDigitsList();
918 Int_t nDigits = GetMultiplicity();
919 Float_t* Energies = GetEnergiesList();
920 Float_t* eFit = new Float_t[nDigits];
922 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
927 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
930 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
931 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
933 fGeom->RelPosInModule(relid, x, z);
935 for(Int_t iMax=0; iMax<Nmax; iMax++)
937 AliPHOSDigit* digitMax = (AliPHOSDigit*)maxAt[iMax];
938 Float_t eMax = maxAtEnergy[iMax];
939 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
940 fGeom->RelPosInModule(relid, xMax, zMax);
941 Float_t dx = xMax - x;
942 Float_t dz = zMax - z;
944 GetReconstructionManager()->AG(eMax,dz,dx,Amp,gx,gy);
951 for(Int_t iMax=0; iMax<Nmax; iMax++)
953 AliPHOSDigit* digitMax = (AliPHOSDigit*)maxAt[iMax];
954 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
955 fGeom->RelPosInModule(relid, xMax, zMax);
956 Float_t eMax = maxAtEnergy[iMax];
958 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
959 newRP->AddDigit(*digitMax,maxAtEnergy[iMax]);
961 //Neighbous ( matrix 3x3 around the local maximum)
962 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
964 AliPHOSDigit* digit = (AliPHOSDigit*)fGetter->Digits()->At( Digits[iDigit] );
965 Float_t eDigit = Energies[iDigit];
966 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
968 fGeom->RelPosInModule(relid, ix, iz);
970 if(AreNeighbours(digitMax,digit))
972 Float_t dx = xMax - ix;
973 Float_t dz = zMax - iz;
974 Float_t single_shower_gain,gxMax,gyMax;
975 GetReconstructionManager()->AG(eMax,dz,dx,single_shower_gain,gxMax,gyMax);
976 Float_t total_gain = eFit[iDigit];
977 Float_t ratio = single_shower_gain/total_gain;
978 cout<<" ratio -> "<<ratio<<endl;
979 eDigit = eDigit*ratio;
980 newRP->AddDigit(*digit,eDigit);
984 newRP->EvalLocalPosition(LogWeight,digits);
985 cout<<"======= Unfold: daughter rec point "<<iMax<<" ================="<<endl;
989 // RemoveFromWorkingPool(this);
998 void AliPHOSEvalRecPoint::PrintPoint(Option_t* opt)
1000 AliPHOSCpvRecPoint::Print(opt);
1003 GetLocalPosition(lpos);
1005 cout<<" Chi2/dof = "<<Chi2Dof()<<endl;
1006 cout<<" Local (x,z) = ("<<lpos.X()<<","<<lpos.Z()<<") in module "<<GetPHOSMod()<<endl;
1007 // if(GetReconstructionManager())
1008 // cout<<" Distance of merger = "<<GetReconstructionManager()->MergeGammasMinDistanceCut()<<endl;
1011 AliPHOSRecManager* AliPHOSEvalRecPoint::GetReconstructionManager() const
1014 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1015 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1016 AliPHOSRecManager* recMng = (AliPHOSRecManager*)wPoolF->FindObject("AliPHOSRecManager");
1018 cout<<" Couldn't find Reconstruction Manager. Exit."<<endl;
1026 AliPHOSRecPoint* AliPHOSEvalRecPoint::Parent()
1028 if(fParent<0) return NULL;
1030 return (AliPHOSRecPoint*)GetFromWorkingPool(fParent);
1034 Float_t AliPHOSEvalRecPoint::Chi2Dof() const
1039 const TObject* AliPHOSEvalRecPoint::GetWorkingPool()
1042 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1043 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1044 TObject* wPool = wPoolF->FindObject("SmartPoints");
1046 cout<<" Couldn't find Working Pool. Exit."<<endl;
1053 void AliPHOSEvalRecPoint::AddToWorkingPool(TObject* obj)
1055 ((TObjArray*)GetWorkingPool())->Add(obj);
1058 TObject* AliPHOSEvalRecPoint::GetFromWorkingPool(Int_t i)
1060 // return fWorkingPool->At(i);
1061 return ((TObjArray*)GetWorkingPool())->At(i);
1064 Int_t AliPHOSEvalRecPoint::InWorkingPool()
1066 return ((TObjArray*)GetWorkingPool())->GetEntries();
1069 void AliPHOSEvalRecPoint::RemoveFromWorkingPool(TObject* obj)
1071 ((TObjArray*)GetWorkingPool())->Remove(obj);
1072 ((TObjArray*)GetWorkingPool())->Compress();
1075 void AliPHOSEvalRecPoint::PrintWorkingPool()
1077 ((TObjArray*)GetWorkingPool())->Print("");