X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=HMPID%2FAliHMPIDRecon.cxx;h=ed3504abf6c7539634f9568849a7040423708a15;hp=3bbaf6d63f57c33b1202795db20f490b362d7f6c;hb=fe536aed16a013c43345ec639e0d4c4c6ddf65d4;hpb=c9e0bd24b76747e7aa81e41b528742cc683ff4e5 diff --git a/HMPID/AliHMPIDRecon.cxx b/HMPID/AliHMPIDRecon.cxx index 3bbaf6d63f5..ed3504abf6c 100644 --- a/HMPID/AliHMPIDRecon.cxx +++ b/HMPID/AliHMPIDRecon.cxx @@ -27,12 +27,14 @@ #include //HoughResponse() #include //CkovAngle() #include //CkovAngle() +#include //CkovAngle() //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ AliHMPIDRecon::AliHMPIDRecon(): - TTask("RichRec","RichPat"), + TNamed("RichRec","RichPat"), fPhotCnt(-1), fPhotFlag(0x0), + fPhotClusIndex(0x0), fPhotCkov(0x0), fPhotPhi(0x0), fPhotWei(0x0), @@ -62,6 +64,7 @@ void AliHMPIDRecon::InitVars(Int_t n) //.. if(n<=0) return; fPhotFlag = new Int_t[n]; + fPhotClusIndex = new Int_t[n]; fPhotCkov = new Double_t[n]; fPhotPhi = new Double_t[n]; fPhotWei = new Double_t[n]; @@ -74,19 +77,19 @@ void AliHMPIDRecon::DeleteVars()const //Delete variables //.. delete [] fPhotFlag; + delete [] fPhotClusIndex; delete [] fPhotCkov; delete [] fPhotPhi; delete [] fPhotWei; } //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ -void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t index,Double_t nmean) +void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t index,Double_t nmean,Float_t xRa,Float_t yRa) { // Pattern recognition method based on Hough transform // Arguments: pTrk - track for which Ckov angle is to be found // pCluLst - list of clusters for this chamber // Returns: - track ckov angle, [rad], - - + const Int_t nMinPhotAcc = 3; // Minimum number of photons required to perform the pattern recognition Int_t nClusTot = pCluLst->GetEntries(); @@ -95,8 +98,8 @@ void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t inde InitVars(nClusTot); - Float_t xRa,yRa,th,ph; - pTrk->GetHMPIDtrk(xRa,yRa,th,ph); //initialize this track: th and ph angles at middle of RAD + Float_t xPc,yPc,th,ph; + pTrk->GetHMPIDtrk(xPc,yPc,th,ph); //initialize this track: th and ph angles at middle of RAD SetTrack(xRa,yRa,th,ph); fParam->SetRefIdx(nmean); @@ -106,8 +109,11 @@ void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t inde fPhotCnt=0; + Int_t nPads = 0; + for (Int_t iClu=0; iCluGetEntriesFast();iClu++){//clusters loop - AliHMPIDCluster *pClu=(AliHMPIDCluster*)pCluLst->UncheckedAt(iClu); //get pointer to current cluster + AliHMPIDCluster *pClu=(AliHMPIDCluster*)pCluLst->UncheckedAt(iClu); //get pointer to current cluster + nPads+=pClu->Size(); if(iClu == index) { // this is the MIP! not a photon candidate: just store mip info mipX = pClu->X(); mipY = pClu->Y(); @@ -116,10 +122,12 @@ void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t inde continue; } chId=pClu->Ch(); + if(pClu->Q()>2*fParam->QCut()) continue; Double_t thetaCer,phiCer; if(FindPhotCkov(pClu->X(),pClu->Y(),thetaCer,phiCer)){ //find ckov angle for this photon candidate fPhotCkov[fPhotCnt]=thetaCer; //actual theta Cerenkov (in TRS) - fPhotPhi [fPhotCnt]=phiCer; //actual phi Cerenkov (in TRS): -pi to come back to "unusual" ref system (X,Y,-Z) + fPhotPhi [fPhotCnt]=phiCer; + fPhotClusIndex[fPhotCnt]=iClu; //actual phi Cerenkov (in TRS): -pi to come back to "unusual" ref system (X,Y,-Z) fPhotCnt++; //increment counter of photon candidates } }//clusters loop @@ -131,25 +139,27 @@ void AliHMPIDRecon::CkovAngle(AliESDtrack *pTrk,TClonesArray *pCluLst,Int_t inde pTrk->SetHMPIDsignal(kNoPhotAccept); //set the appropriate flag return; } - - + fMipPos.Set(mipX,mipY); - - + //PATTERN RECOGNITION STARTED: - Int_t iNrec=FlagPhot(HoughResponse()); //flag photons according to individual theta ckov with respect to most probable + Int_t iNrec=FlagPhot(HoughResponse(),pCluLst,pTrk); //flag photons according to individual theta ckov with respect to most probable + pTrk->SetHMPIDmip(mipX,mipY,mipQ,iNrec); //store mip info if(iNrec<1){ pTrk->SetHMPIDsignal(kNoPhotAccept); //no photon candidates are accepted return; } + + Int_t occupancy = (Int_t)(1000*(nPads/(6.*80.*48.))); + Double_t thetaC = FindRingCkov(pCluLst->GetEntries()); //find the best reconstructed theta Cherenkov // FindRingGeom(thetaC,2); - pTrk->SetHMPIDsignal(thetaC); //store theta Cherenkov - pTrk->SetHMPIDchi2(fCkovSigma2); //store errors squared - + pTrk->SetHMPIDsignal(thetaC+occupancy); //store theta Cherenkov + pTrk->SetHMPIDchi2(fCkovSigma2); //store errors squared + DeleteVars(); }//CkovAngle() //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ @@ -353,7 +363,7 @@ Double_t AliHMPIDRecon::FindRingCkov(Int_t) return weightThetaCerenkov; }//FindCkovRing() //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ -Int_t AliHMPIDRecon::FlagPhot(Double_t ckov) +Int_t AliHMPIDRecon::FlagPhot(Double_t ckov,TClonesArray *pCluLst, AliESDtrack *pTrk) { // Flag photon candidates if their individual ckov angle is inside the window around ckov angle returned by HoughResponse() // Arguments: ckov- value of most probable ckov angle for track as returned by HoughResponse() @@ -374,14 +384,27 @@ Int_t AliHMPIDRecon::FlagPhot(Double_t ckov) Int_t iInsideCnt = 0; //count photons which Theta ckov inside the window for(Int_t i=0;i= tmin && fPhotCkov[i] <= tmax) { - fPhotFlag[i]=2; + if(fPhotCkov[i] >= tmin && fPhotCkov[i] <= tmax) { + fPhotFlag[i]=2; + AddObjectToFriends(pCluLst,i,pTrk); iInsideCnt++; } } + return iInsideCnt; + }//FlagPhot() //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ +void AliHMPIDRecon::AddObjectToFriends(TClonesArray *pCluLst, Int_t photonIndex, AliESDtrack *pTrk) +{ +// Add AliHMPIDcluster object to ESD friends + + AliHMPIDCluster *pClu=(AliHMPIDCluster*)pCluLst->UncheckedAt(fPhotClusIndex[photonIndex]); + AliHMPIDCluster *pClus = new AliHMPIDCluster(*pClu); + pClus->SetChi2(fPhotCkov[photonIndex]); + pTrk->AddCalibObject(pClus); +} +//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ TVector2 AliHMPIDRecon::TracePhot(Double_t ckovThe,Double_t ckovPhi)const { // Trace a single Ckov photon from emission point somewhere in radiator up to photocathode taking into account ref indexes of materials it travereses @@ -464,6 +487,7 @@ Double_t AliHMPIDRecon::HoughResponse() Double_t sumPhots=phots->Integral(bin1,bin2); if(sumPhots<3) continue; // if less then 3 photons don't trust to this ring Double_t sumPhotsw=photsw->Integral(bin1,bin2); + if((Double_t)((i+0.5)*fDTheta)>0.7) continue; resultw->Fill((Double_t)((i+0.5)*fDTheta),sumPhotsw); } // evaluate the "BEST" theta ckov as the maximum value of histogramm @@ -474,3 +498,32 @@ Double_t AliHMPIDRecon::HoughResponse() return (Double_t)(locMax*fDTheta+0.5*fDTheta); //final most probable track theta ckov }//HoughResponse() //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ + Double_t AliHMPIDRecon::FindRingExt(Double_t ckov,Int_t ch,Double_t xPc,Double_t yPc,Double_t thRa,Double_t phRa) +{ +// To find the acceptance of the ring even from external inputs. +// +// + Double_t xRa = xPc - (fParam->RadThick()+fParam->WinThick()+fParam->GapThick())*TMath::Cos(phRa)*TMath::Tan(thRa); //just linear extrapolation back to RAD + Double_t yRa = yPc - (fParam->RadThick()+fParam->WinThick()+fParam->GapThick())*TMath::Sin(phRa)*TMath::Tan(thRa); + + Int_t nStep = 500; + Int_t nPhi = 0; + + Int_t ipc,ipadx,ipady; + + if(ckov>0){ + SetTrack(xRa,yRa,thRa,phRa); + for(Int_t j=0;jIsInDead(pos.X(),pos.Y())) continue; + fParam->Lors2Pad(pos.X(),pos.Y(),ipc,ipadx,ipady); + ipadx+=(ipc%2)*fParam->kPadPcX; + ipady+=(ipc/2)*fParam->kPadPcY; + if(ipadx<0 || ipady>160 || ipady<0 || ipady>144 || ch<0 || ch>6) continue; + if(fParam->IsDeadPad(ipadx,ipady,ch)) continue; + nPhi++; + }//point loop + return ((Double_t)nPhi/(Double_t)nStep); + }//if + return -1; +}