#include <AlgEventSSList.h>
Inheritance diagram for AlgEventSSList:

Public Member Functions | |
| AlgEventSSList () | |
| virtual | ~AlgEventSSList () |
| virtual void | RunAlg (AlgConfig &ac, CandHandle &ch, CandContext &cx) |
| void | BuildEventFromUnassoc (AlgConfig &ac, CandHandle &ch, CandContext &cx, CandSliceHandle *slice) |
| virtual void | Trace (const char *c) const |
| void | SetPrimaryShowers (CandHandle &ch, AlgConfig &ac) |
| void | ReConstructShowers (AlgConfig &ac, CandHandle &ch, CandContext &cx) |
| void | FillRecoList (CandSliceHandle *, CandShowerHandleItr *, CandTrackHandleItr *, TObjArray &) |
| void | AddObjectToEvent (CandRecoHandle *reco, TObjArray *, TObjArray *, Bool_t) |
| void | AddStripToEvent (CandEventHandle *ev, CandShowerHandle *shw, CandSubShowerSRListHandle *subshwlist, CandStripHandle *closeststrip) |
| void | CreatePrimaryShower (AlgConfig &ac, CandContext &cxx, CandEventHandle *closestevent, CandShowerListHandle *showerlist, CandSubShowerSRListHandle *subshowerlist, CandStripHandle *closeststrip) |
| void | FindUnassociated (CandSliceHandleItr &sliceItr, CandEventHandleItr *, TObjArray &unassociated) |
| void | ReFillDist2Map (AlgConfig &ac, CandStripHandle *closeststrip, CandEventHandle *closestevent, TObjArray &unassociated, std::vector< std::map< const CandEventHandle *, Double_t > > &dist2map) |
| void | FillDist2Map (AlgConfig &ac, TObjArray &unassociated, CandHandle &ch, std::vector< std::map< const CandEventHandle *, Double_t > > &dist2map) |
|
|
Definition at line 62 of file AlgEventSSList.cxx. 00063 {
00064 }
|
|
|
Definition at line 67 of file AlgEventSSList.cxx. 00068 {
00069 }
|
|
||||||||||||||||||||
|
Definition at line 1821 of file AlgEventSSList.cxx. References MSG. 01822 {
01823 if (!merge) {
01824 // if this is the first match we have found for this object,
01825 //add object to list
01826 // and point prevlist to this object list
01827 objectlist->Add(reco);
01828 MSG("EventSS",Msg::kVerbose) << " adding to list\n";
01829 }
01830 else {
01831 MSG("EventSS",Msg::kVerbose) << " combining lists\n";
01832 // this object matches more than one event - merge their member
01833 // objects to list prevlist, and add the new object
01834 for (Int_t ireco2=0; ireco2<=objectlist->GetLast();ireco2++)
01835 prevlist->Add(objectlist->At(ireco2));
01836 prevlist->Add(reco);
01837 objectlist->Clear();
01838 }
01839 }
|
|
||||||||||||||||||||
|
Add a strip to the event daughter list, the primary shower, and place the modified cluster on the cluster list Definition at line 966 of file AlgEventSSList.cxx. References CandHandle::AddDaughterLink(), CandShowerSRHandle::AddSubShower(), CandSubShowerSRHandle::DupHandle(), CandSubShowerSRHandle::GetAvgDev(), CandRecoHandle::GetCandSlice(), CandRecoHandle::GetCharge(), CandHandle::GetDaughterIterator(), CandSubShowerSRHandle::GetEnergy(), CandShowerSRHandle::GetLastSubShower(), Calibrator::GetMIP(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandSubShowerSRHandle::GetPlaneView(), CandSubShowerSRHandle::GetSlope(), CandShowerSRHandle::GetSubShower(), CandStripHandle::GetTPos(), CandHandle::GetUidInt(), CandRecoHandle::GetVtxU(), CandRecoHandle::GetVtxV(), CandRecoHandle::GetVtxZ(), CandStripHandle::GetZPos(), Calibrator::Instance(), CandHandle::IsEquivalent(), MSG, CandHandle::RemoveDaughter(), CandShowerSRHandle::RemoveSubShower(), CandRecoHandle::SetCandSlice(), and CandSubShowerSRHandle::SetEnergy(). Referenced by RunAlg(). 00970 {
00971 MSG("EventSS",Msg::kVerbose) << "In AddStripToEvent:" << endl;
00972 closestevent->AddDaughterLink(*closeststrip);
00973 shower->AddDaughterLink(*closeststrip);
00974 if(subshowerlist && shower->InheritsFrom("CandShowerSRHandle")) {
00975 Calibrator& calibrator=Calibrator::Instance();
00976 Int_t foundsubshower=0;
00977 Int_t closestSubShower=-1;
00978 Double_t BestStripDistance = 99999;
00979 Double_t StripDistance = 99999;
00980 CandShowerSRHandle *showerSR =
00981 dynamic_cast<CandShowerSRHandle*> (shower);
00982 MSG("EventSS",Msg::kVerbose) << "Shower UID = "
00983 << shower->GetUidInt() << endl;
00984
00985 for (Int_t isubshower=0;
00986 isubshower<showerSR->GetLastSubShower()+1;
00987 isubshower++) {
00988 CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
00989 (showerSR->GetSubShower(isubshower));
00990 //add strip to cluster of the same view
00991 if (ssl->GetPlaneView()==closeststrip->GetPlaneView()) {
00992 //decide which is most appropriate subshower to add strip to
00993 //based on distance to/spread of subshower:
00994 Double_t tposVtx = 0;
00995 if(ssl->GetPlaneView()==PlaneView::kU ||
00996 ssl->GetPlaneView()==PlaneView::kX)
00997 tposVtx = ssl->GetVtxU();
00998 else if(ssl->GetPlaneView()==PlaneView::kV ||
00999 ssl->GetPlaneView()==PlaneView::kY)
01000 tposVtx = ssl->GetVtxV();
01001 StripDistance = 99999;
01002 if(ssl->GetAvgDev()>0) {
01003 StripDistance =
01004 (TMath::Abs(closeststrip->GetTPos() -
01005 (tposVtx + ssl->GetSlope() *
01006 (closeststrip->GetZPos() -
01007 ssl->GetVtxZ()))) *
01008 TMath::Cos(TMath::ATan(ssl->GetSlope()))) *
01009 calibrator.GetMIP(closeststrip->
01010 GetCharge(CalDigitType::kSigCorr)) /
01011 ssl->GetAvgDev();
01012 }
01013 if(StripDistance>9999){ //bad slope
01014 StripDistance = 9999 +
01015 TMath::Sqrt(TMath::Power(TMath::Abs(closeststrip->GetZPos() -
01016 ssl->GetVtxZ()),2) +
01017 TMath::Power(TMath::Abs(closeststrip->GetTPos() -
01018 tposVtx),2));
01019 }
01020 if(StripDistance<BestStripDistance || foundsubshower==0){
01021 foundsubshower=1;
01022 BestStripDistance = StripDistance;
01023 closestSubShower = isubshower;
01024 }
01025 }
01026 }
01027 if(foundsubshower) {
01028 MSG("EventSS",Msg::kVerbose) << "adding strip to subshower "
01029 << closeststrip->GetPlane()
01030 << endl;
01031 CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
01032 (showerSR->GetSubShower(closestSubShower));
01033 //extract mip -> gev conversion:
01034 Double_t GeVperMIP = 0;
01035 if(ssl->GetEnergy()!=0) GeVperMIP = ssl->GetEnergy()/
01036 Double_t(calibrator.GetMIP(ssl->GetCharge(CalStripType::kSigCorr)));
01037 TIter csshItr(subshowerlist->GetDaughterIterator());
01038 while (CandSubShowerSRHandle *cssh =
01039 dynamic_cast<CandSubShowerSRHandle*>(csshItr())) {
01040 if(cssh->IsEquivalent(ssl) && ssl->GetUidInt()==cssh->GetUidInt()) {
01041 MSG("EventSS",Msg::kVerbose) << "SubShower (old,new) UID "
01042 << ssl->GetUidInt()
01043 << " " << cssh->GetUidInt() << endl;
01044 showerSR->RemoveSubShower(ssl);
01045 CandSubShowerSRHandle *newss = cssh->DupHandle();
01046 newss->SetCandSlice(cssh->GetCandSlice());
01047 newss->AddDaughterLink(*closeststrip);
01048 newss->SetEnergy(GeVperMIP*Double_t(calibrator.GetMIP(closeststrip->
01049 GetCharge(CalDigitType::kSigCorr))) + newss->GetEnergy());
01050 subshowerlist->RemoveDaughter(cssh);
01051 subshowerlist->AddDaughterLink(*newss);
01052 showerSR->AddSubShower(newss);
01053 delete newss;
01054 break;
01055 }
01056 }
01057 }
01058 else {
01059 MSG("EventSS",Msg::kWarning) << "Strip added to shower but not to subshower!"
01060 << " StripDistance = " << StripDistance
01061 << endl;
01062 MSG("EventSS",Msg::kWarning) << "Number of subshowers in shower = "
01063 << showerSR->GetLastSubShower()+1
01064 << endl;
01065 for(int i=0;i<showerSR->GetLastSubShower()+1;i++){
01066 CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
01067 (showerSR->GetSubShower(i));
01068 MSG("EventSS",Msg::kWarning) << "SS " << i << " UID = "
01069 << ssl->GetUidInt() << endl;
01070 }
01071 }
01072 }
01073 }
|
|
||||||||||||||||||||
|
this method is called for slices in which no shower or track was previously found. We see whether there is a 3D cluster in this slice which could be an event Definition at line 1389 of file AlgEventSSList.cxx. References abs(), CandHandle::AddDaughterLink(), CandHandle::FindDaughter(), Registry::Get(), AlgFactory::GetAlgHandle(), CandStripHandle::GetBegTime(), CandHandle::GetCandRecord(), CandContext::GetCandRecord(), CandStripHandle::GetCharge(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), VldContext::GetDetector(), Registry::GetDouble(), AlgFactory::GetInstance(), Registry::GetInt(), CandContext::GetMom(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandStripHandle::GetTPos(), RecMinos::GetVldContext(), CandHandle::GetVldContext(), CandRecoHandle::GetVtxT(), CandRecoHandle::GetVtxU(), CandRecoHandle::GetVtxV(), CandRecoHandle::GetVtxZ(), CandShowerSRHandle::IsUnphysical(), CandEvent::MakeCandidate(), CandShowerSR::MakeCandidate(), CandSubShowerSR::MakeCandidate(), MSG, CandHandle::SetCandRecord(), CandEventHandle::SetCandSlice(), CandRecoHandle::SetCandSlice(), CandSubShowerSRHandle::SetClusterID(), and CandEventHandle::SetPrimaryShower(). Referenced by RunAlg(). 01393 {
01394 assert(cx.GetDataIn());
01395 MSG("EventSS",Msg::kDebug) <<" Building Event - no showers or tracks "
01396 << endl;
01397 if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) {
01398 MSG("EventSS",Msg::kDebug)<<" dataIn inherits from TObjArray " << endl;
01399 return;
01400 }
01401
01402 Double_t pShwXTalkFracCut = ac.GetDouble("ShwXTalkFracCut");
01403 Double_t pShwXTalkDef = ac.GetDouble("ShwXTalkDef");
01404
01405 const CandSliceListHandle *slicelist = 0;
01406 CandShowerListHandle *showerlist = 0;
01407 CandSubShowerSRListHandle *subshowerlist = 0;
01408 const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn());
01409 for (Int_t i=0; i<=cxin->GetLast(); i++) {
01410 TObject *tobj = cxin->At(i);
01411 if (tobj->InheritsFrom("CandSliceListHandle")) {
01412 slicelist = dynamic_cast<CandSliceListHandle*>(tobj);
01413 }
01414 // Retrieve persistant modifiable CandShowerList
01415 if (tobj->InheritsFrom("CandShowerListHandle")) {
01416 showerlist = dynamic_cast<CandShowerListHandle*>(tobj);
01417 }
01418 if (tobj->InheritsFrom("CandSubShowerSRListHandle")) {
01419 subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj);
01420 }
01421 }
01422
01423 if (!showerlist || !slicelist || !subshowerlist) {
01424 MSG("EventSS",Msg::kDebug)
01425 << "Missing either slicelist, showerlist or subshowerlist" << endl;
01426 return;
01427 }
01428 Int_t detectorType = slicelist->GetVldContext()->GetDetector();
01429
01430 CandContext cxx(this,cx.GetMom());
01431 AlgFactory &af = AlgFactory::GetInstance();
01432
01433 const char *pEventAlgConfig = 0;
01434 ac.Get("EventAlgConfig",pEventAlgConfig);
01435 AlgHandle ah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig);
01436
01437 Double_t pHitAssocDPlane = ac.GetInt("HitAssocDPlane");
01438 Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0");
01439 Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1");
01440 Double_t pHitAssocMaxDist2 = ac.GetDouble("HitAssocMaxDist2");
01441 Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm");
01442 Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm");
01443 Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm");
01444 Double_t pShwShwDtpos2 = ac.GetDouble("ShwShwDtpos2");
01445 Double_t pShwShwDz = ac.GetDouble("ShwShwDz");
01446 Double_t pShwShwDt = ac.GetDouble("ShwShwDt");
01447 Double_t pMinShwEFract = ac.GetDouble("MinShwEFract");
01448
01449
01450 const CandRecord *candrec = cx.GetCandRecord();
01451 assert(candrec);
01452 const VldContext *vldcptr = candrec->GetVldContext();
01453 assert(vldcptr);
01454 VldContext vldc = *vldcptr;
01455 UgliGeomHandle ugh(vldc);
01456
01457 // Iterate over strips , finding sets which are 'in time', and within
01458 // a given distance of the nearest neighbor in the set, in both views.
01459 // Each of these sets is then used to construct a shower, which becomes
01460 // the primary shower of a new event. We continue until an iteration
01461 // in which no new event is created.
01462
01463 // fill array of unique unassociated strips
01464 TObjArray unassociated;
01465 TIter stripItr(slice->GetDaughterIterator());
01466 while (CandStripHandle *strip =
01467 dynamic_cast<CandStripHandle*>(stripItr())) {
01468 Bool_t found=false;
01469 for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) {
01470 CandStripHandle *strip0 =
01471 dynamic_cast<CandStripHandle*>(unassociated.At(i));
01472 if (*strip == *strip0) {
01473 found = true;
01474 }
01475 }
01476 if (!found) {
01477 if(!(detectorType==Detector::kNear && strip->GetPlane()>120)){
01478 unassociated.Add(strip);
01479 }
01480 }
01481 }
01482
01483 // loop until no new matches are found in the loop over unassoc hits.
01484 Bool_t found=false;
01485 while (unassociated.GetLast()>0 && !found) {
01486 MSG("EventSS",Msg::kDebug)<< " # unassoc "
01487 << unassociated.GetLast()+1 << endl;
01488 Bool_t firstu=true;
01489 Bool_t firstv=true;
01490 TObjArray neweventu;
01491 TObjArray neweventv;
01492 Bool_t addedstrip = true;
01493 Float_t totcharge=0;
01494 while (addedstrip) {
01495 addedstrip=false;
01496 // loop over unassoc strips.
01497 for (Int_t i=0; i<=unassociated.GetLast(); i++) {
01498 CandStripHandle *strip =
01499 dynamic_cast<CandStripHandle*>(unassociated.At(i));
01500 switch (strip->GetPlaneView()) {
01501 case PlaneView::kU:
01502 // strip is in U plane
01503 // if this is the first strip, add to newevent strip list
01504 if (firstu) {
01505 // if we don't have v strips yet, just add this strip to u list
01506 // seed the event
01507 if (firstv) {
01508 neweventu.Add(strip);
01509 MSG("EventSS",Msg::kDebug)<< " add first " << endl;
01510 totcharge +=strip->GetCharge();
01511 addedstrip = true;
01512 firstu=false;
01513 }
01514 // if we already have one or more v strips check that they are
01515 // in time, and agree in Z position.
01516 else {
01517 Bool_t found2=false;
01518 for (Int_t j=0; j<=neweventv.GetLast() && !found2;
01519 j++) {
01520 CandStripHandle *estrip =
01521 dynamic_cast<CandStripHandle*>(neweventv.At(j));
01522 Double_t dtime =
01523 strip->GetBegTime()-estrip->GetBegTime();
01524 MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9
01525 << " " <<pHitAssocDTime0*1e9
01526 << "/" << pHitAssocDTime1*1e9
01527 << " plane diff "
01528 << abs(strip->GetPlane() -
01529 estrip->GetPlane())
01530 << "/" << pHitAssocDPlane << endl;
01531 if (dtime>pHitAssocDTime0 &&
01532 dtime<pHitAssocDTime1 &&
01533 abs(strip->GetPlane()-estrip->GetPlane()) <=
01534 pHitAssocDPlane) {
01535 MSG("EventSS",Msg::kVerbose)<< " add " << endl;
01536 // match is found, add this strip to array
01537 neweventu.Add(strip);
01538 totcharge +=strip->GetCharge();
01539 addedstrip = true;
01540 firstu=0;
01541 found2=true;
01542 } // if match found
01543 } // for loop over V strips in new event
01544 } // end not first V strip
01545 } // end first U strip
01546 else {
01547 // if we already have u strips, make sure that the one we
01548 // add is adjaccent, using distance criteria based
01549 // on plane, time, and transverse postion separation
01550 Bool_t found2=false;
01551 for (Int_t j=0; j<=neweventu.GetLast() && !found2; j++) {
01552 CandStripHandle *estrip =
01553 dynamic_cast<CandStripHandle*>(neweventu.At(j));
01554 Double_t dtime =
01555 strip->GetBegTime()-estrip->GetBegTime();
01556 MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9
01557 << " " <<pHitAssocDTime0*1e9
01558 << "/" << pHitAssocDTime1*1e9
01559 << " plane diff "
01560 << abs(strip->GetPlane() -
01561 estrip->GetPlane())
01562 << " dtpos " << (strip->GetTPos() -
01563 estrip->GetTPos())
01564 << endl;
01565 if (dtime>pHitAssocDTime0 &&
01566 dtime<pHitAssocDTime1) {
01567 Double_t dplane =
01568 (Double_t)(strip->GetPlane()-estrip->GetPlane());
01569 Double_t dtpos = strip->GetTPos()-estrip->GetTPos();
01570 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01571 (pHitAssocPParm*dtpos*dtpos) +
01572 (pHitAssocTParm*dtime*dtime);
01573 // if dist2 2 parameter meets match criteria
01574 // add strip to array
01575 if (dist2<pHitAssocMaxDist2*4) {
01576 MSG("EventSS",Msg::kVerbose)<< " add " << endl;
01577 neweventu.Add(strip);
01578 totcharge +=strip->GetCharge();
01579 addedstrip = true;
01580 firstu=0;
01581 found2=true;
01582 } // if positions match
01583 } // if in time
01584 } // end loop over event u strips
01585 } // end if already have u strip
01586 break;
01587
01588 case PlaneView::kV:
01589
01590 // if this is the first strip, add to newevent strip list
01591 if (firstv) {
01592 // strip is in V plane
01593 // if this is the first strip, add to newevent strip list
01594 if (firstu) {
01595 neweventv.Add(strip);
01596 totcharge +=strip->GetCharge();
01597 MSG("EventSS",Msg::kDebug)<< " add first" << endl;
01598 addedstrip = true;
01599 firstv = 0;
01600 }
01601 else {
01602 // if we already have one or more u strips check that they are
01603 // in time, and require agreement in Z position.
01604 Bool_t found2=false;
01605 for (Int_t j=0; j<=neweventu.GetLast() && !found2;
01606 j++) {
01607 CandStripHandle *estrip =
01608 dynamic_cast<CandStripHandle*>(neweventu.At(j));
01609 Double_t dtime =
01610 strip->GetBegTime()-estrip->GetBegTime();
01611 MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9
01612 << " " << pHitAssocDTime0*1e9
01613 << "/" << pHitAssocDTime1*1e9
01614 << " plane diff "
01615 << abs(strip->GetPlane() -
01616 estrip->GetPlane())
01617 << "/" << pHitAssocDPlane << endl;
01618 if (dtime>pHitAssocDTime0 &&
01619 dtime<pHitAssocDTime1 &&
01620 abs(strip->GetPlane()-estrip->GetPlane()) <=
01621 pHitAssocDPlane) {
01622 // if dist2 2 parameter meets match criteria,
01623 //add strip to array
01624 MSG("EventSS",Msg::kVerbose) << " add " << endl;
01625 neweventv.Add(strip);
01626 totcharge +=strip->GetCharge();
01627 addedstrip = true;
01628 firstv = false;
01629 found2 = true;
01630 }
01631 } // end loop over u strips in event
01632 } // end if already have u strips
01633 } // end if this is the first v strip
01634 else {
01635 // if we already have v strips, make sure that the one we
01636 // add is adjaccent, using same distance criteria base
01637 // on plane, time, and transverse postion separation
01638 Bool_t found2=false;
01639 for (Int_t j=0; j<=neweventv.GetLast() && !found2; j++) {
01640 CandStripHandle *estrip =
01641 dynamic_cast<CandStripHandle*>(neweventv.At(j));
01642 Double_t dtime = strip->GetBegTime() -
01643 estrip->GetBegTime();
01644 MSG("EventSS",Msg::kVerbose)<< " dtime " << dtime*1e9 << " "
01645 << pHitAssocDTime0*1e9 << "/"
01646 << pHitAssocDTime1*1e9
01647 << " plane diff "
01648 << abs(strip->GetPlane() -
01649 estrip->GetPlane())
01650 << " dtpos " << (strip->GetTPos() -
01651 estrip->GetTPos())
01652 << endl;
01653 if (dtime>pHitAssocDTime0 &&
01654 dtime<pHitAssocDTime1) {
01655 Double_t dplane =
01656 (Double_t)(strip->GetPlane()-estrip->GetPlane());
01657 Double_t dtpos = strip->GetTPos()-estrip->GetTPos();
01658 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01659 (pHitAssocPParm*dtpos*dtpos) +
01660 (pHitAssocTParm*dtime*dtime);
01661 // if dist2 2 parameter meets match criteria,
01662 // add strip to array
01663 if (dist2<pHitAssocMaxDist2*4) {
01664 MSG("EventSS",Msg::kVerbose) << " add " << endl;
01665 neweventv.Add(strip);
01666 totcharge +=strip->GetCharge();
01667 addedstrip = true;
01668 firstv = false;
01669 found2 = true;
01670 } // end if add strip
01671 } //end if in time
01672 } // end loop over v strips
01673 }
01674 break;
01675
01676 default:
01677 break;
01678 }
01679 }
01680
01681 // end of loop over hits - remove new 'associations' and compress
01682 for (Int_t i=0; i<=neweventu.GetLast(); i++) {
01683 CandStripHandle *estrip =
01684 dynamic_cast<CandStripHandle*>(neweventu.At(i));
01685 unassociated.Remove(estrip);
01686 }
01687 for (Int_t i=0; i<=neweventv.GetLast(); i++) {
01688 CandStripHandle *estrip =
01689 dynamic_cast<CandStripHandle*>(neweventv.At(i));
01690 unassociated.Remove(estrip);
01691 }
01692 unassociated.Compress();
01693 }
01694
01695 // if we have hits in both views, build an event from them
01696 // if pulse height>10 pes
01697 MSG("EventSS",Msg::kDebug)<< " considering event formation u/v: "
01698 << !firstu << "/" << !firstv << " charge "
01699 << totcharge << endl;
01700 if (!firstu && !firstv && totcharge>20.0) {
01701 found=true;
01702
01703 //construct CandSubShowerSRs
01704 cxx.SetDataIn(&neweventu);
01705 AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default");
01706 CandSubShowerSRHandle usubshowerhandle =
01707 CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01708 usubshowerhandle.SetCandSlice(slice);
01709 usubshowerhandle.SetClusterID(ClusterType::kUnknown);
01710 subshowerlist->AddDaughterLink(usubshowerhandle);
01711 cxx.SetDataIn(&neweventv);
01712 CandSubShowerSRHandle vsubshowerhandle =
01713 CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01714 vsubshowerhandle.SetCandSlice(slice);
01715 vsubshowerhandle.SetClusterID(ClusterType::kUnknown);
01716 subshowerlist->AddDaughterLink(vsubshowerhandle);
01717 CandHandle *usubshower=subshowerlist->FindDaughter(&usubshowerhandle);
01718 CandHandle *vsubshower=subshowerlist->FindDaughter(&vsubshowerhandle);
01719
01720 //build shower
01721 TObjArray newshower;
01722 newshower.Add(usubshower);
01723 newshower.Add(vsubshower);
01724 cxx.SetDataIn(&newshower);
01725 AlgHandle showerah = af.GetAlgHandle("AlgShowerSS","default");
01726 CandShowerSRHandle showerhandle =
01727 CandShowerSR::MakeCandidate(showerah,cxx);
01728 showerhandle.SetCandSlice(slice);
01729 showerhandle.SetCandRecord(slicelist->GetCandRecord());
01730
01731 Bool_t newShowerOverlapsOld=false;
01732 CandShowerHandleItr showerItr(showerlist->GetDaughterIterator());
01733 while (CandShowerHandle *shower =
01734 dynamic_cast<CandShowerHandle*>(showerItr())){
01735
01736 Double_t dz = showerhandle.GetVtxZ()-shower->GetVtxZ();
01737 Double_t du = showerhandle.GetVtxU()-shower->GetVtxU();
01738 Double_t dv = showerhandle.GetVtxV()-shower->GetVtxV();
01739 Double_t dt = showerhandle.GetVtxT()-shower->GetVtxT();
01740
01741 if( ( (du*du+dv*dv<pShwShwDtpos2 &&fabs(dz)<pShwShwDz) ||
01742 (du*du<pShwShwDtpos2/4. && dv*dv<pShwShwDtpos2*4. &&
01743 fabs(dz)<pShwShwDz/2.) ||
01744 (du*du<pShwShwDtpos2*4 && dv*dv<pShwShwDtpos2/4. &&
01745 fabs(dz)<pShwShwDz/2.) ) &&
01746 fabs(dt)<pShwShwDt ) {
01747 newShowerOverlapsOld = true;
01748 MSG("EventSS",Msg::kDebug)
01749 << " new shower overlaps with previous - don't make new event "
01750 << endl;
01751 break;
01752 }
01753 }
01754
01755 if(!newShowerOverlapsOld){
01756 if(!showerhandle.IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
01757 MSG("EventSS",Msg::kDebug)<< " Building Event! " << endl;
01758 showerlist->AddDaughterLink(showerhandle);
01759 TObjArray recolist;
01760 recolist.Add(&showerhandle);
01761 cxx.SetDataIn(&recolist);
01762 AlgHandle eventah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig);
01763 CandEventHandle eventhandle =
01764 CandEvent::MakeCandidate(eventah,cxx);
01765 eventhandle.SetCandSlice(slice);
01766 eventhandle.SetPrimaryShower(pMinShwEFract,pShwShwDz);
01767 ch.AddDaughterLink(eventhandle);
01768 }
01769 }
01770 }
01771 }
01772 }
|
|
||||||||||||||||||||||||||||
|
Definition at line 1077 of file AlgEventSSList.cxx. References abs(), CandHandle::AddDaughterLink(), CandEventHandle::AddShower(), CandHandle::FindDaughter(), Registry::Get(), AlgFactory::GetAlgHandle(), CandRecoHandle::GetBegPlane(), CandEventHandle::GetCandSlice(), CandStripHandle::GetCharge(), CandHandle::GetDaughterIterator(), Registry::GetDouble(), AlgFactory::GetInstance(), CandStripHandle::GetPlane(), CandSubShowerSRHandle::GetPlaneView(), CandStripHandle::GetPlaneView(), CandEventHandle::GetPrimaryTrack(), CandShowerSR::MakeCandidate(), CandSubShowerSR::MakeCandidate(), MSG, CandRecoHandle::SetCandSlice(), CandSubShowerSRHandle::SetClusterID(), CandContext::SetDataIn(), and CandEventHandle::SetPrimaryShower(). Referenced by RunAlg(). 01082 {
01083 // need to create a shower -- include strips from first plane in each
01084 // view plus this previously unassociated hit to generate a 3D object
01085 MSG("EventSS",Msg::kVerbose) << "In CreatePrimaryShower" << endl;
01086 // Get singleton instance of AlgFactory
01087 const char *pEventAlgConfig = 0;
01088 ac.Get("EventAlgConfig",pEventAlgConfig);
01089
01090 //Double_t shwshwdz = ac.GetDouble("ShwShwDz");
01091 //Double_t minShwEFract = ac.GetDouble("MinShwEFract");
01092 Double_t pMinShwStripPE = ac.GetDouble("MinShwStripPE");
01093 Double_t pMaxNewShwLen = ac.GetDouble("MaxNewShwLen");
01094
01095 AlgFactory &af = AlgFactory::GetInstance();
01096
01097 CandTrackHandle *track = closestevent->GetPrimaryTrack();
01098 if (!showerlist) {
01099 MSG("EventSS",Msg::kError) <<"CandShowerListHandle missing\n" << endl;
01100 return;
01101 }
01102
01103 closestevent->AddDaughterLink(*closeststrip);
01104
01105 // create arrays of strips for each planeview to create
01106 // shower from
01107
01108 MSG("EventSS",Msg::kVerbose) << "creating shower \n";
01109 TObjArray ustriparray;
01110 TObjArray vstriparray;
01111 Int_t uplane=track->GetBegPlane(PlaneView::kU);
01112 Int_t vplane=track->GetBegPlane(PlaneView::kV);
01113 // seed this shower with strips in most upstream track hits
01114 // in each view, along with hits above a charge threshold and
01115 // within pMaxNewShwLen of the track start
01116 TIter stripItr(track->GetDaughterIterator());
01117 while (CandStripHandle *tstrip =
01118 dynamic_cast<CandStripHandle*>(stripItr())) {
01119 if ((tstrip->GetPlane() == uplane) ||
01120 (tstrip->GetPlaneView()==PlaneView::kU &&
01121 tstrip->GetCharge()>pMinShwStripPE &&
01122 abs(tstrip->GetPlane()-uplane)<pMaxNewShwLen)) {
01123 ustriparray.Add(tstrip);
01124 }
01125 if ((tstrip->GetPlane() == vplane) ||
01126 (tstrip->GetPlaneView()==PlaneView::kV &&
01127 tstrip->GetCharge()>pMinShwStripPE &&
01128 abs(tstrip->GetPlane()-vplane)<pMaxNewShwLen)) {
01129 vstriparray.Add(tstrip);
01130 }
01131 }
01132 // add unassoc strip to appropriate array
01133 switch (closeststrip->GetPlaneView()) {
01134 case PlaneView::kU:
01135 ustriparray.Add(closeststrip);
01136 break;
01137 case PlaneView::kV:
01138 vstriparray.Add(closeststrip);
01139 break;
01140 default:
01141 break;
01142 }
01143 // build shower from U/V clusters or subshowers
01144 if(ustriparray.GetLast()>=0 && vstriparray.GetLast()>=0){
01145 if(showerlist->InheritsFrom("CandShowerSRListHandle") && subshowerlist){
01146 MSG("EventSS",Msg::kVerbose) << "U striparray has "
01147 << ustriparray.GetLast()+1 << " entries" << endl;
01148 MSG("EventSS",Msg::kVerbose) << "V striparray has "
01149 << vstriparray.GetLast()+1 << " entries" << endl;
01150
01151 // now construct CandSubShowers
01152 AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default");
01153 cxx.SetDataIn(&ustriparray);
01154 CandSubShowerSRHandle usubshowerhandle =
01155 CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01156 usubshowerhandle.SetCandSlice(closestevent->GetCandSlice());
01157 usubshowerhandle.SetClusterID(ClusterType::kUnknown);
01158 subshowerlist->AddDaughterLink(usubshowerhandle);
01159 MSG("EventSS",Msg::kVerbose) << "U subshower view: "
01160 << usubshowerhandle.GetPlaneView() <<endl;
01161
01162 cxx.SetDataIn(&vstriparray);
01163 CandSubShowerSRHandle vsubshowerhandle =
01164 CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01165 vsubshowerhandle.SetCandSlice(closestevent->GetCandSlice());
01166 vsubshowerhandle.SetClusterID(ClusterType::kUnknown);
01167 subshowerlist->AddDaughterLink(vsubshowerhandle);
01168 MSG("EventSS",Msg::kVerbose) << "V subshower view: "
01169 << vsubshowerhandle.GetPlaneView() <<endl;
01170
01171 //add subshowers to list
01172 CandShowerSRListHandle* showerlistSR =
01173 dynamic_cast<CandShowerSRListHandle*>(showerlist);
01174 TObjArray newshower;
01175 CandSubShowerSRHandle * usubshower = dynamic_cast<CandSubShowerSRHandle*>
01176 (subshowerlist->FindDaughter(&usubshowerhandle));
01177
01178 CandSubShowerSRHandle * vsubshower = dynamic_cast<CandSubShowerSRHandle*>
01179 (subshowerlist->FindDaughter(&vsubshowerhandle));
01180 newshower.Add(usubshower);
01181 newshower.Add(vsubshower);
01182 cxx.SetDataIn(&newshower);
01183
01184 //get ShowerSS alghandle and make shower from subshowers
01185 AlgHandle showerahSS = af.GetAlgHandle("AlgShowerSS","default");
01186 CandShowerSRHandle showerhandleSR = CandShowerSR::MakeCandidate(showerahSS,cxx);
01187 showerhandleSR.SetCandSlice(closestevent->GetCandSlice());
01188 showerlistSR->AddDaughterLink(showerhandleSR);
01189 CandHandle *shw = showerlistSR->FindDaughter(&showerhandleSR);
01190 CandShowerHandle *shwh = dynamic_cast<CandShowerHandle*>(shw);
01191 closestevent->AddShower(shwh);
01192 // Change: cbs 9Aug06
01193 // NEED TO FORCE CREATED SHOWER TO BE PRIMARY
01194 // due to a change in assignment of primary shower, can now create
01195 // a shower here which is *not* called primary by SetPrimaryShower()
01196 // Causes multiple CreatePrimaryShower calls in event builder and
01197 // large numbers of reco'd showers per events.
01198 // SetPrimaryShowers is called at the end of RunAlg so a correct
01199 // assignment will be made later
01200 //closestevent->SetPrimaryShower(minShwEFract,shwshwdz);
01201 closestevent->SetPrimaryShower(shwh);
01202 }
01203 }
01204 }
|
|
||||||||||||||||||||
|
Definition at line 1890 of file AlgEventSSList.cxx. References CandRecoHandle::GetBegPlane(), CandStripHandle::GetBegTime(), CandHandle::GetDaughterIterator(), Registry::GetDouble(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandEventHandle::GetPrimaryShower(), CandEventHandle::GetPrimaryTrack(), CandStripHandle::GetStrip(), CandStripHandle::GetTPos(), CandRecoHandle::GetVtxT(), CandRecoHandle::GetVtxZ(), and MSG. Referenced by RunAlg(). 01890 {
01891 // grab all needed algorithm parameters
01892
01893 Double_t pShwTrkDz = ac.GetDouble("ShwTrkDz");
01894 Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0");
01895 Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1");
01896 Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm");
01897 Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm");
01898 Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm");
01899
01900 // construct the map, first looping over strips
01901 for (Int_t i=0; i<=unassociated.GetLast(); i++) {
01902 CandStripHandle *strip =
01903 dynamic_cast<CandStripHandle*>(unassociated.At(i));
01904 TIter eventItr(ch.GetDaughterIterator());
01905 // inner loop is over events
01906 while (CandEventHandle *event =
01907 dynamic_cast<CandEventHandle*>(eventItr())) {
01908
01909 // calculate separation in space/time
01910 Double_t dtime = (strip->GetBegTime()-event->GetVtxT());
01911 Double_t bestdist2=-1.;
01912 Double_t bestdplane=-1;
01913 Double_t bestdtpos=-1;
01914 Double_t bestdtime=-1;
01915 dist2map[i][event] = 999999999.;
01916 MSG("EventSS",Msg::kVerbose) << i << " " << strip->GetPlane() << " "
01917 << strip->GetStrip() << " stp beg time "
01918 << strip->GetBegTime() << " evt vtx time "
01919 << event->GetVtxT()
01920 << " dtime "
01921 << dtime*1e9 << "/" << pHitAssocDTime0*1e9
01922 << "/" << pHitAssocDTime1*1e9 << endl;
01923 // if there is a rough time match proceed
01924 if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1) {
01925
01926 Float_t vertexsep = 0;
01927 // If event has both shower
01928 // and track calc separation between the two. If this
01929 // separation is greater
01930 // than pSHwTrkDz, we won't bother adding the strip to the
01931 // primary shower, and will instead consider whether
01932 // to create a new shower at track vertex.
01933 if (event->GetPrimaryShower() && event->GetPrimaryTrack())
01934 vertexsep = fabs(event->GetPrimaryShower()->GetVtxZ() -
01935 event->GetPrimaryTrack()->GetVtxZ());
01936
01937 // the event has a shower, and either no track, or track and
01938 // shower vertices are in agreement. In this case, we
01939 // consider adding this strip to this shower, so calc separation.
01940 if (event->GetPrimaryShower() &&
01941 (!event->GetPrimaryTrack() || vertexsep<pShwTrkDz)) {
01942
01943 // create iterator of shower daughter strips in same view
01944 // as unassoc strip
01945 CandStripHandleItr shwstripItr(event->GetPrimaryShower()->
01946 GetDaughterIterator());
01947 CandStripHandleKeyFunc *shwstripKf = shwstripItr.CreateKeyFunc();
01948 shwstripKf->SetFun(CandStripHandle::KeyFromView);
01949 shwstripItr.GetSet()->AdoptSortKeyFunc(shwstripKf);
01950 shwstripKf = 0;
01951 shwstripItr.GetSet()->Slice();
01952 shwstripItr.GetSet()->Slice(strip->GetPlaneView(),
01953 strip->GetPlaneView());
01954
01955 // iterate over the daughter list, and find minimum distane
01956 // between daughter and unassoc strip
01957 while (CandStripHandle *shwstrip =
01958 dynamic_cast<CandStripHandle*>(shwstripItr())) {
01959 Double_t dplane = (Double_t)(strip->GetPlane()-shwstrip->GetPlane());
01960 Double_t dtpos = strip->GetTPos()-shwstrip->GetTPos();
01961 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01962 (pHitAssocPParm*dtpos*dtpos) +
01963 (pHitAssocTParm*dtime*dtime);
01964 if (bestdist2<0. || dist2<bestdist2) {
01965 bestdist2 = dist2;
01966 bestdplane=dplane;
01967 bestdtpos=dtpos;
01968 bestdtime=dtime;
01969 }
01970 }
01971 MSG("EventSS",Msg::kVerbose)
01972 << "primary shower." << " dplane:" << bestdplane <<" dtpos:"
01973 << bestdtpos <<" dtime:" << bestdtime <<" dist2:" << bestdist2<< "\n";
01974 dist2map[i][event] = bestdist2;
01975 } // end if event has shower, etc.
01976 // else if criteria above was not satisfied, consider whether
01977 // strip is near track vertex, in which case we start new shower
01978 else if (event->GetPrimaryTrack()) {
01979
01980 // calc distance to vertex
01981 Double_t dplane = (Double_t)(strip->GetPlane() -
01982 event->GetPrimaryTrack()->
01983 GetBegPlane());
01984 Double_t dtpos=0;
01985 switch (strip->GetPlaneView()) {
01986 case PlaneView::kU:
01987 dtpos=strip->GetTPos()-event->GetPrimaryTrack()->
01988 GetU(event->GetPrimaryTrack()->GetBegPlane());
01989 break;
01990 case PlaneView::kV:
01991 dtpos=strip->GetTPos()-event->GetPrimaryTrack()->
01992 GetV(event->GetPrimaryTrack()->GetBegPlane());
01993 break;
01994 default:
01995 break;
01996 }
01997 Double_t dist2 = ( (pHitAssocZParm*dplane*dplane) +
01998 (pHitAssocPParm*dtpos*dtpos) +
01999 (pHitAssocTParm*dtime*dtime) );
02000 dist2map[i][event] = dist2;
02001 MSG("EventSS",Msg::kVerbose)
02002 << "primary track. begin:" << " dplane:" << dplane
02003 << " dtpos:" << dtpos <<" dtime:" << dtime <<" dist2:"
02004 << dist2 << "\n";
02005 } // end if primary track
02006 } // end if time match OK
02007 } // end loop over events
02008 } // end loop over unassoc strips
02009 }
|
|
||||||||||||||||||||
|
Definition at line 1776 of file AlgEventSSList.cxx. References CandRecoHandle::GetCandSlice(), CandHandle::GetNDaughters(), CandHandle::IsCloneOf(), and MSG. Referenced by RunAlg(). 01777 {
01778
01779 if (showerItr) {
01780 showerItr->Reset();
01781 while (CandShowerHandle *shower =
01782 dynamic_cast<CandShowerHandle*>((*showerItr)())) {
01783 if (shower->GetCandSlice()) {
01784 if (slice->IsCloneOf(*(shower->GetCandSlice()))) {
01785 recolist.Add(shower);
01786 }
01787 }
01788 else {
01789 MSG("EventSS", Msg::kError)
01790 << "Shower doesn't contain a slice. Not added to recolist."
01791 << endl;
01792 }
01793 }
01794 }
01795
01796 if (trackItr) {
01797 trackItr->Reset();
01798 while (CandTrackHandle *track =
01799 dynamic_cast<CandTrackHandle*>((*trackItr)())) {
01800 if (track->GetCandSlice()) {
01801 if (slice->IsCloneOf(*track->GetCandSlice())) {
01802 if(track->GetNDaughters()==0 &&
01803 track->InheritsFrom("CandFitTrackHandle") &&
01804 dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack()){
01805 track = dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack();
01806 MSG("EventSS", Msg::kDebug)<< "Finder Track being used " << endl;
01807 }
01808 recolist.Add(track);
01809 }
01810 }
01811 else {
01812 MSG("EventSS", Msg::kError)
01813 << "Track doesn't contain a slice. Not added to recolist."
01814 << endl;
01815 }
01816 }
01817 }
01818 }
|
|
||||||||||||||||
|
Definition at line 1843 of file AlgEventSSList.cxx. References CandHandle::FindDaughter(), CandEventHandle::GetCandSlice(), CandHandle::GetDaughterIterator(), and CandHandle::IsCloneOf(). Referenced by RunAlg(). 01845 {
01846
01847 while (CandSliceHandle *slice =
01848 dynamic_cast<CandSliceHandle*>(sliceItr())) {
01849 TIter stripItr(slice->GetDaughterIterator());
01850 while (CandStripHandle *strip =
01851 dynamic_cast<CandStripHandle*>(stripItr())) {
01852 Bool_t found=false;
01853 // check for duplicates
01854 for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) {
01855 CandStripHandle *strip0 =
01856 dynamic_cast<CandStripHandle*>(unassociated.At(i));
01857 if (*strip == *strip0) {
01858 found = true;
01859 break;
01860 }
01861 }
01862 // if this strip not already in unnassociated list, then
01863 // loop over objects in this slice, and check against all daughter strips
01864 if (!found) {
01865 if (eventItr) {
01866 eventItr->Reset();
01867 while (CandEventHandle *event =
01868 dynamic_cast<CandEventHandle*>((*eventItr)())) {
01869 if (event->GetCandSlice()) {
01870 if (slice->IsCloneOf(*(event->GetCandSlice()))) {
01871 if (event->FindDaughter(strip)) {
01872 found = true;
01873 break;
01874 }
01875 }
01876 }
01877 }
01878 }
01879 }
01880 // this strip is not in any event, add to unassociated list
01881 if (!found) {
01882 unassociated.Add(strip);
01883 }
01884 }
01885 }
01886 }
|
|
||||||||||||||||
|
This method fills a cluster list for every shower in the snarl, and then reruns the shower algorithm. Following this, it tracks the primary track in each event through the primary shower, and subtracts 1 MIP for each strip passed through. Definition at line 1235 of file AlgEventSSList.cxx. References CandHandle::AddDaughterLink(), CandShowerHandle::CalibrateEnergy(), CandSubShowerSRHandle::DupHandle(), Registry::Get(), AlgHandle::GetAlgConfig(), AlgFactory::GetAlgHandle(), CandEventHandle::GetCandSlice(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), AlgFactory::GetInstance(), CandEventHandle::GetLastShower(), CandShowerSRHandle::GetLastSubShower(), CandEventHandle::GetPrimaryTrack(), CandEventHandle::GetShowerWritable(), CandShowerSRHandle::GetSubShower(), CandHandle::GetUidInt(), CandHandle::IsCloneOf(), CandHandle::IsEqual(), MSG, CandHandle::RemoveDaughter(), AlgHandle::RunAlg(), CandRecoHandle::SetCandSlice(), and CandContext::SetDataIn(). Referenced by RunAlg(). 01238 {
01239 MSG("EventSS",Msg::kVerbose) << "In ReConstructShowers: " << endl;
01240 // Retrieve persistant modifiable CandShowerList
01241 CandShowerListHandle *showerlist = const_cast<CandShowerListHandle*>
01242 (dynamic_cast<const CandShowerListHandle*>(cxx.GetDataIn()));
01243
01244 if(!showerlist->InheritsFrom("CandShowerSRListHandle")) {
01245 MSG("EventSS",Msg::kWarning)
01246 << "ListHandle is not a CandShowerSRListHandle."
01247 << " ShowerList not processed in ReConstructShowers." << endl;
01248 return;
01249 }
01250
01251 // Get singleton instance of AlgFactory
01252 AlgFactory &af = AlgFactory::GetInstance();
01253
01254 // Get an AlgHandle to AlgShowerSR with default AlgConfig
01255 const char *pEventAlgConfig = 0;
01256 ac.Get("EventAlgConfig",pEventAlgConfig);
01257
01258 AlgHandle ahss = af.GetAlgHandle("AlgShowerSS","default");
01259 AlgConfig acss = ahss.GetAlgConfig();
01260 Int_t eventCounter = 0;
01261 TIter evItr(ch.GetDaughterIterator());
01262 while (CandEventHandle *ev =
01263 dynamic_cast<CandEventHandle*>(evItr())) { //event loop
01264 MSG("EventSS",Msg::kVerbose) << "event = " << eventCounter << endl;
01265 for (Int_t ishower=0; ishower<ev->GetLastShower()+1; ishower++) { //shower loop
01266 MSG("EventSS",Msg::kVerbose) << "shower = " << ishower << endl;
01267 CandShowerHandle *shower =
01268 dynamic_cast<CandShowerHandle*>(ev->GetShowerWritable(ishower));
01269 CandShowerSRHandle *showerSR = 0;
01270 if(shower->InheritsFrom("CandShowerSRHandle")) {
01271 showerSR = dynamic_cast<CandShowerSRHandle*>(shower);
01272 }
01273 MSG("EventSS",Msg::kVerbose) << "Shower UID = " << showerSR->GetUidInt() << endl;
01274 //assume subshower chain
01275 if(showerSR) {
01276 TObjArray newshower;
01277 // fill array of subshowers for this shower
01278 for (Int_t isubshower=0;
01279 isubshower<showerSR->GetLastSubShower()+1;
01280 isubshower++) {
01281 MSG("EventSS",Msg::kVerbose) << "subshower = " << isubshower << endl;
01282 const CandSubShowerSRHandle *ssh =
01283 showerSR->GetSubShower(isubshower);
01284 if (ssh) {
01285 // New handle on heap
01286 CandSubShowerSRHandle *sshd = ssh->DupHandle();
01287 // Add new heap-based CandHandle to TObjArray
01288 newshower.Add(sshd);
01289 MSG("EventSS",Msg::kDebug) << "Uid's (old,new) "
01290 << ssh->GetUidInt() << " "
01291 << sshd->GetUidInt() << endl;
01292 }
01293 }
01294 // re-run algorithm for this shower
01295 if (newshower.GetEntriesFast() > 0) {
01296 cxx.SetDataIn(&newshower);
01297 ahss.RunAlg(*showerSR,cxx);
01298 }
01299 else {
01300 MSG("EventSS", Msg::kError)
01301 << "Attempt to pass empty TObjArray newshower to AlgShowerSS"
01302 << " in CandContext. AlgShowerSS::RunAlg() call ignored."
01303 << endl;
01304 }
01305 // Delete heap-based CandHandles in TObjArray
01306 newshower.Delete();
01307 }
01308 else {
01309 MSG("EventSS",Msg::kWarning)
01310 << "Handle is not a CandShowerSRHandle."
01311 << " Shower not processed in ReConstructShowers."
01312 << endl;
01313 continue;
01314 }
01315 MSG("EventSS",Msg::kVerbose) << "Modified shower UID = "
01316 << showerSR->GetUidInt() << endl;
01317 // now loop through strips in this shower, and determine whether primary
01318 // track in this event intercepts this strip. If so, subtract 1 mip
01319 // from shower energy.
01320
01321 // if event has primary track, extrapolate through shower
01322 // and subtract energy loss to recalibrate
01323 CandTrackHandle *primaryTrack = ev->GetPrimaryTrack();
01324 if (primaryTrack) {
01325 shower->CalibrateEnergy(primaryTrack,acss);
01326 }
01327
01328 // Update modified showers in persistent modifiable CandShowerList.
01329 // CandHandle could potentially supply a ReplaceDaughter() function.
01330 // Interface of ReplaceDaughter() would be debatable.
01331 if (showerlist) {
01332 TIter shiter(showerlist->GetDaughterIterator());
01333 CandShowerHandle *target;
01334 Bool_t found = kFALSE;
01335 while (!found &&
01336 (target = dynamic_cast<CandShowerHandle*>(shiter()))) {
01337 if (shower->IsCloneOf(*target)) { // Tests clone or Cand ==
01338 found = kTRUE;
01339
01340 // Replace target only if shower is a modified version of target
01341 if (!(shower->IsEqual(target))) { // CandBase inequality
01342 shower->SetCandSlice(ev->GetCandSlice());
01343 if (!showerlist->RemoveDaughter(target)) { // Failure
01344 MSG("EventSS",Msg::kError)
01345 << endl
01346 << "Failure of ShowerList daughter removal " << endl
01347 << "attempted during replacement by modified Shower."
01348 << endl << "Will result in double counted Shower."
01349 << endl;
01350 }
01351 showerlist->AddDaughterLink(*shower);
01352 }
01353 }
01354 }
01355
01356 // Add new (unfound) shower to persistent modifiable CandShowerList.
01357 if (!found){
01358 shower->SetCandSlice(ev->GetCandSlice());
01359 showerlist->AddDaughterLink(*shower);
01360 }
01361 }
01362 }
01363 eventCounter++;
01364 }
01365 }
|
|
||||||||||||||||||||||||
|
Definition at line 1207 of file AlgEventSSList.cxx. References CandStripHandle::GetBegTime(), Registry::GetDouble(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandStripHandle::GetTPos(), and CandRecoHandle::GetVtxT(). Referenced by RunAlg(). 01208 {
01209 Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm");
01210 Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm");
01211 Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm");
01212
01213 for (Int_t i=0; i<=unassociated.GetLast(); i++) {
01214 CandStripHandle *strip =
01215 dynamic_cast<CandStripHandle*>(unassociated.At(i));
01216 if(strip->GetPlaneView()==closeststrip->GetPlaneView()){
01217 Double_t dtime = strip->GetBegTime()-closestevent->GetVtxT();
01218 Double_t dplane =
01219 (Double_t)(strip->GetPlane()-closeststrip->GetPlane());
01220 Double_t dtpos = strip->GetTPos()-closeststrip->GetTPos();
01221 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01222 (pHitAssocPParm*dtpos*dtpos) +
01223 (pHitAssocTParm*dtime*dtime);
01224 if (dist2<dist2map[i][closestevent]) {
01225 dist2map[i][closestevent] = dist2;
01226 }
01227 }
01228 } // loop over unassoc hits
01229 }
|
|
||||||||||||||||
|
Implements AlgBase. Definition at line 72 of file AlgEventSSList.cxx. References abs(), CandHandle::AddDaughterLink(), AddStripToEvent(), CandShowerSRHandle::BelongsWithShower(), CandTrackHandle::BelongsWithTrack(), CandShowerSRHandle::BelongsWithTrack(), BuildEventFromUnassoc(), CandShowerSRHandle::BuriedTrack(), CreatePrimaryShower(), FillDist2Map(), FillRecoList(), FindUnassociated(), Registry::Get(), AlgFactory::GetAlgHandle(), CandRecoHandle::GetBegPlane(), CandStripHandle::GetBegTime(), CandContext::GetCandRecord(), CandRecoHandle::GetCharge(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), CandRecoHandle::GetDirCosU(), CandRecoHandle::GetDirCosV(), Registry::GetDouble(), CandRecoHandle::GetEndPlane(), AlgFactory::GetInstance(), Registry::GetInt(), CandEventHandle::GetLastShower(), CandContext::GetMom(), CandHandle::GetNDaughters(), CandRecoHandle::GetNStrip(), CandStripHandle::GetPlane(), CandEventHandle::GetPrimaryShower(), CandEventHandle::GetPrimaryTrack(), CandEventHandle::GetShower(), CandEventHandle::GetShowerWritable(), CandStripHandle::GetTPos(), RecMinos::GetVldContext(), CandRecoHandle::GetVtxPlane(), CandRecoHandle::GetVtxT(), CandHandle::IsCloneOf(), CandTrackHandle::IsUnphysical(), CandShowerSRHandle::IsUnphysical(), CandEvent::MakeCandidate(), max, min, MSG, ChopHelp::nchop, PITCHINMETRES, ChopHelp::Print(), ReConstructShowers(), ReFillDist2Map(), CandEventHandle::SetCandSlice(), and SetPrimaryShowers(). 00073 {
00074
00075 assert(cx.GetDataIn());
00076 if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) return;
00077
00078 CandContext cxx(this,cx.GetMom());
00079 AlgFactory &af = AlgFactory::GetInstance();
00080 const char *pEventAlgConfig = 0;
00081 ac.Get("EventAlgConfig",pEventAlgConfig);
00082 AlgHandle ah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig);
00083 // grab all needed algorithm parameters
00084
00085 //generic parameters:
00086 Double_t pTrkTrkDtpos2 = ac.GetDouble("TrkTrkDtpos2");
00087 Double_t pTrkTrkDz = ac.GetDouble("TrkTrkDz");
00088 Double_t pTrkTrkDt = ac.GetDouble("TrkTrkDt");
00089 Double_t pShwTrkDtpos2 = ac.GetDouble("ShwTrkDtpos2");
00090 Double_t pShwTrkDz = ac.GetDouble("ShwTrkDz");
00091 Double_t pShwTrkDt = ac.GetDouble("ShwTrkDt");
00092 Double_t pShwShwDtpos2 = ac.GetDouble("ShwShwDtpos2");
00093 Double_t pShwShwDz = ac.GetDouble("ShwShwDz");
00094 Double_t pShwShwDt = ac.GetDouble("ShwShwDt");
00095 Double_t pMaxNewShwLen = ac.GetDouble("MaxNewShwLen");
00096 Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0");
00097 Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1");
00098 Double_t pHitAssocMaxDist2 = ac.GetDouble("HitAssocMaxDist2");
00099 Double_t pMinTrackIsolation = ac.GetDouble("MinTrackIsolation");
00100
00101 //SS specific parameters
00102 Double_t pHardBuriedFrac = ac.GetDouble("HardBuriedFrac");
00103 Double_t pSoftBuriedFrac = ac.GetDouble("SoftBuriedFrac");
00104 Double_t pDCosUVCut = ac.GetDouble("DCosUVCut");
00105 Double_t pTrkGapFracCut = ac.GetDouble("TrkGapFracCut");
00106 Double_t pTrkAsymCut = ac.GetDouble("TrkAsymCut");
00107 Double_t pTrkXTalkFracCut = ac.GetDouble("TrkXTalkFracCut");
00108 Double_t pTrkXTalkDef = ac.GetDouble("TrkXTalkDef");
00109 Double_t pShwXTalkFracCut = ac.GetDouble("ShwXTalkFracCut");
00110 Double_t pShwXTalkDef = ac.GetDouble("ShwXTalkDef");
00111 Int_t useChopper = ac.GetInt("UseChopper");
00112
00113 MSG("EventSS",Msg::kDebug) << pHardBuriedFrac << " " << pSoftBuriedFrac << " "
00114 << pDCosUVCut << " " << pTrkGapFracCut << " "
00115 << pTrkAsymCut << " " << pTrkXTalkFracCut << " "
00116 << pTrkXTalkDef << " " << pShwXTalkFracCut << " "
00117 << pShwXTalkDef << endl;
00118
00119 // grab all needed input lists
00120 const CandRecord *candrec = cx.GetCandRecord();
00121 assert(candrec);
00122 const VldContext *vldcptr = candrec->GetVldContext();
00123 assert(vldcptr);
00124 VldContext vldc = *vldcptr;
00125 const CandSliceListHandle *slicelist = 0;
00126 CandShowerListHandle *showerlist = 0;
00127 const CandTrackListHandle *tracklist = 0;
00128 CandSubShowerSRListHandle *subshowerlist = 0;
00129 const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn());
00130 for (Int_t i=0; i<=cxin->GetLast(); i++) {
00131 TObject *tobj = cxin->At(i);
00132 if (tobj->InheritsFrom("CandSliceListHandle")) {
00133 slicelist = dynamic_cast<CandSliceListHandle*>(tobj);
00134 }
00135 if (tobj->InheritsFrom("CandShowerListHandle")) {
00136 showerlist = dynamic_cast<CandShowerListHandle*>(tobj);
00137 }
00138 if (tobj->InheritsFrom("CandTrackListHandle")) {
00139 tracklist = dynamic_cast<CandTrackListHandle*>(tobj);
00140 }
00141 if (tobj->InheritsFrom("CandSubShowerSRListHandle")) {
00142 subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj);
00143 }
00144 }
00145 if (!slicelist) {
00146 MSG("EventSS",Msg::kError) << "CandSliceListHandle missing\n";
00147 return;
00148 }
00149
00150 CandShowerHandleItr * showerItr=0;
00151 // create iterators for track, shower, and cluster lists, keying on slice
00152 if (showerlist && showerlist->GetNDaughters()>0) {
00153 showerItr = new CandShowerHandleItr(showerlist->GetDaughterIterator());
00154 CandShowerHandleKeyFunc *showerKf = showerItr->CreateKeyFunc();
00155 showerKf->SetFun(CandShowerHandle::KeyFromSlice);
00156 showerItr->GetSet()->AdoptSortKeyFunc(showerKf);
00157 showerKf = 0;
00158 }
00159
00160 CandTrackHandleItr * trackItr = 0;
00161 if (tracklist && tracklist->GetNDaughters()>0) {
00162 trackItr = new CandTrackHandleItr(tracklist->GetDaughterIterator());
00163 CandTrackHandleKeyFunc *trackKf = trackItr->CreateKeyFunc();
00164 trackKf->SetFun(CandTrackHandle::KeyFromSlice);
00165 trackItr->GetSet()->AdoptSortKeyFunc(trackKf);
00166 trackKf = 0;
00167 }
00168
00169 CandSubShowerSRHandleItr * subshowerItr = 0;
00170 if (subshowerlist && subshowerlist->GetNDaughters()>0) {
00171 subshowerItr = new CandSubShowerSRHandleItr(subshowerlist->
00172 GetDaughterIterator());
00173 CandSubShowerSRHandleKeyFunc *subshowerKf = subshowerItr->CreateKeyFunc();
00174 subshowerKf->SetFun(CandSubShowerSRHandle::KeyFromSlice);
00175 subshowerItr->GetSet()->AdoptSortKeyFunc(subshowerKf);
00176 subshowerKf = 0;
00177 }
00178
00179 Int_t nslice=0;
00180
00181 // loop over slices
00182 CandSliceHandleItr sliceItr(slicelist->GetDaughterIterator());
00183 while (CandSliceHandle *slice =
00184 dynamic_cast<CandSliceHandle*>(sliceItr())) { //loop over slices
00185
00186 ++nslice;
00187 MSG("EventSS", Msg::kDebug)
00188 << " ****** Slice " << nslice << " *******"
00189 << endl;
00190
00191 //TObjArray newevent;
00192 TObjArray recolist;
00193 TObjArray eventlist;
00194
00195 //construct list of all showers and tracks in this slice
00196 FillRecoList(slice,showerItr,trackItr,recolist);
00197
00198 //first try to remove tracks buried in showers
00199 std::vector<Bool_t> removeMe(recolist.GetLast()+1);
00200 for (Int_t i=0;i<=recolist.GetLast();i++) {
00201 removeMe[i] = false;
00202 if(recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00203 CandTrackHandle *track =
00204 dynamic_cast<CandTrackHandle*>(recolist.At(i));
00205 Bool_t removeTrack = false;
00206 for (Int_t j=0;j<=recolist.GetLast();j++) {
00207 if (recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00208 CandShowerSRHandle * shower =
00209 dynamic_cast<CandShowerSRHandle*>(recolist.At(j));
00210 if(shower->BuriedTrack(track,pMinTrackIsolation)) {
00211 removeTrack = true;
00212 break;
00213 }
00214 }
00215 }
00216 if(removeTrack){
00217 MSG("EventSS",Msg::kDebug)
00218 << " REMOVING TRACK" << endl;
00219 removeMe[i] = true;
00220 }
00221 }
00222 }
00223 for (Int_t i=0;i<=recolist.GetLast();i++) {
00224 if(removeMe[i]) recolist.RemoveAt(i);
00225 }
00226 recolist.Compress();
00227
00228 MSG("EventSS",Msg::kDebug) << "Number of reco objects in slice: "
00229 << recolist.GetLast()+1 << endl;
00230 //now make matrix to hold match results between all objects.
00231 //this also includes buried tracks that have not been removed
00232 //NOTE: matching requires both views match;
00233 // burial just counts up matching planes in either view;
00234 // One view being significantly buried may indicate a bad track, which
00235 // should be associated with the shower, not used to make a new event
00236 std::vector<std::vector<Int_t> > objectMatch(recolist.GetLast()+1,0);
00237 //initialize:
00238 for(int i=0;i<=recolist.GetLast();i++){
00239 for(int j=0;j<=recolist.GetLast();j++) {
00240 objectMatch[i].push_back(0);
00241 }
00242 }
00243
00244 //use this vector to keep track of objects already assigned
00245 std::vector<Bool_t> objectUsed(recolist.GetLast()+1,0);
00246
00247 //loop over all objects
00248 for(Int_t i=0;i<=recolist.GetLast();i++) { //loop over i
00249 objectUsed[i] = false; //initialization
00250
00251 CandShowerSRHandle *shower=0;
00252 CandTrackHandle *track=0;
00253 if (recolist.At(i)->InheritsFrom("CandShowerSRHandle"))
00254 shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00255 if (recolist.At(i)->InheritsFrom("CandTrackHandle"))
00256 track = dynamic_cast<CandTrackHandle*>(recolist.At(i));
00257
00258 //loop over all later objects only (since matching matrix is symmetric)
00259 for(Int_t j=i+1;j<=recolist.GetLast();j++) { //loop over j
00260
00261 MSG("EventSS",Msg::kDebug) << "Checking objectMatch between "
00262 << i << " and " << j << endl;
00263
00264 CandShowerSRHandle *shower2=0;
00265 CandTrackHandle *track2=0;
00266 if (recolist.At(j)->InheritsFrom("CandShowerSRHandle"))
00267 shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(j));
00268 if (recolist.At(j)->InheritsFrom("CandTrackHandle"))
00269 track2 = dynamic_cast<CandTrackHandle*>(recolist.At(j));
00270
00271 if(shower) {
00272 if(shower2) {
00273 if(!shower->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef) &&
00274 !shower2->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00275 if(shower->BelongsWithShower(shower2,ac,vldcptr,pShwShwDtpos2,
00276 pShwShwDz,pShwShwDt)) objectMatch[i][j] = 1;
00277 }
00278 objectMatch[j][i] = objectMatch[i][j];
00279 }
00280 else if(track2) {
00281 if(shower->BelongsWithTrack(track2,ac,vldcptr,pShwTrkDtpos2,
00282 pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1;
00283
00284 if(objectMatch[i][j]!=1) {
00285 //check track isn't buried:
00286 Int_t *stats = new Int_t[8];
00287 shower->BuriedTrack(track2,1000,stats); //1000 forces all calculations
00288 if(stats[3]>0) { //at least one buried track planes
00289 Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00290 MSG("EventSS",Msg::kDebug)
00291 << "BuriedTrack Stats: " << "\n"
00292 << "Number of track planes = " << nplanes << "\n"
00293 << "nSharedPlanes = " << stats[0] << "\n"
00294 << "nMissingIsoPlanes = " << stats[1] << "\n"
00295 << "nMissingSharedPlanes = " << stats[2] << "\n"
00296 << "nBuriedPlanes = " << stats[3] << "\n"
00297 << "nPhysBuriedPlanes = " << stats[4] << "\n"
00298 << "nTrkLike = " << stats[5] << "\n"
00299 << "nTrkLikeTopol = " << stats[6] << "\n"
00300 << "nXTalk = " << stats[7] << endl;
00301 nplanes -= (stats[1] + stats[2]); //subtract missing track planes
00302 //if more than half of the track planes are in the shower core: match
00303 if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac)
00304 objectMatch[i][j] = -1;
00305 //else if more than half the track planes are anywhere in
00306 //the shower or are x-talk like and are not "trklike": match
00307 else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) >
00308 pHardBuriedFrac) objectMatch[i][j] = -1;
00309 //check whether track passes through shower at high angle in either view
00310 //expect significant number of buried planes
00311 else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) &&
00312 ( TMath::Abs(track2->GetDirCosU())>pDCosUVCut ||
00313 TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) {
00314 objectMatch[i][j] = -1;
00315 }
00316 }
00317 delete [] stats;
00318 }
00319 objectMatch[j][i] = objectMatch[i][j];
00320 }
00321 }
00322
00323 if(track) {
00324 if(shower2) {
00325 if(shower2->BelongsWithTrack(track,ac,vldcptr,pShwTrkDtpos2,
00326 pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1;
00327 if(objectMatch[i][j]!=1) {
00328 //check track isn't buried:
00329 Int_t *stats = new Int_t[8];
00330 shower2->BuriedTrack(track,1000,stats); //1000 forces all calculations
00331 if(stats[3]>0) { //at least one buried track planes
00332 Int_t nplanes = track->GetEndPlane() - track->GetBegPlane() + 1;
00333 MSG("EventSS",Msg::kDebug)
00334 << "BuriedTrack Stats: " << "\n"
00335 << "Number of track planes = " << nplanes << "\n"
00336 << "nSharedPlanes = " << stats[0] << "\n"
00337 << "nMissingIsoPlanes = " << stats[1] << "\n"
00338 << "nMissingSharedPlanes = " << stats[2] << "\n"
00339 << "nBuriedPlanes = " << stats[3] << "\n"
00340 << "nPhysBuriedPlanes = " << stats[4] << "\n"
00341 << "nTrkLike = " << stats[5] << "\n"
00342 << "nTrkLikeTopol = " << stats[6] << "\n"
00343 << "nXTalk = " << stats[7] << endl;
00344 nplanes -= (stats[1] + stats[2]); //subtract missing track planes
00345 //if more than half of the track planes are in the shower core: match
00346 if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac)
00347 objectMatch[i][j] = -1;
00348 //else if more than half the track planes are anywhere in
00349 //the shower or are x-talk like and are not "trklike": match
00350 else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) >
00351 pHardBuriedFrac) objectMatch[i][j] = -1;
00352 //check whether track passes through shower at high angle in either view
00353 //expect significant number of buried planes
00354 else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) &&
00355 (TMath::Abs(track2->GetDirCosU())>pDCosUVCut ||
00356 TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) {
00357 objectMatch[i][j] = -1;
00358 }
00359 }
00360 delete [] stats;
00361 }
00362 objectMatch[j][i] = objectMatch[i][j];
00363 }
00364 else if(track2) {
00365 if(track->BelongsWithTrack(track2,ac,vldcptr,pTrkTrkDtpos2,
00366 pTrkTrkDz,pTrkTrkDt)) objectMatch[i][j] = 1;
00367 objectMatch[j][i] = objectMatch[i][j];
00368 }
00369 }
00370 MSG("EventSS",Msg::kDebug) << "Match for objects " << i << " and "
00371 << j << " : " << objectMatch[i][j] << endl;
00372 }
00373 }
00374
00375 MSG("EventSS",Msg::kDebug) << "Starting to form events" << endl;
00376
00377 //start to make events:
00378 Int_t nRecoObjects = recolist.GetLast()+1;
00379 //get total charge from each object
00380 std::vector<Float_t> ph(nRecoObjects,0);
00381 for(int i=0;i<nRecoObjects;i++){
00382 if (recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {
00383 CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00384 ph[i] = shower->GetCharge(CalStripType::kMIP);
00385 }
00386 if (recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00387 CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(i));
00388 //ph[i] = track->GetCharge(CalStripType::kMIP)+1e7; //force tracks to come first
00389 //force tracks to come first, order tracks by number of strips
00390 ph[i] = track->GetNStrip()+1e7;
00391 }
00392 }
00393
00394 Int_t nEvents = 0;
00395 Int_t nUnusedObjects = nRecoObjects;
00396 while(nUnusedObjects>0) { //loop until all reco objects are in events
00397
00398 //find largest object:
00399 Float_t maxPH = 0;
00400 Int_t maxPHIndex = -1;
00401 for(int i=0;i<nRecoObjects;i++){
00402 if(ph[i]>maxPH) {
00403 maxPH = ph[i];
00404 maxPHIndex = i;
00405 }
00406 }
00407
00408 //add largest objects first:
00409 TObjArray *event = new TObjArray;
00410 event->Add(recolist.At(maxPHIndex));
00411 ph[maxPHIndex] = 0;
00412 objectUsed[maxPHIndex] = true;
00413 nUnusedObjects-=1;
00414
00415 MSG("EventSS",Msg::kDebug) << "Added object " << maxPHIndex
00416 << " to event " << nEvents << endl;
00417
00418 if(recolist.At(maxPHIndex)->InheritsFrom("CandTrackHandle")) {
00419 MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a track" << endl;
00420
00421 //maxPH object is a track
00422 CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(maxPHIndex));
00423
00424 //now add objects clearly associated with largest object:
00425 for(int i=0;i<nRecoObjects;i++){
00426 if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) {
00427 if(recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00428 MSG("EventSS",Msg::kDebug) << "Track match. Adding object " << i
00429 << " to event " << nEvents << endl;
00430 //if this match is a track, add it to event
00431 event->Add(recolist.At(i));
00432 ph[i] = 0;
00433 objectUsed[i] = true;
00434 nUnusedObjects-=1;
00435 } // end if for CandTrackHandle
00436 else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {
00437 MSG("EventSS",Msg::kDebug) << "Shower match. Testing whether to add object "
00438 << i << " to event " << nEvents << endl;
00439 CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00440 //if this match is a shower, check where shower is relative to track vertex
00441 //if(TMath::Abs(shower2->GetVtxZ()-track->GetVtxZ())<pShwTrkDz/2.){
00442 if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) -
00443 track->GetBegPlane(PlaneView::kU)) <
00444 pShwTrkDz/(PITCHINMETRES*2.) ||
00445 TMath::Abs(shower2->GetBegPlane(PlaneView::kV) -
00446 track->GetBegPlane(PlaneView::kV)) <
00447 pShwTrkDz/(PITCHINMETRES*2.)){
00448 MSG("EventSS",Msg::kDebug) << "Shower " << i << " is close to track vertex. "
00449 << "Adding to event " << nEvents << endl;
00450 event->Add(recolist.At(i));
00451 ph[i] = 0;
00452 objectUsed[i] = true;
00453 nUnusedObjects-=1;
00454 //look for objects associated with this vertex shower
00455 for(int j=0;j<nRecoObjects;j++){
00456 if(j!=maxPHIndex && objectMatch[i][j]!=0 && !objectUsed[j]) { //got a match
00457 if(objectMatch[maxPHIndex][j]==0) { // j does not match maxPHIndex
00458 if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00459 MSG("EventSS",Msg::kDebug) << "Shower " << i
00460 << " has another track match. "
00461 << "Test how buried track " << j
00462 << " is in shower" << endl;
00463 if(objectMatch[i][j]==-1){
00464 event->Add(recolist.At(j));
00465 ph[j] = 0;
00466 objectUsed[j] = true;
00467 nUnusedObjects-=1;
00468 }
00469 else {
00470 CandTrackHandle *track2 =
00471 dynamic_cast<CandTrackHandle*>(recolist.At(j));
00472 //check if track2 is buried
00473 Int_t *stats = new Int_t[8];
00474 shower2->BuriedTrack(track2,1000,stats);
00475 if(stats[3]>0){
00476 Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00477 MSG("EventSS",Msg::kDebug)
00478 << "BuriedTrack Stats: " << "\n"
00479 << "Number of track planes = " << nplanes << "\n"
00480 << "nSharedPlanes = " << stats[0] << "\n"
00481 << "nMissingIsoPlanes = " << stats[1] << "\n"
00482 << "nMissingSharedPlanes = " << stats[2] << "\n"
00483 << "nBuriedPlanes = " << stats[3] << "\n"
00484 << "nPhysBuriedPlanes = " << stats[4] << "\n"
00485 << "nTrkLike = " << stats[5] << "\n"
00486 << "nTrkLikeTopol = " << stats[6] << "\n"
00487 << "nXTalk = " << stats[7] << endl;
00488 nplanes -= (stats[1] + stats[2]);
00489 if( (Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac) ||
00490 (Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) >
00491 pHardBuriedFrac) ||
00492 (stats[3]>=Int_t(nplanes*pSoftBuriedFrac) &&
00493 (TMath::Abs(track2->GetDirCosU())>pDCosUVCut ||
00494 TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) {
00495 MSG("EventSS",Msg::kDebug) << "Adding track " << j
00496 << " to event " << nEvents << endl;
00497 event->Add(recolist.At(j));
00498 ph[j] = 0;
00499 objectUsed[j] = true;
00500 nUnusedObjects-=1;
00501 }
00502 }
00503 delete [] stats;
00504 }
00505 }
00506 else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00507 MSG("EventSS",Msg::kDebug) << "Shower " << i
00508 << " has another shower match. "
00509 << "Test whether shower " << j
00510 << " has a better association" << endl;
00511 Bool_t isAssocTrack = false; //is this shower assoc with a track
00512 for(int k=0;k<nRecoObjects;k++) {
00513 //check for match
00514 if(k!=maxPHIndex && k!=i && objectMatch[j][k]!=0){
00515 //check for track
00516 if(recolist.At(k)->InheritsFrom("CandTrackHandle"))
00517 //check it does not match with shower2
00518 if(objectMatch[i][k]==0) isAssocTrack = true;
00519 }
00520 }
00521 if(!isAssocTrack) { //no track associated
00522 MSG("EventSS",Msg::kDebug) << "Shower " << j
00523 << " has no track association. "
00524 << "Add to event " << nEvents << endl;
00525 event->Add(recolist.At(j));
00526 ph[j] = 0;
00527 objectUsed[j] = true;
00528 nUnusedObjects-=1;
00529 }
00530 }
00531 }
00532 else { //j matches maxPHIndex, so add it to event
00533 MSG("EventSS",Msg::kDebug) << "Shower " << j
00534 << " matches orginal track " << maxPHIndex
00535 << " Add to event " << nEvents << endl;
00536 event->Add(recolist.At(j));
00537 ph[j] = 0;
00538 objectUsed[j] = true;
00539 nUnusedObjects-=1;
00540 }
00541 }
00542 }
00543 } //end if for vertex shower
00544 else { //if not a vertex shower, check other associations
00545 MSG("EventSS",Msg::kDebug) << "Testing whether shower " << i
00546 << " has better associations. " << endl;
00547 Bool_t otherTrackVertexMatch = false;
00548 Bool_t otherTrackMatch = false;
00549 Bool_t otherShowerMatch = false;
00550 for(int j=0;j<nRecoObjects;j++){
00551 if(j!=maxPHIndex && objectMatch[i][j]!=0 &&
00552 !objectUsed[j] && objectMatch[maxPHIndex][j]==0) {
00553 if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00554 CandTrackHandle *track2 =
00555 dynamic_cast<CandTrackHandle*>(recolist.At(j));
00556 if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) -
00557 track->GetBegPlane(PlaneView::kU)) <
00558 pShwTrkDz/(PITCHINMETRES*2.) ||
00559 TMath::Abs(shower2->GetBegPlane(PlaneView::kV) -
00560 track->GetBegPlane(PlaneView::kV)) <
00561 pShwTrkDz/(PITCHINMETRES*2.) )
00562 otherTrackVertexMatch = true;
00563 else {
00564 //first test whether track2 is buried
00565 //if so then can still associate shower to track
00566 if(objectMatch[i][j]!=-1) {
00567 //object match is not flagged as buried
00568 //perform check:
00569 Bool_t track2Buried = false;
00570 Int_t *stats2 = new Int_t[8];
00571 shower2->BuriedTrack(track2,1000,stats2);
00572 if(stats2[3]>0){
00573 Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00574 MSG("EventSS",Msg::kDebug)
00575 << "BuriedTrack Stats: " << "\n"
00576 << "Number of track planes = " << nplanes << "\n"
00577 << "nSharedPlanes = " << stats2[0] << "\n"
00578 << "nMissingIsoPlanes = " << stats2[1] << "\n"
00579 << "nMissingSharedPlanes = " << stats2[2] << "\n"
00580 << "nBuriedPlanes = " << stats2[3] << "\n"
00581 << "nPhysBuriedPlanes = " << stats2[4] << "\n"
00582 << "nTrkLike = " << stats2[5] << "\n"
00583 << "nTrkLikeTopol = " << stats2[6] << "\n"
00584 << "nXTalk = " << stats2[7] << endl;
00585 nplanes -= (stats2[1] + stats2[2]);
00586 if( (Float_t(stats2[4])/Float_t(nplanes) > pHardBuriedFrac) ||
00587 (Float_t(stats2[3]-stats2[5]-stats2[6])/Float_t(nplanes) >
00588 pHardBuriedFrac) ||
00589 (stats2[3]>=Int_t(nplanes*pSoftBuriedFrac) &&
00590 (TMath::Abs(track2->GetDirCosU())>pDCosUVCut ||
00591 TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) {
00592 track2Buried = true;
00593 objectMatch[i][j]=-1;
00594 objectMatch[j][i]=-1;
00595 }
00596 }
00597 if(!track2Buried){
00598 //test whether track2 matches better than track
00599 Int_t *stats = new Int_t[8];
00600 shower2->BuriedTrack(track,1000,stats);
00601 if( ( stats2[3] - stats2[5] - stats2[6] ) >=
00602 ( stats[3] - stats[5] - stats[6] ) &&
00603 ( stats2[4] >= stats[4] ) ) {
00604 otherTrackMatch = true;
00605 }
00606 delete [] stats;
00607 }
00608 delete [] stats2;
00609 }
00610 }
00611 }
00612 else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00613 Int_t *stats = new Int_t[8];
00614 shower2->BuriedTrack(track,1000,stats);
00615 if(stats[4]==0) { //how convincing is first track/shower match?
00616 otherShowerMatch = true;
00617 }
00618 delete [] stats;
00619 }
00620 }
00621 }
00622
00623 //if there are no other viable matches:
00624 if(otherTrackVertexMatch || otherShowerMatch ||
00625 otherTrackMatch) {
00626 MSG("EventSS",Msg::kDebug) << "Shower " << i << " has better association: "
00627 << " otherTrackVertexMatch = "
00628 << otherTrackVertexMatch
00629 << " otherShowerMatch = "
00630 << otherShowerMatch
00631 << " otherTrackMatch = "
00632 << otherTrackMatch << endl;
00633 continue;
00634 }
00635 //then add this shower to this track
00636 MSG("EventSS",Msg::kDebug) << "Shower " << i << " has no better association. "
00637 << "Add to event " << nEvents << endl;
00638 event->Add(recolist.At(i));
00639 ph[i] = 0;
00640 objectUsed[i] = true;
00641 nUnusedObjects-=1;
00642 } //end if for non-vertex shower
00643 } //end if for CandShowerSRHandle
00644 } //end if for valid association
00645 } //end for loop over all objects
00646 } //end if maxPH object is a track
00647
00648 else if(recolist.At(maxPHIndex)->InheritsFrom("CandShowerSRHandle")) {
00649 MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a shower" << endl;
00650 MSG("EventSS",Msg::kDebug) << "Check whether shower " << maxPHIndex
00651 << " has other associations" << endl;
00652 //maxPH object is a shower
00653 //CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(maxPHIndex));
00654 //now add objects clearly associated with largest object:
00655 for(int i=0;i<nRecoObjects;i++){
00656 if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) {
00657 if(recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00658 MSG("EventSS",Msg::kWarning) << "Arrgh - Logic flaw in eventSR. "
00659 << "Seeing tracks in event building, "
00660 << "when all tracks should already be "
00661 << "assigned!" << endl;
00662 } //end if CandTrackHandle
00663 else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {
00664 //check for two more levels of associated showers
00665 //CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00666 MSG("EventSS",Msg::kDebug) << "Adding shower " << i << " to event "
00667 << nEvents << endl;
00668 MSG("EventSS",Msg::kDebug) << "Check whether shower " << i
00669 << " has other associations" << endl;
00670 event->Add(recolist.At(i));
00671 ph[i] = 0;
00672 objectUsed[i] = true;
00673 nUnusedObjects-=1;
00674 for(int j=0;j<nRecoObjects;j++){
00675 if(j!=maxPHIndex && objectMatch[i][j]!=0 &&
00676 !objectUsed[j] && objectMatch[maxPHIndex][j]==0) {
00677 if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00678 MSG("EventSS",Msg::kWarning) << "Arrgh again - Logic flaw in eventSR. "
00679 << "Seeing tracks in event building, "
00680 << "when all tracks should already be "
00681 << "assigned!" << endl;
00682 }
00683 else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00684 MSG("EventSS",Msg::kDebug) << "Adding shower " << j << " to event "
00685 << nEvents << endl;
00686 event->Add(recolist.At(j));
00687 ph[j] = 0;
00688 objectUsed[j] = true;
00689 nUnusedObjects-=1;
00690 }
00691 }
00692 }
00693 } //end if CandShowerSRHandle
00694 } //end if valid association
00695 } //end loop over all objects
00696 } //end if maxPH object is a shower
00697 //add event to eventlist
00698 eventlist.Add(event);
00699 nEvents+=1;
00700 }
00701
00702 // remove all empty events from event list
00703
00704 for (Int_t i=0; i<=eventlist.GetLast(); i++) {
00705 TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i));
00706 if (objectlist->GetLast()<0) eventlist.Remove(objectlist);
00707 else if(objectlist->GetLast()==0) { //consider single shower/track events
00708 if(objectlist->At(0)->InheritsFrom("CandShowerSRHandle")) {
00709 MSG("EventSS",Msg::kDebug) << "Have a single candshower event " << endl;
00710 CandShowerSRHandle *showerSR =
00711 dynamic_cast<CandShowerSRHandle*> (objectlist->At(0));
00712 //if shower seems unphysical, remove it
00713 if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00714 eventlist.Remove(objectlist);
00715 MSG("EventSS",Msg::kDebug) << "REMOVING SHOWER-ONLY EVENT" << endl;
00716 continue;
00717 }
00718 }
00719 else if(objectlist->At(0)->InheritsFrom("CandTrackHandle")) {
00720 CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(0));
00721 if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) {
00722 //if track seems unphysical, remove it
00723
00724 eventlist.Remove(objectlist);
00725 MSG("EventSS",Msg::kDebug) << "REMOVING TRACK-ONLY EVENT" << endl;
00726 continue;
00727 }
00728 }
00729 }
00730 else if(objectlist->GetLast()==1) { //multiple shower/track events
00731 for(int ind_loop = 0;ind_loop<2;ind_loop++){
00732 Int_t ind_loop1 = ind_loop;
00733 Int_t ind_loop2 = ind_loop+1; if(ind_loop2==2) ind_loop2 = 0;
00734 if(objectlist->At(ind_loop1)->InheritsFrom("CandShowerSRHandle")){
00735 CandShowerSRHandle *showerSR =
00736 dynamic_cast<CandShowerSRHandle*> (objectlist->At(ind_loop1));
00737 if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00738 if(objectlist->At(ind_loop2)->InheritsFrom("CandTrackHandle")){
00739 CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(ind_loop2));
00740 if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) {
00741 eventlist.Remove(objectlist);
00742 MSG("EventSS",Msg::kDebug) << "REMOVING UNPHYSICAL SHOWER+TRACK EVENT" << endl;
00743 continue;
00744 }
00745 }
00746 }
00747 }
00748 }
00749 }
00750 if(useChopper==1) {
00751 std::vector<CandHandle> evtVec;
00752 for(int i=0;i<=objectlist->GetLast();i++){
00753 CandHandle *event = dynamic_cast<CandHandle*>(objectlist->At(i));
00754 evtVec.push_back(*event);
00755 }
00756 ChopHelper chopper(cx.GetMom());
00757 ChopHelp *chop = chopper.GetChopHelp(evtVec);
00758 if(chop->nchop==1) continue;
00759 else {
00760 chop->Print();
00761 //do something with chop info
00762 //for(int i=0;i<ncandidates;i++){
00763 //}
00764 }
00765 }
00766 }
00767 eventlist.Compress();
00768
00769
00770 // now make CandEvents for each event we have found
00771 Bool_t buildEvent=false;
00772 for (Int_t i=0; i<=eventlist.GetLast(); i++) {
00773 TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i));
00774 MSG("EventSS",Msg::kDebug) << "making event of "
00775 << objectlist->GetLast()+1
00776 << " objects \n";
00777 cxx.SetDataIn(objectlist);
00778 buildEvent=true;
00779 CandEventHandle eventhandle = CandEvent::MakeCandidate(ah,cxx);
00780 MSG("EventSS",Msg::kDebug) << " # of strips: "
00781 << eventhandle.GetNDaughters() << "\n";
00782 eventhandle.SetCandSlice(slice);
00783 ch.AddDaughterLink(eventhandle);
00784 delete eventlist.At(i);
00785 }
00786 if (!buildEvent) BuildEventFromUnassoc(ac,ch,cx,slice);
00787 } // end loop over slice list
00788
00789 delete trackItr;
00790 delete showerItr;
00791 delete subshowerItr;
00792
00793 MSG("EventSS",Msg::kDebug) << "starting unassociated hits assignment \n";
00794
00795 // find unassociated hits, and place them in tobjarray
00796 sliceItr.Reset();
00797 TObjArray unassociated;
00798 CandEventListHandle &eventlist = dynamic_cast<CandEventListHandle &>(ch);
00799 CandEventHandleItr * eventItr2 = 0;
00800 if (eventlist.GetNDaughters()>0) {
00801 eventItr2 = new CandEventHandleItr(eventlist.GetDaughterIterator());
00802 }
00803 FindUnassociated(sliceItr,eventItr2,unassociated);
00804 delete eventItr2;
00805 SetPrimaryShowers(ch,ac);
00806
00807 // now loop over unnassociated hits, and check for adjacency with all
00808 // showers, and track vertices (in the case that no shower exists).
00809 // If a match is found in the former case, add strip to shower.
00810 // In the latter case, create a new shower.
00811
00812 // this map contains the separation between each event/unassociated hit
00813 vector<map<const CandEventHandle*, Double_t> > dist2map(unassociated.GetLast()+1);
00814 FillDist2Map(ac,unassociated,ch,dist2map);
00815
00816 // initialize a vector which holds a Bool_t
00817 // indicating whether a hit has already been used
00818 vector<Bool_t> alreadyUsed(unassociated.GetLast()+1);
00819 for (Int_t i=0; i<=unassociated.GetLast(); i++) alreadyUsed[i] = false;
00820
00821 // loop over unnassociated hits, adding them to events.
00822 // After each loop the affected distance map entries are recalculated,
00823 // and the loop is repeated. This continues until a loop adds no strips.
00824 Bool_t found=true;
00825 Int_t n_found=0;
00826 while (found) {
00827
00828 //found indicates whether a strip was added to an event in this loop
00829 found=false;
00830
00831 CandStripHandle *closeststrip= 0;
00832 CandEventHandle *closestevent= 0;
00833 Double_t closestdist2 = 9999999.;
00834 Int_t closesti=-1;
00835
00836 // loop over unassoc strips
00837 for (Int_t i=0; i<=unassociated.GetLast(); i++) {
00838 if(alreadyUsed[i]) continue;
00839
00840 CandStripHandle *strip =
00841 dynamic_cast<CandStripHandle*>(unassociated.At(i));
00842
00843 CandEventHandle *bestevent = 0;
00844 Double_t bestdist2 = 9999999.;
00845 // loop over events, and find event which is closest to strip
00846 // and in time
00847 TIter eventItr(ch.GetDaughterIterator());
00848 while (CandEventHandle *event =
00849 dynamic_cast<CandEventHandle*>(eventItr())) {
00850 Double_t dtime = strip->GetBegTime()-event->GetVtxT();
00851 if ( dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1 &&
00852 (!bestevent || dist2map[i][event]<bestdist2) ) {
00853 bestevent = event;
00854 bestdist2 = dist2map[i][event];
00855 }
00856 }
00857 if (bestevent && (!closeststrip || bestdist2<closestdist2)) {
00858 closeststrip = strip;
00859 closesti=i;
00860 closestdist2 = bestdist2;
00861 closestevent = bestevent;
00862 }
00863 }
00864 MSG("EventSS",Msg::kVerbose)<< "closest strip " << closestdist2
00865 << "/"<<pHitAssocMaxDist2<< endl;
00866
00867 // if distance to closest event less than maximum allowable add to event
00868 if(closeststrip && closestevent && closestdist2<pHitAssocMaxDist2){
00869 MSG("EventSS",Msg::kVerbose) << "closest strip "
00870 << closeststrip->GetPlane() << " "
00871 << closeststrip->GetTPos()
00872 << " " << closestdist2 << "\n";
00873 alreadyUsed[closesti] = true;
00874 found = true;
00875 ++n_found;
00876
00877 // now either add this strip to existing primary shower, or
00878 // create a new shower at the track vertex if one does not exist
00879 // Get PrimaryShower CandHandle from fShowerList for modification.
00880 CandShowerHandle *shower = 0;
00881 if (CandShowerHandle *primsh = closestevent->GetPrimaryShower()) {
00882 for (Int_t ishower=0;
00883 !shower && ishower<closestevent->GetLastShower()+1;
00884 ishower++) {
00885 const CandShowerHandle *target = closestevent->GetShower(ishower);
00886 if (primsh->IsCloneOf(*target)) // Tests clone or Cand ==
00887 shower = closestevent->GetShowerWritable(ishower);
00888 }
00889 }
00890
00891 Float_t minShwPlane=0;
00892 Float_t maxShwPlane=0;
00893 if(shower){
00894 minShwPlane=min(shower->GetVtxPlane(),shower->GetEndPlane());
00895 maxShwPlane=max(shower->GetVtxPlane(),shower->GetEndPlane());
00896 }
00897 Float_t vertexsep = 0;
00898 CandTrackHandle *track = closestevent->GetPrimaryTrack();
00899 if (shower && track) {
00900 vertexsep = abs(shower->GetBegPlane() - track->GetBegPlane());
00901 vertexsep *= PITCHINMETRES;
00902 MSG("EventSS",Msg::kVerbose) << "Found track and primary shower. \n"
00903 << "Vertex separation (limit) = "
00904 << vertexsep << "("
00905 << pShwTrkDz << ")" << " \n";
00906 }
00907 // if event has shower and either no track, or shower is close
00908 // to track vertex, add this hit to the shower which was
00909 // compared in position to this strip earler.
00910 if (shower && (!track || vertexsep<pShwTrkDz)) {
00911 // NOTE: CandShowerHandle *shower, used to modify CandShower, is owned
00912 // by fShowerList.
00913 // add new strip to this shower, and the appropriate cluster
00914 MSG("EventSS",Msg::kVerbose) << "Shower min/max plane (limit): "
00915 << minShwPlane << "/" << maxShwPlane
00916 << " (" << pMaxNewShwLen << ")" << "\n";
00917 if(closeststrip->GetPlane()>minShwPlane-pMaxNewShwLen &&
00918 closeststrip->GetPlane()<maxShwPlane+pMaxNewShwLen ){
00919 MSG("EventSS",Msg::kVerbose) << "adding strip to shower \n";
00920 AddStripToEvent(closestevent,shower,subshowerlist,closeststrip);
00921 }
00922 }
00923 // no shower consistent with track vertex, and this strip
00924 // is close to vertex, so start new shower.
00925 else if (track) CreatePrimaryShower(ac,cxx,closestevent,showerlist,
00926 subshowerlist,closeststrip);
00927
00928 // need to recalculate distances bewteen remaining
00929 // unassoc hits and this modified event, if the added strip
00930 // results in a smaller separation, replace value in dist2map
00931 ReFillDist2Map(ac,closeststrip,closestevent,unassociated,dist2map);
00932 } // if adding strip to event in form of new shower
00933 } // if found (if an unassoc hit was matched to event in last loop over them
00934
00935 // because we have added strips to previously constructed showers, their
00936 // runalg method should be run on them again. In the future, an
00937 // optimization would be to do this only on showers which have been
00938 // touched.
00939 MSG("EventSS",Msg::kVerbose) << "added " << n_found << " of "
00940 << unassociated.GetLast()+1
00941 << " unassociated hits " << endl;
00942
00943 cxx.SetDataIn(showerlist); // Pass along persistable CandShowerList
00944 ReConstructShowers(ac,ch,cxx);
00945
00946 //Final step: use the chopper to combine events across slices:
00947 if(useChopper==1) {
00948 std::vector<CandHandle> evtVec;
00949 CandEventHandleItr evitr(ch.GetDaughterIterator());
00950 while(CandEventHandle *event = dynamic_cast<CandEventHandle*>(evitr())){
00951 evtVec.push_back(*event);
00952 }
00953 ChopHelper chopper(cx.GetMom());
00954 ChopHelp *chop = chopper.GetChopHelp(evtVec);
00955 chop->Print();
00956 }
00957
00958 // For each event, set primary shower based on
00959 // distance from vertex and energy
00960 SetPrimaryShowers(ch,ac);
00961 }
|
|
||||||||||||
|
Redetermine primary shower for each event. If the largest shower is greater than shwshwdz from the track vertex, and the closest shower has energy greater than minshwEfract call the closest shower the primary, otherwise the largest shower is the primary Definition at line 1370 of file AlgEventSSList.cxx. References CandHandle::GetDaughterIterator(), Registry::GetDouble(), and CandEventHandle::SetPrimaryShower(). Referenced by RunAlg(). 01371 {
01372
01373 // This method loops over events in the eventlist, finding the largest
01374 // shower and setting that to the primary shower in the event.
01375
01376 Double_t pShwShwDz = ac.GetDouble("ShwShwDz");
01377 Double_t pMinShwEFract = ac.GetDouble("MinShwEFract");
01378
01379 TIter evItr(ch.GetDaughterIterator());
01380 while (CandEventHandle *ev =
01381 dynamic_cast<CandEventHandle*>(evItr())) {
01382 ev->SetPrimaryShower(pMinShwEFract,pShwShwDz);
01383 }
01384 }
|
|
|
Reimplemented from AlgBase. Definition at line 2012 of file AlgEventSSList.cxx. 02013 {
02014 }
|
1.3.9.1