// 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();
if(nClusTot>fParam->MultCut()) fIsWEIGHT = kTRUE; // offset to take into account bkg in reconstruction
else fIsWEIGHT = kFALSE;
}
}
}//clusters loop
+
fMipPos.Set(mipX,mipY);
- if(fPhotCnt<=3) pTrk->SetHMPIDsignal(kNoPhotAccept); //no reconstruction with <=3 photon candidates
+
+ if(fPhotCnt<=nMinPhotAcc) { //no reconstruction with <=3 photon candidates
+ pTrk->SetHMPIDsignal(kNoPhotAccept); //set the appropriate flag
+ pTrk->SetHMPIDmip(mipX,mipY,mipQ,fPhotCnt); //store mip info
+ return;
+ }
+
+ if(mipId==-1) {pTrk->SetHMPIDsignal(kMipQdcCut); return;} //no clusters with QDC more the threshold at all
+ if(dMin>fParam->DistCut()) {pTrk->SetHMPIDsignal(kMipDistCut); return;} //closest cluster with enough charge is still too far from intersection
+
+//PATTERN RECOGNITION STARTED:
+
Int_t iNrec=FlagPhot(HoughResponse()); //flag photons according to individual theta ckov with respect to most probable
pTrk->SetHMPIDmip(mipX,mipY,mipQ,iNrec); //store mip info
- if(mipId==-1) {pTrk->SetHMPIDsignal(kMipQdcCut); return;} //no clusters with QDC more the threshold at all
- if(dMin>fParam->DistCut()) {pTrk->SetHMPIDsignal(kMipDistCut); return;} //closest cluster with enough charge is still too far from intersection
pTrk->SetHMPIDcluIdx(chId,mipId); //set index of cluster
+
if(iNrec<1){
pTrk->SetHMPIDsignal(kNoPhotAccept); //no photon candidates are accepted
}
{fPc.Set(xPc,yPc);} //set track impact to PC
void SetMip (Double_t xmip,Double_t ymip )
{fMipPos.Set(xmip,ymip);} //set track impact to PC
- enum ETrackingFlags {kMipDistCut=-9,kMipQdcCut=-5,kNoPhotAccept=-11};
+ enum ETrackingFlags {kNotPerformed=-20,kMipDistCut=-9,kMipQdcCut=-5,kNoPhotAccept=-11};
//
protected:
Int_t fPhotCnt; // counter of photons candidate
return kTRUE;
}//GetTrackPoint()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Int_t AliHMPIDTracker::IntTrkCha(AliESDtrack *pTrk,Float_t &xPc,Float_t &yPc)
+Int_t AliHMPIDTracker::IntTrkCha(AliESDtrack *pTrk,Float_t &xPc,Float_t &yPc,Float_t &xRa,Float_t &yRa,Float_t &theta,Float_t &phi)
{
// Static method to find intersection in between given track and HMPID chambers
// Arguments: pTrk- ESD track; xPc,yPc- track intersection with PC in LORS [cm]
// Returns: intersected chamber ID or -1
AliHMPIDParam *pParam=AliHMPIDParam::Instance();
- Float_t xRa=0,yRa=0,theta=0,phi=0; //track intersection at PC and angles at RAD, LORS
for(Int_t i=AliHMPIDParam::kMinCh;i<=AliHMPIDParam::kMaxCh;i++){ //chambers loop
Double_t p1[3],n1[3]; pParam->Norm(i,n1); pParam->Point(i,p1,AliHMPIDParam::kRad); //point & norm for middle of radiator plane
Double_t p2[3],n2[3]; pParam->Norm(i,n2); pParam->Point(i,p2,AliHMPIDParam::kPc); //point & norm for entrance to PC plane
pParam->Mars2LorsVec(i,n1,theta,phi); //track angles at RAD
pParam->Mars2Lors (i,p1,xRa,yRa); //TRKxRAD position
pParam->Mars2Lors (i,p2,xPc,yPc); //TRKxPC position
- if(AliHMPIDParam::IsInside(xPc,yPc,pParam->DistCut())==kFALSE) continue; //not in active area
- pTrk->SetHMPIDtrk (xRa,yRa,theta,phi); //store track intersection info
- pTrk->SetHMPIDcluIdx (i,0);
- return i;
+ if(AliHMPIDParam::IsInside(xPc,yPc,pParam->DistCut())==kTRUE) return i; //return intersected chamber
} //chambers loop
- return -1; //no intersection with HMPID chambers
+ return -1; //no intersection with HMPID chambers
}//IntTrkCha()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Int_t AliHMPIDTracker::LoadClusters(TTree *pCluTree)
// Arguments: pEsd- pointer ESD; pClu- pointer to clusters for all chambers; pNmean - pointer to all function Nmean=f(time)
// Returns: error code, 0 if no errors
AliHMPIDRecon recon; //instance of reconstruction class, nothing important in ctor
- Float_t xPc,yPc;
+ Float_t xPc,yPc,xRa,yRa,theta,phi;
for(Int_t iTrk=0;iTrk<pEsd->GetNumberOfTracks();iTrk++){ //ESD tracks loop
AliESDtrack *pTrk = pEsd->GetTrack(iTrk); //get reconstructed track
- Int_t cham=IntTrkCha(pTrk,xPc,yPc); //get chamber intersected by this track
- if(cham<0) continue; //no intersection at all, go after next track
+ Int_t cham=IntTrkCha(pTrk,xPc,yPc,xRa,yRa,theta,phi); //get chamber intersected by this track
+ if(cham<0) { //no intersection at all, go after next track
+ pTrk->SetHMPIDtrk(0,0,0,0); //no intersection found
+ pTrk->SetHMPIDcluIdx (99,99999); //chamber not found, mip not yet considered
+ pTrk->SetHMPIDsignal(AliHMPIDRecon::kNotPerformed); //ring reconstruction not yet performed
+ continue;
+ }
+ pTrk->SetHMPIDtrk(xRa,yRa,theta,phi); //store initial infos
Double_t nmean=((TF1*)pNmean->At(3*cham))->Eval(pEsd->GetTimeStamp()); //C6F14 Nmean for this chamber
Int_t hvsec = AliHMPIDParam::InHVSector(xPc,yPc);
Double_t qthre=((TF1*)pQthre->At(6*cham+hvsec))->Eval(pEsd->GetTimeStamp());
Int_t RefitInward (AliESDEvent * ) {return 0;} //pure virtual from AliTracker
void UnloadClusters ( ) { } //pure virtual from AliTracker
//private part
- static Int_t IntTrkCha (AliESDtrack *pTrk,Float_t &xPc,Float_t &yPc ); //find track-PC intersection, retuns chamber ID
+ static Int_t IntTrkCha (AliESDtrack *pTrk,Float_t &xPc,Float_t &yPc,Float_t &xRa,Float_t &yRa,Float_t &theta,Float_t &phi);//find track-PC intersection, retuns chamber ID
static Int_t Recon (AliESDEvent *pEsd,TObjArray *pCluAll,TObjArray *pNmean=0,TObjArray *pQthre=0);//do actual job, returns status code
static Int_t ReconHiddenTrk(Int_t iCh,Int_t iHVsec,AliESDtrack *pTrk,TClonesArray *pClus,TObjArray *pNmean, TObjArray *pQthre);//do actual job with Hidden Track Algorithm
if(mtid>=0) continue; // only primaries
if(pTrack->GetPDG()->Charge()==0) continue;
AliESDtrack trk(pTrack);
- Float_t xPc,yPc;
- Int_t iCh=AliHMPIDTracker::IntTrkCha(&trk,xPc,yPc); //get chamber intersected by this track
- if(iCh<0) continue; //no intersection at all, go after next track
+ Float_t xPc,yPc,xRa,yRa,thRa,phRa;
+ Int_t iCh=AliHMPIDTracker::IntTrkCha(&trk,xPc,yPc,xRa,yRa,thRa,phRa); //get chamber intersected by this track
+ if(iCh<0) {
+ trk.SetHMPIDtrk(0,0,0,0); //no intersection found
+ trk.SetHMPIDcluIdx (99,99999); //chamber not found, mip not yet considered
+ trk.SetHMPIDsignal(AliHMPIDRecon::kNotPerformed); //ring reconstruction not yet performed
+ continue; //no intersection at all, go after next track
+ }
// Printf(" particle analyzed %d with mtid %d is %s hitting chamber %d",i,mtid,pTrack->GetName(),iCh);
+ trk.SetHMPIDcluIdx (iCh,99999); //chamber not found, mip not yet considered
+ trk.SetHMPIDtrk(xRa,yRa,thRa,phRa); //store initial infos
pEsd->AddTrack(&trk);
AliHMPIDTracker::Recon(pEsd,pH->CluLst(),pNmean,pQthre);
}// track loop
//find the chamber that intersects HMPID
AliESDtrack trk(pTrack);
- Float_t xPc,yPc;
- Int_t iCh=AliHMPIDTracker::IntTrkCha(&trk,xPc,yPc); //get chamber intersected by this track
+ Float_t xPc,yPc,xRa,yRa,thRa,phRa;
+ Int_t iCh=AliHMPIDTracker::IntTrkCha(&trk,xPc,yPc,xRa,yRa,thRa,phRa); //get chamber intersected by this track
if(iCh<0) continue; //no intersection at all, go after next track
+ trk.SetHMPIDtrk(xRa,yRa,thRa,phRa); //store initial infos
Float_t radX,radY,thetaTrk,phiTrk;
trk.GetHMPIDtrk(radX,radY,thetaTrk,phiTrk);
// Printf("simulated track theta %f phi %f",thetaTrk*rd,phiTrk*rd);
TObjArray *pClus = pH->CluLst();
- if(AliHMPIDTracker::ReconHiddenTrk(iCh,&trk,(TClonesArray *)pClus->At(iCh),pNmean,pQthre)!=0) continue;
+ if(AliHMPIDTracker::ReconHiddenTrk(iCh,0,&trk,(TClonesArray *)pClus->At(iCh),pNmean,pQthre)!=0) continue;
pEsd->AddTrack(&trk);
Double_t thetaCerSim = TMath::ACos(pTrack->Energy()/(1.292*pTrack->P()));
if(fDetBG->GetButton(kTRD )->GetState()) det+="TRD ";
if(fDetBG->GetButton(kTOF )->GetState()) det+="TOF ";
if(!fVerBG->GetButton(kNo)->GetState()) det+="HMPID ";
- char *sBatchName="Hbatch";
- FILE *fp=fopen(Form("%s.C",sBatchName),"w"); if(!fp){Info("CreateBatch","Cannot open output file: %s.C",sBatchName);return;}
+ char *sBatchName="sim";
+ FILE *fp=fopen(Form("%s.C",sBatchName),"w"); if(!fp){Info("CreateSim","Cannot open output file: %s.C",sBatchName);return;}
- fprintf(fp,"void %s(Int_t iNevt,Bool_t isDbg,char *sCfg)\n{\n",sBatchName);
+ fprintf(fp,"void %s(Int_t iNevt=1,Bool_t isDbg=kFALSE,char *sCfg=\"Config.C\")\n{\n",sBatchName);
fprintf(fp," gSystem->Exec(\"rm -rf hlt hough gphysi* fort* ZZZ* raw*\"); //remove garbage\n");
fprintf(fp," gBenchmark->Start(\"ALICE\"); TDatime time; //start benchmarking\n\n");
else if(fRawBG->GetButton(kRoo)->GetState()) fprintf(fp," pSim->SetWriteRawData(\"%s\",\"raw.root\"); //raw data as ROOT\n",det.Data());
fprintf(fp," pSim->SetRunHLT(\"\"); //no HLT stuff\n");
- fprintf(fp," pSim->Run(iNevt); //run iNevt events\n delete pSim;\n\n");
+ fprintf(fp," pSim->SetQA(kFALSE); //no QA\n");
+ fprintf(fp," pSim->Run(iNevt); //run iNevt events\n delete pSim;\n\n");
}//sim section
+ fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <sim.C>: Start time: \";time.Print();\n");
+ fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <sim.C>: Stop time: \";time.Set(); time.Print();\n");
+ fprintf(fp," gBenchmark->Show(\"ALICE\");\n");
+
+ fprintf(fp," gSystem->Exec(\"aliroot rec.C\");\n");
+ fprintf(fp," gSystem->Exec(\"touch ZZZ______finished_______SSS\");\n}\n");
+ fclose(fp);
+ char *sBatchName="rec";
+ FILE *fp=fopen(Form("%s.C",sBatchName),"w"); if(!fp){Info("CreateRec","Cannot open output file: %s.C",sBatchName);return;}
+
+ fprintf(fp,"void %s()\n{\n",sBatchName);
+ fprintf(fp," gSystem->Exec(\"rm -rf RRR* \"); //remove garbage\n");
if(fRecB->GetState()){
fprintf(fp," AliReconstruction *pRec=new AliReconstruction;\n");
+ fprintf(fp," gBenchmark->Start(\"ALICE\"); TDatime time; //start benchmarking\n\n");
//---------------------------------------------
if (fTrkBG->GetButton(kRecoPar)->GetState())
fprintf(fp," pRec->Run();delete pRec;\n\n");
}//rec part
//benchmarks
- fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <my/Batch.C>: Start time: \";time.Print();\n");
- fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <my/Batch.C>: Stop time: \";time.Set(); time.Print();\n");
+ fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <rec.C>: Start time: \";time.Print();\n");
+ fprintf(fp," cout<<\"!!!!!!!!!!!!Info in <rec.C>: Stop time: \";time.Set(); time.Print();\n");
fprintf(fp," gBenchmark->Show(\"ALICE\");\n");
- fprintf(fp," gSystem->Exec(\"touch ZZZ______finished_______ZZZ\");\n}\n");
+ fprintf(fp," gSystem->Exec(\"touch ZZZ______finished_______RRR\");\n}\n");
fclose(fp);
}//WriteBatch()
//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
AliHMPIDRecon rec;
for(Int_t iTrk=0;iTrk<pEsd->GetNumberOfTracks();iTrk++){//tracks loop to collect cerenkov rings and intersection points
AliESDtrack *pTrk=pEsd->GetTrack(iTrk); Int_t ch=pTrk->GetHMPIDcluIdx(); //get track and chamber intersected by it
- if(ch<0) continue; //this track does not intersect any chamber
- Float_t thRa,phRa,xRa,yRa; pTrk->GetHMPIDtrk(xRa,yRa,thRa,phRa); //get info on current track
- ch/=1000000;
- Float_t xPc=0,yPc=0; AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc); //find again intersection of track with PC--> it is not stored in ESD!
+ ch/=1000000;
+ if(ch<AliHMPIDParam::AliHMPIDParam::kMinCh||ch>AliHMPIDParam::kMaxCh) continue;//this track does not intersect any chamber
+ Float_t xPc,yPc,xRa,yRa,thRa,phRa;
+ Int_t chamb = AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc,xRa,yRa,thRa,phRa); //find again intersection of track with PC--> it is not stored in ESD!
+ if(ch!=chamb){(" CHAMBER MISMATCH: in ESDTrack chamber %i - in IntTrkCha chamber %i",ch,chamb);End();}
Int_t npTrk = fRenTxC[ch]->SetNextPoint(xPc,yPc); //add this intersection point
Float_t ckov=pTrk->GetHMPIDsignal(); //get ckov angle stored for this track
if(ckov>0){
else if (symbol==kTrack || symbol==kRing) {
if(symbol==kRing) index = b->GetUniqueID();
AliESDtrack *pTrk=fEsd->GetTrack(index);
- Float_t thRa,phRa,xRa,yRa; pTrk->GetHMPIDtrk(xRa,yRa,thRa,phRa);
- Float_t xPc,yPc; AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc);
- text0="";text0.Append(Form("TRACK n.%d: x %6.2f y %6.2f at PC plane",index,xPc,yPc));
+ Float_t xPc,yPc,thRa,phRa,xRa,yRa;
+ Int_t ch = AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc,xRa,yRa,thRa,phRa);
+ text0="";text0.Append(Form("TRACK n.%d: x %6.2f y %6.2f at PC plane in chamber %i",index,xPc,yPc,ch));
text2="";text2.Append(Form("p = %7.2f GeV/c",pTrk->GetP()));
Float_t ckov=pTrk->GetHMPIDsignal();
Double_t prob[5];
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff