#include <AlgCalDetPID.h>
Inheritance diagram for AlgCalDetPID:
Public Member Functions | |
| AlgCalDetPID () | |
| virtual | ~AlgCalDetPID () |
| virtual void | RunAlg (AlgConfig &ac, CandHandle &ch, CandContext &cx) |
| virtual void | Trace (const char *c) const |
Private Member Functions | |
| bool | HitTimeAnal (const CandDigitListHandle *cdlh, Double_t toftime, Float_t lowwin, Float_t highwin, Int_t &nbefore, Int_t &nafter, Int_t &nout, Int_t &nin, Float_t &chi2) const |
| bool | IsAnOverlap (AlgConfig &ac, Int_t nbefore, Int_t nafter, Int_t nout, Int_t nin) |
| bool | CERInTimeWin (AlgConfig &ac, const CandCalDetSIHandle *csih, Double_t abstime, Double_t toftime, Float_t win0low, Float_t win0high, Float_t win1low, Float_t win1high, Float_t win2low, Float_t win2high, bool &cer0, bool &cer1, bool &cer2) |
| Float_t | OverlapLikelihood (TH1F *ehfar, TH1F *ehnear) const |
Private Attributes | |
| Bool_t | fFirstEvent |
| Int_t | fNAllowedBefore |
| Int_t | fNAllowedAfter |
| Int_t | fNAllowedOut |
| Int_t | fNRequiredIn |
| Int_t | fTimeFileNum |
| TFile * | fTimeFile |
| TH1F * | fThist [2] |
|
|
Definition at line 46 of file AlgCalDetPID.cxx. 00046 : 00047 fFirstEvent(kTRUE), 00048 fNAllowedBefore(1), 00049 fNAllowedAfter(1), 00050 fNAllowedOut(2), 00051 fNRequiredIn(1), 00052 fTimeFileNum(0) 00053 { 00054 fTimeFile=0; 00055 fThist[0]=0; 00056 fThist[1]=0; 00057 00058 }
|
|
|
Definition at line 61 of file AlgCalDetPID.cxx. References fTimeFile. 00062 {
00063 if(fTimeFile){
00064 if(fTimeFile->IsOpen()){
00065 fTimeFile->Close();
00066 delete fTimeFile;
00067 fTimeFile = 0;
00068 }
00069 }
00070
00071 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 601 of file AlgCalDetPID.cxx. References CandCalDetSIHandle::GetKovADC1(), CandCalDetSIHandle::GetKovADC2(), CandCalDetSIHandle::GetKovADC3(), CandCalDetSIHandle::GetKovTimeStamp1(), CandCalDetSIHandle::GetKovTimeStamp2(), CandCalDetSIHandle::GetKovTimeStamp3(), and MSG. Referenced by RunAlg(). 00607 {
00608 // returns in cerX true if counter X passes the selection
00609 // returns cer0&cer1&cer2
00610
00611 cer0=cer1=cer2=false;
00612 if(csih==0){
00613 MSG("CalDetPID",Msg::kError)
00614 <<"CERInTimeWin got a null CandCalDetSIHandle!"
00615 <<"Returning false!"<<endl;
00616 return false;
00617 }
00618 else{
00619 // in principle could use algorithm config to avoid considering
00620 // one of the counters for some particles
00621
00622 // candcaldetsihandle starts counting cer counters at 1
00623 const Double_t tdc_convert=1.5625;
00624
00625
00626 Int_t cer0adc=csih->GetKovADC2(); // yes, I really mean this
00627 Float_t cer0time= toftime-
00628 (csih->GetKovTimeStamp2()*tdc_convert -abstime/Munits::ns); // too
00629 Int_t cer1adc=csih->GetKovADC1(); // too
00630 Float_t cer1time= toftime -
00631 (csih->GetKovTimeStamp1()*tdc_convert -abstime/Munits::ns);//too
00632 Int_t cer2adc=csih->GetKovADC3();
00633 Float_t cer2time= toftime -
00634 (csih->GetKovTimeStamp3()*tdc_convert-abstime/Munits::ns);
00635
00636 MSG("CalDetPID",Msg::kDebug)<<"toftime - (kov2*tdc_convert -abstime/Munits::ns)= "
00637 <<toftime<<"- ("
00638 <<csih->GetKovTimeStamp2()*tdc_convert
00639 <<" - "<<abstime/Munits::ns<<")"
00640 <<" = "<<toftime -(csih->GetKovTimeStamp2()*tdc_convert -abstime/Munits::ns)
00641 <<endl;
00642
00643 MSG("CalDetPID",Msg::kDebug)<<"toftime - (kov1*tdc_convert -abstime/Munits::ns)= "
00644 <<toftime<<"- ("
00645 <<csih->GetKovTimeStamp1()*tdc_convert
00646 <<" - "<<abstime/Munits::ns<<")"
00647 <<" = "<<toftime -(csih->GetKovTimeStamp1()*tdc_convert -abstime/Munits::ns)
00648 <<endl;
00649
00650 MSG("CalDetPID",Msg::kDebug)<<"toftime - (kov3*tdc_convert -abstime/Munits::ns)= "
00651 <<toftime<<"- ("
00652 <<csih->GetKovTimeStamp3()*tdc_convert
00653 <<" - "<<abstime/Munits::ns<<")"
00654 <<" = "<<toftime -(csih->GetKovTimeStamp3()*tdc_convert -abstime/Munits::ns)
00655 <<endl;
00656 const float nullwin=100000;
00657 // set cerX true if counter X passes the selection
00658 if(cer0adc>0){
00659 if(win0low==-nullwin && win0high==nullwin) cer0=true;
00660 else if((cer0time<win0low) || (cer0time>win0high)) cer0=false;
00661 else cer0=true;
00662 }
00663 else cer0=true;
00664
00665 if(cer1adc>0){
00666 if(win1low==-nullwin && win1high==nullwin) cer1=true;
00667 else if((cer1time<win1low) || (cer1time>win1high)) cer1=false;
00668 else cer1=true;
00669 }
00670 else cer1=true;
00671
00672 if(cer2adc>0){
00673 if(win2low==-nullwin && win2high==nullwin) cer2=true;
00674 else if((cer2time<win2low) || (cer2time>win2high)) cer2=false;
00675 else cer2=true;
00676 }
00677 else cer2=true;
00678
00679 return (cer0&cer1&cer2);
00680 }
00681
00682 }
|
|
||||||||||||||||||||||||||||||||||||||||
|
Definition at line 487 of file AlgCalDetPID.cxx. References fThist, CandHandle::GetDaughterIterator(), MSG, and OverlapLikelihood(). Referenced by RunAlg(). 00492 {
00493 // returns false if there is an error
00494
00495 // zero counters
00496 nbefore=0; nafter=0; nout=0; nin=0;
00497
00498 //make histograms to hold timing of current event
00499 int NBFAR = fThist[0]->GetNbinsX();
00500 Axis_t LOWFAR = fThist[0]->GetBinLowEdge(1);
00501 Axis_t HIFAR = fThist[0]->GetBinLowEdge(NBFAR)+
00502 fThist[0]->GetBinWidth(NBFAR);
00503
00504 int NBNEAR = fThist[1]->GetNbinsX();
00505 Axis_t LOWNEAR = fThist[1]->GetBinLowEdge(1);
00506 Axis_t HINEAR = fThist[1]->GetBinLowEdge(NBNEAR)+
00507 fThist[1]->GetBinWidth(NBNEAR);
00508
00509 MSG("CalDetPID",Msg::kDebug)<<"making event thists "
00510 <<" NBFAR "<<NBFAR
00511 <<" LOWFAR "<<LOWFAR
00512 <<" HIFAR "<<HIFAR
00513 <<" NBNEAR "<<NBNEAR
00514 <<" LOWNEAR "<<LOWNEAR
00515 <<" HINEAR "<<HINEAR<<endl;
00516
00517 TH1F ehfar("ehfar","ehfar",NBFAR,LOWFAR,HIFAR);
00518 TH1F ehnear("ehnear","ehnear",NBNEAR,LOWNEAR,HINEAR);
00519
00520 // check digit list
00521 CandDigitHandleItr hi(cdlh->GetDaughterIterator());
00522 if(!hi.IsValid()){
00523 MSG("CalDetPID",Msg::kError)<<"DaugherIter not valid"<<endl;
00524 return false;
00525 }
00526 //begin looping over hits
00527 for(hi.Reset();hi.IsValid();hi.Next()){
00528 // to check that it isn't a cosmic counter or unknown
00529 const PlaneView::PlaneView_t pv=
00530 (*hi)->GetPlexSEIdAltL().GetPlaneView();
00531 // to check that it IS connected to scint
00532 const PlexHandle ph(*((*hi)->GetVldContext()));
00533 const ReadoutType::Readout_t rt=
00534 ph.GetReadoutType((*hi)->GetChannelId());
00535 const ElecType::Elec_t et = (*hi)->GetChannelId().GetElecType();
00536 /*
00537 // the logic used to be
00538 if((pv==PlaneView::kA)||(pv==PlaneView::kB)
00539 ||(pv==PlaneView::kUnknown)) {
00540 */
00541 /* changing the logic so that it will consider near detector hits
00542 if( !(pv==PlaneView::kU || pv==PlaneView::kV)
00543 || rt!=ReadoutType::kScintStrip || (et!=ElecType::kVA)){
00544 // don't do anything with this hit
00545 continue;
00546 }
00547 */
00548 if( !(pv==PlaneView::kU || pv==PlaneView::kV)
00549 || rt!=ReadoutType::kScintStrip){
00550 // don't do anything with this hit
00551 continue;
00552 }
00553
00554 Float_t time = (*hi)->GetSubtractedTime(CalTimeType::kT0);
00555 time/=Munits::nanosecond;
00556 MSG("CalDetPID",Msg::kDebug)<<"hittime - toftime= "
00557 <<time<<"-"<<toftime
00558 <<" = "<<time-toftime<<endl;
00559 time-=toftime;
00560 if(et==ElecType::kVA){
00561 ehfar.Fill(time);
00562 if(time<lowwin){nbefore++; nout++;}
00563 else if(time>highwin){nafter++; nout++;}
00564 else nin++;
00565 }
00566 else if(et==ElecType::kQIE){
00567 ehnear.Fill(time);
00568 }
00569 }
00570 chi2 = OverlapLikelihood(&ehfar,&ehnear);
00571 return true;
00572 }
|
|
||||||||||||||||||||||||
|
Definition at line 574 of file AlgCalDetPID.cxx. References fNAllowedAfter, fNAllowedBefore, fNAllowedOut, fNRequiredIn, and Registry::Get(). Referenced by RunAlg(). 00577 {
00578 // returns false if there is no overlap
00579 // returns true if there is an overlap
00580 // get the number of hits allowed out of the time interval from
00581 // AlgConfig
00582
00583 Int_t tmpi;
00584 if(ac.Get("NAllowedBefore", tmpi)) fNAllowedBefore=tmpi;
00585 if(ac.Get("NAllowedAfter", tmpi)) fNAllowedAfter=tmpi;
00586 if(ac.Get("NAllowedOut", tmpi)) fNAllowedOut=tmpi;
00587 if(ac.Get("NRequiredIn", tmpi)) fNRequiredIn=tmpi;
00588
00589 if(nout<=fNAllowedOut&&
00590 nbefore<=fNAllowedBefore&&
00591 nafter<=fNAllowedAfter&&
00592 nin>=fNRequiredIn){
00593 return false;
00594 }
00595 else{
00596 return true;
00597 }
00598
00599 }
|
|
||||||||||||
|
Definition at line 684 of file AlgCalDetPID.cxx. Referenced by HitTimeAnal(). 00685 {
00686 //now loop over the histogram bins to calculate the chi2
00687 Int_t nhits=0;
00688 Float_t flsum=0.;
00689
00690 TH1F *eventhist[2];
00691 eventhist[0] = ehfar;
00692 eventhist[1] = ehnear;
00693
00694 for(int i=0;i<2;i++){
00695 const int nevents = (const int)eventhist[i]->Integral();
00696 nhits+=nevents;
00697 for(int j=1;j<=eventhist[i]->GetNbinsX();j++){
00698 if(eventhist[i]->GetBinContent(j)>0){
00699 int x = (int)eventhist[i]->GetBinContent(j);
00700 float mean = fThist[i]->GetBinContent(j)*nevents;
00701 MSG("CalDetPID",Msg::kDebug)<<" crate "<<i
00702 <<" at time "<<eventhist[i]->GetBinCenter(j)
00703 <<" x "<<x<<" mean "<<mean<<endl;
00704
00705 if(mean==0){
00706 continue;
00707 }
00708
00709 double lnprob = 1.*x*log(mean)-mean;
00710 double lnnorm = 1.*x*log(1.*x)-1.*x;
00711
00712 flsum+=(lnprob-lnnorm);
00713
00714 MSG("CalDetPID",Msg::kDebug)<<"i "<<i
00715 <<" ln prob "<<lnprob
00716 <<" ln norm "<<lnnorm
00717 <<" fl sum "<<flsum<<endl;
00718 }
00719 }
00720 MSG("CalDetPID",Msg::kDebug)<<"i "<<i
00721 <<" flsum "<<flsum<<endl;
00722 }
00723 if(nhits>0){
00724 flsum = -2*(flsum/(1.*nhits));
00725 }
00726 else{
00727 flsum = 0;
00728 }
00729 return flsum;
00730 }
|
|
||||||||||||||||
|
Implements AlgBase. Definition at line 74 of file AlgCalDetPID.cxx. References abs(), CalDetParticleType::AsString(), CERInTimeWin(), fFirstEvent, CandRecord::FindCandHandle(), RecDataRecord< T >::FindComponent(), fThist, fTimeFile, fTimeFileNum, Registry::Get(), CandDigitListHandle::GetAbsTime(), CandContext::GetCandRecord(), CalDetCERRange::GetCER0High(), CalDetCERRange::GetCER0Low(), CalDetCERRange::GetCER1High(), CalDetCERRange::GetCER1Low(), CalDetCERRange::GetCER2High(), CalDetCERRange::GetCER2Low(), CandContext::GetDataIn(), MomNavigator::GetFragment(), CandCalDetSIHandle::GetKovADC1(), CandCalDetSIHandle::GetKovADC2(), CandCalDetSIHandle::GetKovADC3(), CandContext::GetMom(), DbiResultPtr< T >::GetNumRows(), CalDetCERTimeWin::GetParticleType(), CalDetOverlapWin::GetParticleType(), CalDetCERRange::GetParticleType(), CalDetTOFRange::GetParticleType(), DbiResultPtr< T >::GetRow(), VldContext::GetSimFlag(), CalDetTOFRange::GetTDC0High(), CalDetTOFRange::GetTDC0Low(), CalDetTOFRange::GetTDC2High(), CalDetTOFRange::GetTDC2Low(), CalDetTOFRange::GetTDC2MinusTDC0High(), CalDetTOFRange::GetTDC2MinusTDC0Low(), CandCalDetSIHandle::GetTofTDC0(), CandCalDetSIHandle::GetTofTDC1(), CandCalDetSIHandle::GetTofTDC2(), CandCalDetSIHandle::GetTofTimeStamp(), DbiResultPtr< T >::GetValidityRec(), RecMinos::GetVldContext(), CalDetCERTimeWin::GetWin0High(), CalDetCERTimeWin::GetWin0Low(), CalDetCERTimeWin::GetWin1High(), CalDetCERTimeWin::GetWin1Low(), CalDetCERTimeWin::GetWin2High(), CalDetCERTimeWin::GetWin2Low(), CalDetOverlapWin::GetWinHigh(), CalDetOverlapWin::GetWinLow(), HitTimeAnal(), IsAnOverlap(), MSG, CandCalDetPIDHandle::SetInCERTime(), CandCalDetPIDHandle::SetInCERTimeBits(), CandCalDetPIDHandle::SetNoOverlap(), CandCalDetPIDHandle::SetNoOverlapBits(), CandCalDetPIDHandle::SetOLChi2(), and CandCalDetPIDHandle::SetPIDType(). 00076 {
00077 bool ok=true;
00078
00079 MSG("CalDetPID", Msg::kDebug)
00080 << "Starting AlgCalDetPID::RunAlg()" << endl;
00081
00082 if(fFirstEvent){
00083 //get filename from alg config
00084 Int_t tfn;
00085 if(ac.Get("TimeFileNum", tfn)){
00086 fTimeFileNum = tfn;
00087 MSG("CalDetPID", Msg::kDebug)<<"got TimeFileNum from config "
00088 <<fTimeFileNum<<endl;
00089 }
00090 else{
00091 fTimeFileNum = 0;
00092 MSG("CalDetPID", Msg::kError)<<"Could not get TimeFileName from config."
00093 <<" Using default "<<fTimeFileNum
00094 <<" instead"<<endl;
00095 }
00096 string fTimeFileName;
00097 switch(fTimeFileNum){
00098 case 0: //far-far readout t11 2002
00099 fTimeFileName = "timinghists-2002-t11.root";
00100 break;
00101 case 1: //far-far readout t7 2002
00102 fTimeFileName = "timinghists-2002-t7.root";
00103 break;
00104 case 2: //far-near readout 2003
00105 fTimeFileName = "timinghists-2003-t7-nf.root";
00106 break;
00107 case 3: //near readout t7 2003
00108 fTimeFileName = "timinghists-2003-t7-n.root";
00109 break;
00110 case 4: //near readout t11 2003
00111 fTimeFileName = "timinghists-2003-t11.root";
00112 break;
00113 default:
00114 fTimeFileName = "timinghists.root";
00115 break;
00116 }
00117
00118 // open input root file
00119 string fname(getenv("SRT_PRIVATE_CONTEXT"));
00120 fname+=("/CalDetPID/data/"+fTimeFileName);
00121
00122 fTimeFile = new TFile(fname.c_str());
00123 if(!fTimeFile->IsOpen()){
00124 MSG("CalDetPID",Msg::kInfo)<<"Could not find "<<fname<<" in the base release. Will try the test release."<<endl;
00125 string fname_save=fname;
00126 fname=getenv("SRT_PUBLIC_CONTEXT");
00127 fname+="/CalDetPID/data/"+fTimeFileName;
00128 //crashes if it get's here because now fTimeFile goes out of scope
00129 //tv 6-16-2004
00130 // TFile *fTimeFile = new TFile(fname.c_str());
00131 fTimeFile = new TFile(fname.c_str());
00132 if(!fTimeFile->IsOpen()){
00133 MSG("CalDetPID",Msg::kError)<<"Could not open timing hist file"
00134 <<fname<<endl
00135 <<"Also tried "<<fname_save<<endl;
00136 return;
00137 }
00138 }
00139 MSG("CalDetPID",Msg::kInfo)<<"Opened timing hist file "<<fname<<endl;
00140 // fTimeFile->Print();
00141 for(int i=0;i<2;i++){
00142 char hn[50];
00143 sprintf(hn,"all_thist_c%d",i);
00144 fThist[i] = (TH1F *)fTimeFile->Get(hn);
00145 if(fThist[i]!=0){
00146 MSG("CalDetPID",Msg::kDebug)<<"Found histogram "<<hn<<endl;
00147 }
00148 else{
00149 MSG("CalDetPID",Msg::kError)<<"Could not find histogram "<<hn<<endl;
00150 }
00151 }
00152 fFirstEvent=kFALSE;
00153 }
00154
00155
00156 assert(ch.InheritsFrom("CandCalDetPIDHandle"));
00157 CandCalDetPIDHandle& cpidh = static_cast<CandCalDetPIDHandle&>(ch);
00158
00159 assert(cx.GetDataIn());
00160 const CandRecord* crec = cx.GetCandRecord();
00161 // Get VldContext
00162 const VldContext vc = *(crec->GetVldContext());
00163 // is this data?
00164 const bool isdata = (vc.GetSimFlag()&SimFlag::kData) ? true : false;
00165 // is this mc?
00166 const bool ismc = ((vc.GetSimFlag()&SimFlag::kMC) ||
00167 (vc.GetSimFlag()&SimFlag::kReroot)) ? true : false;
00168
00169 if(isdata){
00170 // for data we must have a caldetsi candidate
00171 const CandCalDetSIHandle *csih =
00172 dynamic_cast<const CandCalDetSIHandle *> (cx.GetDataIn());
00173 assert(csih);
00174 // we must also have some digits
00175 const CandDigitListHandle* cdlh=
00176 dynamic_cast<const CandDigitListHandle*>
00177 (crec->FindCandHandle("CandDigitListHandle","canddigitlist"));
00178 if(cdlh == 0) {//failure
00179 MSG("CalDetPID",Msg::kFatal)<<"cdlh==0"<<endl;
00180 assert(0);
00181 }
00182
00183 const Double_t toftime=
00184 (csih->GetTofTimeStamp()*1.5625 - cdlh->GetAbsTime()/Munits::ns);
00185
00186 // talk to the db
00187
00188 DbiResultPtr<CalDetTOFRange> rp_tof(vc,0);
00189 if(!rp_tof.GetValidityRec()){
00190 MSG("CalDetPID", Msg::kError)
00191 <<"Could not look up selection ranges "
00192 <<"from CalDetTOFRange table!\n"
00193 <<"Bad DbiResultPtr!?"
00194 <<endl;
00195 ok=false;
00196 }
00197
00198 DbiResultPtr<CalDetCERRange> rp_cer(vc,0);
00199 if(!rp_cer.GetValidityRec()){
00200 MSG("CalDetPID", Msg::kError)
00201 <<"Could not look up selection ranges "
00202 <<"from CalDetCERRange table!\n"
00203 <<"Bad DbiResultPtr!?"
00204 <<endl;
00205 ok=false;
00206 }
00207
00208
00209
00210 // The algorithm:
00211 // Look through each entry in the db.
00212 // There should not be many of them.
00213 // One entry may correspond to one or more types of particles.
00214 // The PARTICLETYPE entry in the table rows will be used to set high
00215 // the corresponding bits in case the data is in the appropriate range.
00216 // This can be done seperately for tof and cer with the results anded
00217 // at the end.
00218
00219 // do id by tof
00220 UInt_t tof_pid=0;
00221 if(ok){
00222 // here is the tof data
00223 const Int_t tdc0=csih->GetTofTDC0();
00224 const Int_t tdc1=csih->GetTofTDC1();
00225 const Int_t tdc2=csih->GetTofTDC2();
00226 MSG("CalDetPID", Msg::kDebug)
00227 <<"TOF 0 1 2 = "<<tdc0<<" "<<tdc1<<" "<<tdc2<<endl;
00228 // loop over rows
00229 const int nrows_tof = rp_tof.GetNumRows();
00230 for(int i=0; i<nrows_tof; i++){
00231 // get a row
00232 const CalDetTOFRange* tr = rp_tof.GetRow(i);
00233 if(!tr){
00234 // something is really hosed
00235 MSG("CalDetPID", Msg::kFatal)
00236 <<"Could not look up row "<<i
00237 <<" in CalDetTOFRange table!"
00238 <<endl;
00239 ok=false;
00240 assert(0);
00241 }
00242 else{
00243 // for now, only test on tdc2 and tdc0
00244 // tdc1 was never used...
00245 // when we ran middle and ds counters the middle
00246 // was plugged into tdc2
00247 if((tdc0 < tr->GetTDC0High())
00248 && (tdc0 > tr->GetTDC0Low())
00249 && (tdc2 < tr->GetTDC2High())
00250 && (tdc2 > tr->GetTDC2Low())
00251 && (tdc2-tdc0 < tr->GetTDC2MinusTDC0High())
00252 && (tdc2-tdc0 > tr->GetTDC2MinusTDC0Low())
00253 )
00254 tof_pid|=(tr->GetParticleType());
00255 }
00256 }
00257 }
00258
00259 // do id by cer
00260 UInt_t cer_pid=0;
00261 if(ok){
00262 // here is the cer data
00263 const Int_t cer0=csih->GetKovADC2();
00264 const Int_t cer1=csih->GetKovADC1();
00265 const Int_t cer2=csih->GetKovADC3();
00266 MSG("CalDetPID", Msg::kDebug)
00267 <<"CER 0 1 2 = "<<cer0<<" "<<cer1<<" "<<cer2<<endl;
00268 // loop over rows
00269 const int nrows_cer = rp_cer.GetNumRows();
00270 for(int i=0; i<nrows_cer; i++){
00271 const CalDetCERRange* cr = rp_cer.GetRow(i);
00272 if(!cr){
00273 // something is really hosed
00274 MSG("CalDetPID", Msg::kFatal)
00275 <<"Could not look up row "<<i
00276 <<" in CalDetCERRange table!"
00277 <<endl;
00278 ok=false;
00279 assert(0);
00280 }
00281 else{
00282 if((cer0<cr->GetCER0High())
00283 &&(cer0>cr->GetCER0Low())
00284 &&(cer1<cr->GetCER1High())
00285 &&(cer1>cr->GetCER1Low())
00286 &&(cer2<cr->GetCER2High())
00287 &&(cer2>cr->GetCER2Low())
00288 )
00289 cer_pid|=(cr->GetParticleType());
00290 }
00291 }
00292 }
00293
00294 MSG("CalDetPID", Msg::kDebug)
00295 <<"tof_pid = "<<hex<<tof_pid<<dec
00296 <<", Identified via TOF as "
00297 <<CalDetParticleType::AsString(tof_pid)<<endl;
00298
00299 MSG("CalDetPID", Msg::kDebug)
00300 <<"cer_pid = "<<hex<<cer_pid<<dec
00301 <<", Identified via CER as "
00302 <<CalDetParticleType::AsString(cer_pid)<<endl;
00303
00304 // Now we have the basic pid information.
00305 // Look to see if this is an overlapper.
00306
00307 // For things that are not overlappers, set the appropriate bits high
00308 // otherwise do not set them
00309 // If we cannot talk to the db then cannot_run_overlap=true;
00310 UInt_t overlap_pid=0;
00311 bool cannot_run_overlap=false;
00312 Float_t chi2 = -1.;
00313 DbiResultPtr<CalDetOverlapWin> rp_cow(vc,0);
00314 if(!rp_cow.GetValidityRec()){
00315 MSG("CalDetPID", Msg::kWarning)
00316 <<"Could not look up selection ranges "
00317 <<"from CalDetOverlapWin table!\n"
00318 <<"Bad DbiResultPtr!?"
00319 <<endl;
00320 // but maybe it could live with this?
00321 cannot_run_overlap=true;
00322 ok=false;
00323 }
00324 else{
00325 const int nrows_cow=rp_cow.GetNumRows();
00326 for(int i=0; i<nrows_cow; i++){
00327 const CalDetOverlapWin* cow = rp_cow.GetRow(i);
00328 if(!cow){
00329 // something is really hosed
00330 MSG("CalDetPID", Msg::kFatal)
00331 <<"Could not look up row "<<i
00332 <<" in CalDetOverlapWin table!"
00333 <<endl;
00334 ok=false;
00335 assert(0);
00336 }
00337 else{
00338 const Float_t lowwin = cow->GetWinLow();
00339 const Float_t highwin = cow->GetWinHigh();
00340 Int_t nbefore,nafter,nout,nin;
00341 HitTimeAnal(cdlh,toftime,lowwin,highwin,
00342 nbefore,nafter,nout,nin,chi2);
00343 bool overlap = IsAnOverlap(ac,nbefore,nafter,
00344 nout,nin);
00345 if(!overlap){
00346 overlap_pid|=(cow->GetParticleType());
00347 }
00348 }
00349 }
00350 }
00351
00352 //
00353 // now work on certimewin
00354 UInt_t ctw_pid=0;
00355 bool cannot_run_ctw=false;
00356 DbiResultPtr<CalDetCERTimeWin> rp_ctw(vc,0);
00357 if(!rp_ctw.GetValidityRec()){
00358 MSG("CalDetPID", Msg::kError)
00359 <<"Could not look up selection ranges "
00360 <<"from CalDetCERTimeWin table!\n"
00361 <<"Bad DbiResultPtr!?"
00362 <<endl;
00363 // but maybe it could live with this?
00364 cannot_run_ctw=true;
00365 ok=false;
00366 }
00367 else{
00368 const int nrows_ctw=rp_ctw.GetNumRows();
00369 for(int i=0; i<nrows_ctw; i++){
00370 const CalDetCERTimeWin* ctw = rp_ctw.GetRow(i);
00371 if(!ctw){
00372 // something is really hosed
00373 MSG("CalDetPID", Msg::kFatal)
00374 <<"Could not look up row "<<i
00375 <<" in CalDetCERTimeWin table!"
00376 <<endl;
00377 ok=false;
00378 assert(0);
00379 }
00380 else{
00381 const Float_t win0low = ctw->GetWin0Low();
00382 const Float_t win0high = ctw->GetWin0High();
00383 const Float_t win1low = ctw->GetWin1Low();
00384 const Float_t win1high = ctw->GetWin1High();
00385 const Float_t win2low = ctw->GetWin2Low();
00386 const Float_t win2high = ctw->GetWin2High();
00387 bool cer0=false;
00388 bool cer1=false;
00389 bool cer2=false;
00390 MSG("CalDetPID", Msg::kDebug)
00391 <<(*ctw)<<endl;
00392 const bool in_win = CERInTimeWin(ac, csih,
00393 cdlh->GetAbsTime(),
00394 toftime,
00395 win0low,win0high,
00396 win1low,win1high,
00397 win2low,win2high,
00398 cer0, cer1, cer2);
00399 if(in_win){
00400 ctw_pid|=(ctw->GetParticleType());
00401 }
00402 }
00403 }
00404 }
00405
00406 // final step:
00407 // and the tof and cer results together
00408 const UInt_t pid_result = (tof_pid&cer_pid);
00409 MSG("CalDetPID", Msg::kDebug)
00410 <<"result = "<<hex<<pid_result<<dec
00411 <<", Identified as "
00412 <<CalDetParticleType::AsString(pid_result)<<endl;
00413 cpidh.SetPIDType(pid_result);
00414
00415 MSG("CalDetPID", Msg::kDebug)
00416 <<"overlap_pid = "<<hex<<overlap_pid<<dec
00417 <<", Non overlapping for "
00418 <<CalDetParticleType::AsString(overlap_pid)<<endl;
00419
00420 cpidh.SetNoOverlapBits(overlap_pid);
00421
00422 MSG("CalDetPID", Msg::kDebug)
00423 <<"ctw_pid = "<<hex<<ctw_pid<<dec
00424 <<", Cer in time window for "
00425 <<CalDetParticleType::AsString(ctw_pid)<<endl;
00426
00427 cpidh.SetInCERTimeBits(ctw_pid);
00428
00429 Bool_t nooverlap;
00430 if((pid_result&overlap_pid)==pid_result) nooverlap=kTRUE;
00431 else nooverlap=kFALSE;
00432 MSG("CalDetPID", Msg::kDebug)
00433 <<"NoOverlap = "<<( (nooverlap) ? "true": "false")<<endl;
00434 cpidh.SetNoOverlap(nooverlap);
00435
00436 Bool_t incertime;
00437 if((pid_result&ctw_pid)==pid_result) incertime=kTRUE;
00438 else incertime=kFALSE;
00439 MSG("CalDetPID", Msg::kDebug)
00440 <<"InCERTime = "<<( (incertime) ? "true": "false")<<endl;
00441 cpidh.SetInCERTime(incertime);
00442
00443 cpidh.SetOLChi2(chi2);
00444 }
00445 else if(ismc){
00446 cpidh.SetPIDType(CalDetParticleType::kUnknown);
00447 cpidh.SetNoOverlapBits(0);
00448 cpidh.SetInCERTimeBits(0);
00449 cpidh.SetNoOverlap(kTRUE);
00450 cpidh.SetInCERTime(kTRUE);
00451 cpidh.SetOLChi2(0.);
00452 // look for the primary in stdhep
00453 const SimSnarlRecord* simrec = dynamic_cast<const SimSnarlRecord*>
00454 (cx.GetMom()->GetFragment("SimSnarlRecord"));
00455 if(simrec){
00456 const TClonesArray* stdhep = dynamic_cast<const TClonesArray*>
00457 (simrec->FindComponent("TClonesArray", "StdHep"));
00458 if(stdhep){
00459 // take the first in the list
00460 const TParticle* part = dynamic_cast<const TParticle*>
00461 (stdhep->At(0));
00462 if(part){
00463 Int_t ipdg = part->GetPdgCode();
00464 if(abs(ipdg)==11)
00465 cpidh.SetPIDType(CalDetParticleType::kElectron);
00466 else if(abs(ipdg)==13)
00467 cpidh.SetPIDType(CalDetParticleType::kMuon);
00468 else if(abs(ipdg)==211)
00469 cpidh.SetPIDType(CalDetParticleType::kPion);
00470 else if(abs(ipdg)==321)
00471 cpidh.SetPIDType(CalDetParticleType::kKaon);
00472 else if(abs(ipdg)==2212)
00473 cpidh.SetPIDType(CalDetParticleType::kProton);
00474 }
00475 }
00476 }
00477 }
00478
00479
00480 }
|
|
|
Reimplemented from AlgBase. Definition at line 483 of file AlgCalDetPID.cxx. 00484 {
00485 }
|
|
|
Definition at line 45 of file AlgCalDetPID.h. Referenced by RunAlg(). |
|
|
Definition at line 47 of file AlgCalDetPID.h. Referenced by IsAnOverlap(). |
|
|
Definition at line 46 of file AlgCalDetPID.h. Referenced by IsAnOverlap(). |
|
|
Definition at line 48 of file AlgCalDetPID.h. Referenced by IsAnOverlap(). |
|
|
Definition at line 49 of file AlgCalDetPID.h. Referenced by IsAnOverlap(). |
|
|
Definition at line 52 of file AlgCalDetPID.h. Referenced by HitTimeAnal(), OverlapLikelihood(), and RunAlg(). |
|
|
Definition at line 51 of file AlgCalDetPID.h. Referenced by RunAlg(), and ~AlgCalDetPID(). |
|
|
Definition at line 50 of file AlgCalDetPID.h. Referenced by RunAlg(). |
1.3.9.1