1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 //*-- Author: Boris Polichtchouk, IHEP
19 //////////////////////////////////////////////////////////////////////////////
20 // Reconstructed point operations for the Clusterization class for IHEP reconstruction.
21 // Performs clusterization (collects neighbouring active cells)
22 // It differs from AliPHOSClusterizerv1 in neighbour definition only
24 // ---- ROOT system ---
26 #include <TDirectory.h>
32 // --- AliRoot header files ---
34 #include "AliConfig.h"
35 #include "AliPHOSDigit.h"
36 #include "AliPHOSClusterizer.h"
37 #include "AliPHOSEvalRecPoint.h"
39 #include "AliPHOSLoader.h"
40 #include "AliPHOSRecCpvManager.h"
41 #include "AliPHOSRecEmcManager.h"
42 #include "AliPHOSDigitizer.h"
43 #include "AliPHOSGeometry.h"
45 // --- Standard library ---
48 ClassImp(AliPHOSEvalRecPoint)
50 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint() :
55 fEventFolderName(AliConfig::GetDefaultEventFolderName())
60 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Bool_t cpv, AliPHOSEvalRecPoint* parent) :
72 TObjArray* wPool = (TObjArray*)GetWorkingPool();
74 Fatal("ctor", "Couldn't find working pool") ;
77 fParent = wPool->IndexOf((TObject*)parent);
78 fChi2Dof = parent->Chi2Dof();
89 // Add itself to working pool
90 AddToWorkingPool((TObject*)this);
94 AliPHOSEvalRecPoint::AliPHOSEvalRecPoint(Int_t i, Bool_t cpv) :
99 fEventFolderName(AliConfig::GetDefaultEventFolderName())
102 AliPHOSEmcRecPoint* rp=0;
104 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
107 rp = (AliPHOSCpvRecPoint*)fLoader->CpvRecPoints()->At(i);
112 rp = (AliPHOSEmcRecPoint*)fLoader->EmcRecPoints()->At(i);
117 Int_t* digits = rp->GetDigitsList();
118 Float_t* energies = rp->GetEnergiesList();
119 Int_t nDigits = rp->GetMultiplicity();
121 for(Int_t iDigit=0; iDigit<nDigits; iDigit++) {
122 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
123 Float_t eDigit = energies[iDigit];
124 this->AddDigit(*digit,eDigit);
128 rp->GetLocalPosition(locpos);
131 // Add itself to working pool
132 AddToWorkingPool((TObject*)this);
136 AliPHOSClusterizer* AliPHOSEvalRecPoint::GetClusterizer()
138 // returns clusterizer task
139 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
140 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
141 AliPHOSClusterizer* clu = (AliPHOSClusterizer*)wPoolF->FindObject("PHOS:clu-v1");
143 Fatal("GetClusterizer", "Couldn't find Clusterizer") ;
149 Bool_t AliPHOSEvalRecPoint::TooClose(AliPHOSRecPoint* pt) const
151 // check if a rec.point pt is too close to this one
154 pt->GetLocalPosition(herPos);
155 this->GetLocalPosition(myPos);
156 Float_t dx = herPos.X() - myPos.X();
157 Float_t dz = herPos.Z() - myPos.Z();
158 Float_t dr = TMath::Sqrt(dx*dx + dz*dz);
160 if(dr<GetReconstructionManager()->MergeGammasMinDistanceCut())
167 Bool_t AliPHOSEvalRecPoint::NeedToSplit() const
169 // rec.point needs to be split
173 void AliPHOSEvalRecPoint::DeleteParent()
179 void AliPHOSEvalRecPoint::UpdateWorkingPool()
181 // update pool of rec.points
182 Int_t i; //loop variable
184 for(i=0; i<this->InWorkingPool(); i++) {
185 AliPHOSEvalRecPoint* parent = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
187 Int_t nChild = parent->HasChild(children);
188 for(Int_t iChild=0; iChild<nChild; iChild++)
190 ((AliPHOSEvalRecPoint*)children.At(iChild))->DeleteParent();
194 RemoveFromWorkingPool(parent);
200 for(i=0; i<this->InWorkingPool(); i++) {
201 AliPHOSEvalRecPoint* weak = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
202 if (weak->KillWeakPoint()) delete weak;
205 for(i=0; i<this->InWorkingPool(); i++) {
206 AliPHOSEvalRecPoint* close = (AliPHOSEvalRecPoint*)GetFromWorkingPool(i);
207 close->MergeClosePoint();
210 for(i=0; i<this->InWorkingPool(); i++)
211 ((AliPHOSEvalRecPoint*)AliPHOSEvalRecPoint::GetFromWorkingPool(i))->SetIndexInList(i);
215 Int_t AliPHOSEvalRecPoint::HasChild(TObjArray& children)
217 // returns the number of children
218 for(Int_t iChild=0; iChild<InWorkingPool(); iChild++)
220 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)GetFromWorkingPool(iChild);
221 TObject* parent = (TObject*)child->Parent();
222 TObject* me = (TObject*)this;
223 if(parent==me) children.Add(child);
226 return children.GetEntries();
229 void AliPHOSEvalRecPoint::Init()
232 AliPHOSClusterizer* clusterizer = GetClusterizer();
234 Fatal("Init", "Cannot get clusterizer") ;
237 TClonesArray* digits = AliPHOSLoader::GetPHOSLoader(fEventFolderName)->Digits();
241 logWeight = clusterizer->GetEmcLogWeight();
244 logWeight = clusterizer->GetCpvLogWeight();
247 TVector3 vtx(0.,0.,0.) ;
248 EvalLocalPosition(logWeight,vtx,digits); // evaluate initial position
252 void AliPHOSEvalRecPoint::MakeJob()
254 // Reconstruction algoritm implementation.
256 AliPHOSRecManager* recMng = GetReconstructionManager();
263 Int_t nChild = HasChild(children);
267 if(Chi2Dof()>recMng->OneGamChisqCut()) SplitMergedShowers();
270 for(Int_t iChild=0; iChild<nChild; iChild++) {
271 AliPHOSEvalRecPoint* child = (AliPHOSEvalRecPoint*)children.At(iChild);
272 child->EvaluatePosition();
274 if(child->Chi2Dof()>recMng->OneGamChisqCut())
275 child->SplitMergedShowers();
280 void AliPHOSEvalRecPoint::InitTwoGam(Float_t* gamma1, Float_t* gamma2)
282 //Compute start values for two gamma fit algorithm.
283 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed "gammas".
286 TVector3 lpos; // start point choosing.
287 GetLocalPosition(lpos);
288 Float_t xx = lpos.Z();
289 Float_t yy = lpos.X();
290 Float_t e = GetEnergy();
292 AliInfo(Form("(x,z,e)[old] = (%f, %f, %f)", yy, xx, e)) ;
299 AliPHOSDigit * digit ;
300 Int_t nDigits = GetMultiplicity();
301 Int_t * digits = GetDigitsList();
302 Float_t * energies = GetEnergiesList();
314 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
315 const AliPHOSGeometry* fGeom = fLoader->PHOSGeometry();
317 Int_t iDigit; //loop variable
319 for(iDigit = 0 ; iDigit < nDigits ; iDigit ++)
321 digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
322 eDigit = energies[iDigit];
323 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
324 fGeom->RelPosInModule(relid, iy, ix);
326 Float_t dx = ix - xx;
327 Float_t dy = iy - yy;
334 sqr = TMath::Sqrt(4.*exy*exy + (exx-eyy)*(exx-eyy));
335 Float_t euu = (exx+eyy+sqr)/2.;
336 if(sqr>1.e-10) cos2fi = (exx - eyy)/sqr;
337 Float_t cosfi = TMath::Sqrt((1.+cos2fi)/2.);
338 Float_t sinfi = TMath::Sqrt((1.-cos2fi)/2.);
339 if(exy<0) sinfi = -sinfi;
342 for(iDigit = 0 ; iDigit < nDigits ; iDigit ++)
344 digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
345 eDigit = energies[iDigit];
346 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
347 fGeom->RelPosInModule(relid, iy, ix);
349 Float_t dx = ix - xx;
350 Float_t dy = iy - yy;
351 Float_t du = dx*cosfi + dy*sinfi;
352 eu3 += eDigit*du*du*du;
355 Float_t c = e*eu3*eu3/(euu*euu*euu)/2.;
357 if(eu3<0) sign = -1.;
358 Float_t r = 1.+ c + sign*TMath::Sqrt(2.*c + c*c);
360 if(TMath::Abs(r-1.)>0.1) u = eu3/euu*(r+1.)/(r-1.);
361 if(TMath::Abs(r-1.)<0.1) u = TMath::Sqrt(sqr/e/r)*(1.+r);
363 Float_t e2c = e/(1.+r);
365 Float_t u1 = -u/(1.+r);
368 Float_t x1c = xx+u1*cosfi;
369 Float_t y1c = yy+u1*sinfi;
370 Float_t x2c = xx+u2*cosfi;
371 Float_t y2c = yy+u2*sinfi;
373 // printf("e1c -> %f\n",e1c);
374 // printf("x1c -> %f\n",x1c);
375 // printf("y1c -> %f\n",y1c);
376 // printf("e2c -> %f\n",e2c);
377 // printf("x2c -> %f\n",x2c);
378 // printf("y2c -> %f\n",y2c);
390 void AliPHOSEvalRecPoint::TwoGam(Float_t* gamma1, Float_t* gamma2)
392 //Fitting algorithm for the two very closed gammas
393 //that merged into the cluster with one maximum.
394 // Starting points gamma1 and gamma2 must be provided before ( see InitTwoGam)
395 //Set chisquare of the fit.
398 Float_t st0 = GetReconstructionManager()->TwoGamInitialStep();
399 Float_t chmin = GetReconstructionManager()->TwoGamChisqMin();
400 Float_t emin = GetReconstructionManager()->TwoGamEmin();
401 Float_t stpmin = GetReconstructionManager()->TwoGamStepMin();
402 Int_t nIter = GetReconstructionManager()->TwoGamNumOfIterations(); // Number of iterations.
404 Float_t chisq = 100.; //Initial chisquare.
406 Int_t nadc = GetMultiplicity();
410 Int_t dof = nadc - 5;
412 Float_t chstop = chmin*dof;
422 Float_t ee1=0,xx1=0,yy1=0,ee2=0,xx2=0,yy2=0,cosi;
424 Float_t e1c = gamma1[0];
425 Float_t y1c = gamma1[1];
426 Float_t x1c = gamma1[2];
428 Float_t e2c = gamma2[0];
429 Float_t y2c = gamma2[1];
430 Float_t x2c = gamma2[2];
432 Float_t e = GetEnergy();
435 AliPHOSDigit * digit ;
436 Int_t nDigits = GetMultiplicity();
437 Int_t * digits = GetDigitsList();
438 Float_t * energies = GetEnergiesList();
444 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
445 const AliPHOSGeometry* fGeom = fLoader->PHOSGeometry();
447 for(Int_t iter=0; iter<nIter; iter++)
456 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
458 digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
459 eDigit = energies[iDigit];
460 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
461 fGeom->RelPosInModule(relid, iy, ix);
466 Float_t dx1 = x1c - ix;
467 Float_t dy1 = y1c - iy;
469 GetReconstructionManager()->AG(e1c,dx1,dy1,a1,gx1,gy1);
471 Float_t dx2 = x2c - ix;
472 Float_t dy2 = y2c - iy;
474 GetReconstructionManager()->AG(e2c,dx2,dy2,a2,gx2,gy2);
477 // Float_t D = Const*a*(1. - a/e);
480 // // D = 9.; // ????
482 // Float_t da = a - eDigit;
483 // chisqc += da*da/D;
484 // Float_t dd = da/D;
485 // dd = dd*(2.-dd*Const*(1.-2.*a/e));
488 chisqc += GetReconstructionManager()->TwoGamChi2(a,eDigit,e,dd);
494 grec += (a1/e1c - a2/e2c)*e*dd;
498 Float_t grc = TMath::Sqrt(grx1c*grx1c + gry1c*gry1c + grx2c*grx2c + gry2c*gry2c);
499 if(grc<1.e-10) grc=1.e-10;
501 Float_t sc = 1. + chisqc/ch;
506 cosi = (grx1*grx1c + gry1*gry1c + grx2*grx2c + gry2*gry2c + gre*grec)/gr/grc;
507 st = st*sc/(1.4 - cosi);
529 Float_t step = st*gr;
531 AliInfo(Form("Iteration %d dof %d chisq/dof %f chstop/dof %f step %d stpmin %d",
532 iter, dof, ch/dof, chstop/dof, step, stpmin)) ;
538 Float_t de = st*gre*e;
542 if( (e1c>emin) && (e2c>emin) )
558 // if(ch>chisq*(nadc-2)-delch)
578 message = " (x,z,e)[1 fit] = (%f, %f, %f)\n" ;
579 message = " (x,z,e)[2 fit] = (%f, %f, %f)\n" ;
580 AliInfo(Form(message.Data(),
582 x2New, z2New, e2New)) ;
588 void AliPHOSEvalRecPoint::UnfoldTwoMergedPoints(Float_t* gamma1, Float_t* gamma2)
590 //Unfold TWO merged rec. points in the case when cluster has only one maximum,
591 //but it's fitting to the one gamma shower is too bad.
592 // Use TwoGam() to estimate the positions and energies of merged points.
598 Int_t* digits = GetDigitsList();
599 Int_t nDigits = GetMultiplicity();
600 Float_t* energies = GetEnergiesList();
601 Float_t* eFit = new Float_t[nDigits];
603 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
604 const AliPHOSGeometry* fGeom = fLoader->PHOSGeometry();
606 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
608 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
610 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
612 fGeom->RelPosInModule(relid, x, z);
615 for(Int_t iMax=0; iMax<nMax; iMax++)
622 Float_t eMax = gamma[0];
623 Float_t xMax = gamma[1];
624 Float_t zMax = gamma[2];
626 Float_t dx = xMax - x;
627 Float_t dz = zMax - z;
629 GetReconstructionManager()->AG(eMax,dz,dx,amp,gx,gy);
637 for( Int_t iMax=0; iMax<nMax; iMax++)
644 Float_t eMax = gamma[0];
645 Float_t xMax = gamma[1];
646 Float_t zMax = gamma[2];
648 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
649 TVector3 newpos(xMax,0,zMax);
650 newRP->SetLocalPosition(newpos);
652 for( Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
654 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
655 Float_t eDigit = energies[iDigit];
657 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
659 fGeom->RelPosInModule(relid, ix, iz);
661 Float_t dx = xMax - ix;
662 Float_t dz = zMax - iz;
663 Float_t singleShowerGain,gxMax,gyMax;
664 GetReconstructionManager()->AG(eMax,dz,dx,singleShowerGain,gxMax,gyMax);
665 Float_t totalGain = eFit[iDigit];
666 Float_t ratio = singleShowerGain/totalGain;
667 eDigit = eDigit*ratio;
668 newRP->AddDigit(*digit,eDigit);
670 AliInfo(Form("======= Split: daughter rec point %d =================",
682 void AliPHOSEvalRecPoint::EvaluatePosition()
684 // One gamma fit algorithm.
685 // Set chisq/dof of the fit.
686 // gamma1[0], gamma1[1], gamma1[2] are Energy,x,z of the reconstructed gamma, respectively.
688 Int_t nDigits = GetMultiplicity();
692 Int_t nIter = GetReconstructionManager()->OneGamNumOfIterations(); // number of iterations
693 Float_t st0 = GetReconstructionManager()->OneGamInitialStep(); // initial step
694 // const Float_t stpMin = 0.005;
695 Float_t stpMin = GetReconstructionManager()->OneGamStepMin();
696 Float_t chmin = GetReconstructionManager()->OneGamChisqMin();
704 GetLocalPosition(locpos);
705 Float_t e = GetEnergy();
706 Float_t xc = locpos.Z();
707 Float_t yc = locpos.X();
708 Float_t dx,dy,gx,gy,grxc,gryc;
710 Float_t chisq = 1.e+20;
714 Int_t dof = GetMultiplicity() - 2;
717 Float_t chstop = chmin*dof;
718 Float_t cosi,x1=0,y1=0;
721 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
722 const AliPHOSGeometry* fGeom = fLoader->PHOSGeometry();
724 for(Int_t iter=0; iter<nIter; iter++)
730 for(Int_t iDigit = 0 ; iDigit < nDigits ; iDigit ++)
733 Float_t* energies = GetEnergiesList();
734 Int_t* digits = GetDigitsList();
735 digit = (AliPHOSDigit*)fLoader->Digits()->At( digits[iDigit] );
736 eDigit = energies[iDigit];
737 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
738 fGeom->RelPosInModule(relid, iy, ix);
747 GetReconstructionManager()->AG(e,dx,dy,a,gx,gy);
749 AliInfo(Form(" (ix iy xc yc dx dy) %f %f %f %f %f %f",
750 ix, iy, xc, yc, dx, dy)) ;
751 Float_t chi2dg = GetReconstructionManager()->OneGamChi2(a,eDigit,e,dd);
753 // Exclude digit with too large chisquare.
754 if(chi2dg > 10) { continue; }
762 Float_t grc = TMath::Sqrt(grxc*grxc + gryc*gryc);
763 if(grc<1.e-10) grc=1.e-10;
764 Float_t sc = 1. + chisqc/chisq;
765 AliInfo(Form(" chisq: %f", chisq)) ;
769 cosi = (grx*grxc + gry*gryc)/gr/grc;
770 st = st*sc/(1.4 - cosi);
782 Float_t step = st*gr;
784 AliInfo(Form(" Iteration %d dof %d chisq/dof %f chstop/dof %f step %d stpMin %d",
785 iter, dof, chisq/dof, chisq/dof, chstop/dof, step, stpMin)) ;
810 TVector3 newpos(gamma1[1],0,gamma1[2]);
811 //SetLocalPosition(newpos);
817 // RP->GetLocalPosition(pos);
826 Bool_t AliPHOSEvalRecPoint::KillWeakPoint()
828 // Remove this point from further procession
829 // if it's energy is too small.
831 Float_t thr0 = GetReconstructionManager()->KillGamMinEnergy();
833 if(GetEnergy()<thr0) {
834 AliInfo(Form("+++++++ Killing this rec point ++++++++++")) ;
835 RemoveFromWorkingPool(this);
843 void AliPHOSEvalRecPoint::SplitMergedShowers()
845 // Split two very close rec. points.
847 if(GetMultiplicity()<2)
853 InitTwoGam(gamma1,gamma2); // start points for fit
854 TwoGam(gamma1,gamma2); // evaluate the positions and energies
855 UnfoldTwoMergedPoints(gamma1,gamma2); // make the job
860 void AliPHOSEvalRecPoint::MergeClosePoint()
862 // merge rec.points if they are too close
863 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
864 // AliPHOSDigitizer* digitizer = fLoader->Digitizer();
865 // Float_t fPedestal = digitizer->GetPedestal(); // YVK 30.09.2001
866 // Float_t fSlope = digitizer->GetSlope();
868 for (Int_t i=0;i<InWorkingPool(); i++)
870 TObject* obj = GetFromWorkingPool(i);
871 if(!((TObject*)this)->IsEqual(obj))
873 AliPHOSRecPoint* rp = (AliPHOSRecPoint*)obj;
874 if(GetPHOSMod() == rp->GetPHOSMod())
878 AliInfo(Form("+++++++ Merging point 1: ++++++")) ;
880 AliInfo(Form("+++++++ and point 2: ++++++++++")) ;
881 ((AliPHOSEvalRecPoint*)rp)->Print();
883 //merge two rec. points
886 this->GetLocalPosition(lpos1);
887 rp->GetLocalPosition(lpos2);
889 Int_t* digits = rp->GetDigitsList();
890 Float_t dE = rp->GetEnergy()/(rp->GetEnergy()+this->GetEnergy());
891 Float_t newX = lpos1.X()*dE + lpos2.X()*(1.-dE);
892 Float_t newZ = lpos1.Z()*dE + lpos2.Z()*(1.-dE);
893 Float_t newE = rp->GetEnergy()+this->GetEnergy();
894 Float_t* energies = ((AliPHOSEmcRecPoint*)rp)->GetEnergiesList();
896 for(Int_t iDigit=0; iDigit<rp->GetDigitsMultiplicity(); iDigit++)
898 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At(digits[iDigit]);
899 Float_t eDigit = energies[iDigit];
900 this->AddDigit(*digit,eDigit);
903 TVector3 newpos(newX,0,newZ);
906 RemoveFromWorkingPool(rp);
909 AliInfo(Form("++++++ Resulting point: ++++++++")) ;
920 Int_t AliPHOSEvalRecPoint::UnfoldLocalMaxima()
922 // Make unfolding in the reconstruction point with several local maxima.
923 // Return the number of local maxima.
925 // if multiplicity less then 2 - nothing to unfold
926 if(GetMultiplicity()<2) return 1;
928 AliPHOSDigit * maxAt[1000];
929 Float_t maxAtEnergy[1000];
930 Float_t locMaxCut, logWeight;
935 AliPHOSClusterizer* clusterizer = GetClusterizer();
937 AliFatal(Form("Cannot get clusterizer")) ;
941 locMaxCut = clusterizer->GetEmcLocalMaxCut();
942 logWeight = clusterizer->GetEmcLogWeight();
945 locMaxCut = clusterizer->GetCpvLocalMaxCut();
946 logWeight = clusterizer->GetCpvLogWeight();
949 AliPHOSLoader* fLoader = AliPHOSLoader::GetPHOSLoader(fEventFolderName);
950 const AliPHOSGeometry* fGeom = fLoader->PHOSGeometry();
951 TClonesArray* digits = fLoader->Digits();
953 // if number of local maxima less then 2 - nothing to unfold
954 Int_t nMax = GetNumberOfLocalMax(maxAt,maxAtEnergy,locMaxCut,digits);
955 if(nMax<2) return nMax;
957 Int_t* digitsList = GetDigitsList();
958 Int_t nDigits = GetMultiplicity();
959 Float_t* energies = GetEnergiesList();
960 Float_t* eFit = new Float_t[nDigits];
962 Int_t iDigit; //loop variable
964 for(iDigit=0; iDigit<nDigits; iDigit++)
969 for(iDigit=0; iDigit<nDigits; iDigit++)
972 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At( digitsList[iDigit] );
973 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
975 fGeom->RelPosInModule(relid, x, z);
977 for(Int_t iMax=0; iMax<nMax; iMax++)
979 AliPHOSDigit* digitMax = maxAt[iMax];
980 Float_t eMax = maxAtEnergy[iMax];
981 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
982 fGeom->RelPosInModule(relid, xMax, zMax);
983 Float_t dx = xMax - x;
984 Float_t dz = zMax - z;
986 GetReconstructionManager()->AG(eMax,dz,dx,amp,gx,gy);
993 for(Int_t iMax=0; iMax<nMax; iMax++)
995 AliPHOSDigit* digitMax = maxAt[iMax];
996 fGeom->AbsToRelNumbering(digitMax->GetId(), relid) ;
997 fGeom->RelPosInModule(relid, xMax, zMax);
998 Float_t eMax = maxAtEnergy[iMax];
1000 AliPHOSEvalRecPoint* newRP = new AliPHOSEvalRecPoint(IsCPV(),this);
1001 newRP->AddDigit(*digitMax,maxAtEnergy[iMax]);
1003 //Neighbous ( matrix 3x3 around the local maximum)
1004 for(Int_t iDigit=0; iDigit<nDigits; iDigit++)
1006 AliPHOSDigit* digit = (AliPHOSDigit*)fLoader->Digits()->At( digitsList[iDigit] );
1007 Float_t eDigit = energies[iDigit];
1008 fGeom->AbsToRelNumbering(digit->GetId(), relid) ;
1010 fGeom->RelPosInModule(relid, ix, iz);
1012 if(AreNeighbours(digitMax,digit))
1014 Float_t dx = xMax - ix;
1015 Float_t dz = zMax - iz;
1016 Float_t singleShowerGain,gxMax,gyMax;
1017 GetReconstructionManager()->AG(eMax,dz,dx,singleShowerGain,gxMax,gyMax);
1018 Float_t totalGain = eFit[iDigit];
1019 Float_t ratio = singleShowerGain/totalGain;
1020 AliInfo(Form(" ratio -> %f", ratio)) ;
1021 eDigit = eDigit*ratio;
1022 newRP->AddDigit(*digit,eDigit);
1026 TVector3 vtx(0.,0.,0.) ;
1027 newRP->EvalLocalPosition(logWeight,vtx,digits);
1028 AliInfo(Form("======= Unfold: daughter rec point %d =================",
1033 // RemoveFromWorkingPool(this);
1041 void AliPHOSEvalRecPoint::PrintPoint()
1043 // print rec.point to stdout
1044 AliPHOSCpvRecPoint::Print();
1047 GetLocalPosition(lpos);
1050 message = " Chi2/dof = %f" ;
1051 message += " Local (x,z) = (%f, %f) in module %d" ;
1052 AliInfo(Form(message.Data(), Chi2Dof(), lpos.X(), lpos.Z(), GetPHOSMod())) ;
1056 AliPHOSRecManager* AliPHOSEvalRecPoint::GetReconstructionManager() const
1058 // returns reconstruction manager
1059 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1060 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1061 AliPHOSRecManager* recMng = (AliPHOSRecManager*)wPoolF->FindObject("AliPHOSRecManager");
1063 AliFatal(Form("Couldn't find Reconstruction Manager")) ;
1070 AliPHOSRecPoint* AliPHOSEvalRecPoint::Parent()
1073 if(fParent<0) return NULL;
1075 return (AliPHOSRecPoint*)GetFromWorkingPool(fParent);
1079 Float_t AliPHOSEvalRecPoint::Chi2Dof() const
1081 // returns a chi^2 per degree of freedom
1085 const TObject* AliPHOSEvalRecPoint::GetWorkingPool()
1087 // returns a pool of rec.points
1088 TFolder* aliceF = (TFolder*)gROOT->FindObjectAny("YSAlice");
1089 TFolder* wPoolF = (TFolder*)aliceF->FindObject("WhiteBoard/RecPoints/PHOS/SmP");
1090 TObject* wPool = wPoolF->FindObject("SmartPoints");
1092 AliFatal(Form("Couldn't find Working Pool")) ;
1098 void AliPHOSEvalRecPoint::AddToWorkingPool(TObject* obj)
1100 // add a rec.point to a pool
1101 ((TObjArray*)GetWorkingPool())->Add(obj);
1104 TObject* AliPHOSEvalRecPoint::GetFromWorkingPool(Int_t i)
1106 // return a rec.point from a pool at an index i
1107 // return fWorkingPool->At(i);
1108 return ((TObjArray*)GetWorkingPool())->At(i);
1111 Int_t AliPHOSEvalRecPoint::InWorkingPool()
1113 // return the number of rec.points in a pool
1114 return ((TObjArray*)GetWorkingPool())->GetEntries();
1117 void AliPHOSEvalRecPoint::RemoveFromWorkingPool(TObject* obj)
1119 // remove a rec.point from a pool
1120 ((TObjArray*)GetWorkingPool())->Remove(obj);
1121 ((TObjArray*)GetWorkingPool())->Compress();
1124 void AliPHOSEvalRecPoint::PrintWorkingPool()
1126 // print pool of rec.points to stdout
1127 ((TObjArray*)GetWorkingPool())->Print();