AlgEventSSList Class Reference

#include <AlgEventSSList.h>

Inheritance diagram for AlgEventSSList:

List of all members.

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 (AlgConfig &ac, 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)

---

Constructor & Destructor Documentation

AlgEventSSList::AlgEventSSList (   )

Definition at line 62 of file AlgEventSSList.cxx. { }

AlgEventSSList::~AlgEventSSList (   )  [virtual]

Definition at line 67 of file AlgEventSSList.cxx. { }

---

Member Function Documentation

void AlgEventSSList::AddObjectToEvent (  CandRecoHandle *  reco, TObjArray *  , TObjArray *  , Bool_t  )

Definition at line 1836 of file AlgEventSSList.cxx. References MSG. { if (!merge) { // if this is the first match we have found for this object, //add object to list // and point prevlist to this object list objectlist->Add(reco); MSG("EventSS",Msg::kVerbose) << " adding to list\n"; } else { MSG("EventSS",Msg::kVerbose) << " combining lists\n"; // this object matches more than one event - merge their member // objects to list prevlist, and add the new object for (Int_t ireco2=0; ireco2<=objectlist->GetLast();ireco2++) prevlist->Add(objectlist->At(ireco2)); prevlist->Add(reco); objectlist->Clear(); } }

void AlgEventSSList::AddStripToEvent (  CandEventHandle *  closestevent, CandShowerHandle *  shower, CandSubShowerSRListHandle *  subshowerlist, CandStripHandle *  closeststrip )

Add a strip to the event daughter list, the primary shower, and place the modified cluster on the cluster list Definition at line 976 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(). { MSG("EventSS",Msg::kVerbose) << "In AddStripToEvent:" << endl; closestevent->AddDaughterLink(*closeststrip); shower->AddDaughterLink(*closeststrip); if(subshowerlist && shower->InheritsFrom("CandShowerSRHandle")) { Calibrator& calibrator=Calibrator::Instance(); Int_t foundsubshower=0; Int_t closestSubShower=-1; Double_t BestStripDistance = 99999; Double_t StripDistance = 99999; CandShowerSRHandle *showerSR = dynamic_cast<CandShowerSRHandle*> (shower); MSG("EventSS",Msg::kVerbose) << "Shower UID = " << shower->GetUidInt() << endl; for (Int_t isubshower=0; isubshower<showerSR->GetLastSubShower()+1; isubshower++) { CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*> (showerSR->GetSubShower(isubshower)); //add strip to cluster of the same view if (ssl->GetPlaneView()==closeststrip->GetPlaneView()) { //decide which is most appropriate subshower to add strip to //based on distance to/spread of subshower: Double_t tposVtx = 0; if(ssl->GetPlaneView()==PlaneView::kU || ssl->GetPlaneView()==PlaneView::kX) tposVtx = ssl->GetVtxU(); else if(ssl->GetPlaneView()==PlaneView::kV || ssl->GetPlaneView()==PlaneView::kY) tposVtx = ssl->GetVtxV(); StripDistance = 99999; if(ssl->GetAvgDev()>0) { StripDistance = (TMath::Abs(closeststrip->GetTPos() - (tposVtx + ssl->GetSlope() * (closeststrip->GetZPos() - ssl->GetVtxZ()))) * TMath::Cos(TMath::ATan(ssl->GetSlope()))) * calibrator.GetMIP(closeststrip-> GetCharge(CalDigitType::kSigCorr)) / ssl->GetAvgDev(); } if(StripDistance>9999){ //bad slope StripDistance = 9999 + TMath::Sqrt(TMath::Power(TMath::Abs(closeststrip->GetZPos() - ssl->GetVtxZ()),2) + TMath::Power(TMath::Abs(closeststrip->GetTPos() - tposVtx),2)); } if(StripDistance<BestStripDistance || foundsubshower==0){ foundsubshower=1; BestStripDistance = StripDistance; closestSubShower = isubshower; } } } if(foundsubshower) { MSG("EventSS",Msg::kVerbose) << "adding strip to subshower " << closeststrip->GetPlane() << endl; CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*> (showerSR->GetSubShower(closestSubShower)); //extract mip -> gev conversion: Double_t GeVperMIP = 0; if(ssl->GetEnergy()!=0) GeVperMIP = ssl->GetEnergy()/ Double_t(calibrator.GetMIP(ssl->GetCharge(CalStripType::kSigCorr))); TIter csshItr(subshowerlist->GetDaughterIterator()); while (CandSubShowerSRHandle *cssh = dynamic_cast<CandSubShowerSRHandle*>(csshItr())) { if(cssh->IsEquivalent(ssl) && ssl->GetUidInt()==cssh->GetUidInt()) { MSG("EventSS",Msg::kVerbose) << "SubShower (old,new) UID " << ssl->GetUidInt() << " " << cssh->GetUidInt() << endl; showerSR->RemoveSubShower(ssl); CandSubShowerSRHandle *newss = cssh->DupHandle(); newss->SetCandSlice(cssh->GetCandSlice()); newss->AddDaughterLink(*closeststrip); newss->SetEnergy(GeVperMIP*Double_t(calibrator.GetMIP(closeststrip-> GetCharge(CalDigitType::kSigCorr))) + newss->GetEnergy()); subshowerlist->RemoveDaughter(cssh); subshowerlist->AddDaughterLink(*newss); showerSR->AddSubShower(newss); delete newss; break; } } } else { MSG("EventSS",Msg::kWarning) << "Strip added to shower but not to subshower!" << " StripDistance = " << StripDistance << endl; MSG("EventSS",Msg::kWarning) << "Number of subshowers in shower = " << showerSR->GetLastSubShower()+1 << endl; for(int i=0;i<showerSR->GetLastSubShower()+1;i++){ CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*> (showerSR->GetSubShower(i)); MSG("EventSS",Msg::kWarning) << "SS " << i << " UID = " << ssl->GetUidInt() << endl; } } } }

void AlgEventSSList::BuildEventFromUnassoc (  AlgConfig &  ac, CandHandle &  ch, CandContext &  cx, CandSliceHandle *  slice )

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 1399 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(). { assert(cx.GetDataIn()); MSG("EventSS",Msg::kDebug) <<" Building Event - no showers or tracks " << endl; if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) { MSG("EventSS",Msg::kDebug)<<" dataIn inherits from TObjArray " << endl; return; } Double_t pShwXTalkFracCut = ac.GetDouble("ShwXTalkFracCut"); Double_t pShwXTalkDef = ac.GetDouble("ShwXTalkDef"); const CandSliceListHandle *slicelist = 0; CandShowerListHandle *showerlist = 0; CandSubShowerSRListHandle *subshowerlist = 0; const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn()); for (Int_t i=0; i<=cxin->GetLast(); i++) { TObject *tobj = cxin->At(i); if (tobj->InheritsFrom("CandSliceListHandle")) { slicelist = dynamic_cast<CandSliceListHandle*>(tobj); } // Retrieve persistant modifiable CandShowerList if (tobj->InheritsFrom("CandShowerListHandle")) { showerlist = dynamic_cast<CandShowerListHandle*>(tobj); } if (tobj->InheritsFrom("CandSubShowerSRListHandle")) { subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj); } } if (!showerlist || !slicelist || !subshowerlist) { MSG("EventSS",Msg::kDebug) << "Missing either slicelist, showerlist or subshowerlist" << endl; return; } Int_t detectorType = slicelist->GetVldContext()->GetDetector(); CandContext cxx(this,cx.GetMom()); AlgFactory &af = AlgFactory::GetInstance(); const char *pEventAlgConfig = 0; ac.Get("EventAlgConfig",pEventAlgConfig); AlgHandle ah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig); Double_t pHitAssocDPlane = ac.GetInt("HitAssocDPlane"); Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0"); Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1"); Double_t pHitAssocMaxDist2 = ac.GetDouble("HitAssocMaxDist2"); Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm"); Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm"); Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm"); Double_t pShwShwDtpos2 = ac.GetDouble("ShwShwDtpos2"); Double_t pShwShwDz = ac.GetDouble("ShwShwDz"); Double_t pShwShwDt = ac.GetDouble("ShwShwDt"); Double_t pMinShwEFract = ac.GetDouble("MinShwEFract"); Double_t pMinUnassocStripPE = ac.GetDouble("MinUnassocStripPE"); const CandRecord *candrec = cx.GetCandRecord(); assert(candrec); const VldContext *vldcptr = candrec->GetVldContext(); assert(vldcptr); VldContext vldc = *vldcptr; UgliGeomHandle ugh(vldc); // Iterate over strips , finding sets which are 'in time', and within // a given distance of the nearest neighbor in the set, in both views. // Each of these sets is then used to construct a shower, which becomes // the primary shower of a new event. We continue until an iteration // in which no new event is created. // fill array of unique unassociated strips TObjArray unassociated; TIter stripItr(slice->GetDaughterIterator()); while (CandStripHandle *strip = dynamic_cast<CandStripHandle*>(stripItr())) { // Ignore strips below the required minimum PE if (strip->GetCharge(CalDigitType::kPE)<pMinUnassocStripPE) continue; Bool_t found=false; for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) { CandStripHandle *strip0 = dynamic_cast<CandStripHandle*>(unassociated.At(i)); if (*strip == *strip0) { found = true; } } if (!found) { if(!(detectorType==Detector::kNear && strip->GetPlane()>120)){ unassociated.Add(strip); } } } // loop until no new matches are found in the loop over unassoc hits. Bool_t found=false; while (unassociated.GetLast()>0 && !found) { MSG("EventSS",Msg::kDebug)<< " # unassoc " << unassociated.GetLast()+1 << endl; Bool_t firstu=true; Bool_t firstv=true; TObjArray neweventu; TObjArray neweventv; Bool_t addedstrip = true; Float_t totcharge=0; while (addedstrip) { addedstrip=false; // loop over unassoc strips. for (Int_t i=0; i<=unassociated.GetLast(); i++) { CandStripHandle *strip = dynamic_cast<CandStripHandle*>(unassociated.At(i)); switch (strip->GetPlaneView()) { case PlaneView::kU: // strip is in U plane // if this is the first strip, add to newevent strip list if (firstu) { // if we don't have v strips yet, just add this strip to u list // seed the event if (firstv) { neweventu.Add(strip); MSG("EventSS",Msg::kDebug)<< " add first " << endl; totcharge +=strip->GetCharge(); addedstrip = true; firstu=false; } // if we already have one or more v strips check that they are // in time, and agree in Z position. else { Bool_t found2=false; for (Int_t j=0; j<=neweventv.GetLast() && !found2; j++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventv.At(j)); Double_t dtime = strip->GetBegTime()-estrip->GetBegTime(); MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 << " " <<pHitAssocDTime0*1e9 << "/" << pHitAssocDTime1*1e9 << " plane diff " << abs(strip->GetPlane() - estrip->GetPlane()) << "/" << pHitAssocDPlane << endl; if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1 && abs(strip->GetPlane()-estrip->GetPlane()) <= pHitAssocDPlane) { MSG("EventSS",Msg::kVerbose)<< " add " << endl; // match is found, add this strip to array neweventu.Add(strip); totcharge +=strip->GetCharge(); addedstrip = true; firstu=0; found2=true; } // if match found } // for loop over V strips in new event } // end not first V strip } // end first U strip else { // if we already have u strips, make sure that the one we // add is adjaccent, using distance criteria based // on plane, time, and transverse postion separation Bool_t found2=false; for (Int_t j=0; j<=neweventu.GetLast() && !found2; j++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventu.At(j)); Double_t dtime = strip->GetBegTime()-estrip->GetBegTime(); MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 << " " <<pHitAssocDTime0*1e9 << "/" << pHitAssocDTime1*1e9 << " plane diff " << abs(strip->GetPlane() - estrip->GetPlane()) << " dtpos " << (strip->GetTPos() - estrip->GetTPos()) << endl; if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1) { Double_t dplane = (Double_t)(strip->GetPlane()-estrip->GetPlane()); Double_t dtpos = strip->GetTPos()-estrip->GetTPos(); Double_t dist2 = (pHitAssocZParm*dplane*dplane) + (pHitAssocPParm*dtpos*dtpos) + (pHitAssocTParm*dtime*dtime); // if dist2 2 parameter meets match criteria // add strip to array if (dist2<pHitAssocMaxDist2*4) { MSG("EventSS",Msg::kVerbose)<< " add " << endl; neweventu.Add(strip); totcharge +=strip->GetCharge(); addedstrip = true; firstu=0; found2=true; } // if positions match } // if in time } // end loop over event u strips } // end if already have u strip break; case PlaneView::kV: // if this is the first strip, add to newevent strip list if (firstv) { // strip is in V plane // if this is the first strip, add to newevent strip list if (firstu) { neweventv.Add(strip); totcharge +=strip->GetCharge(); MSG("EventSS",Msg::kDebug)<< " add first" << endl; addedstrip = true; firstv = 0; } else { // if we already have one or more u strips check that they are // in time, and require agreement in Z position. Bool_t found2=false; for (Int_t j=0; j<=neweventu.GetLast() && !found2; j++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventu.At(j)); Double_t dtime = strip->GetBegTime()-estrip->GetBegTime(); MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 << " " << pHitAssocDTime0*1e9 << "/" << pHitAssocDTime1*1e9 << " plane diff " << abs(strip->GetPlane() - estrip->GetPlane()) << "/" << pHitAssocDPlane << endl; if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1 && abs(strip->GetPlane()-estrip->GetPlane()) <= pHitAssocDPlane) { // if dist2 2 parameter meets match criteria, //add strip to array MSG("EventSS",Msg::kVerbose) << " add " << endl; neweventv.Add(strip); totcharge +=strip->GetCharge(); addedstrip = true; firstv = false; found2 = true; } } // end loop over u strips in event } // end if already have u strips } // end if this is the first v strip else { // if we already have v strips, make sure that the one we // add is adjaccent, using same distance criteria base // on plane, time, and transverse postion separation Bool_t found2=false; for (Int_t j=0; j<=neweventv.GetLast() && !found2; j++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventv.At(j)); Double_t dtime = strip->GetBegTime() - estrip->GetBegTime(); MSG("EventSS",Msg::kVerbose)<< " dtime " << dtime*1e9 << " " << pHitAssocDTime0*1e9 << "/" << pHitAssocDTime1*1e9 << " plane diff " << abs(strip->GetPlane() - estrip->GetPlane()) << " dtpos " << (strip->GetTPos() - estrip->GetTPos()) << endl; if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1) { Double_t dplane = (Double_t)(strip->GetPlane()-estrip->GetPlane()); Double_t dtpos = strip->GetTPos()-estrip->GetTPos(); Double_t dist2 = (pHitAssocZParm*dplane*dplane) + (pHitAssocPParm*dtpos*dtpos) + (pHitAssocTParm*dtime*dtime); // if dist2 2 parameter meets match criteria, // add strip to array if (dist2<pHitAssocMaxDist2*4) { MSG("EventSS",Msg::kVerbose) << " add " << endl; neweventv.Add(strip); totcharge +=strip->GetCharge(); addedstrip = true; firstv = false; found2 = true; } // end if add strip } //end if in time } // end loop over v strips } break; default: break; } } // end of loop over hits - remove new 'associations' and compress for (Int_t i=0; i<=neweventu.GetLast(); i++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventu.At(i)); unassociated.Remove(estrip); } for (Int_t i=0; i<=neweventv.GetLast(); i++) { CandStripHandle *estrip = dynamic_cast<CandStripHandle*>(neweventv.At(i)); unassociated.Remove(estrip); } unassociated.Compress(); } // if we have hits in both views, build an event from them // if pulse height>10 pes MSG("EventSS",Msg::kDebug)<< " considering event formation u/v: " << !firstu << "/" << !firstv << " charge " << totcharge << endl; if (!firstu && !firstv && totcharge>20.0) { found=true; //construct CandSubShowerSRs cxx.SetDataIn(&neweventu); AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default"); CandSubShowerSRHandle usubshowerhandle = CandSubShowerSR::MakeCandidate(subshowerah,cxx); usubshowerhandle.SetCandSlice(slice); usubshowerhandle.SetClusterID(ClusterType::kUnknown); subshowerlist->AddDaughterLink(usubshowerhandle); cxx.SetDataIn(&neweventv); CandSubShowerSRHandle vsubshowerhandle = CandSubShowerSR::MakeCandidate(subshowerah,cxx); vsubshowerhandle.SetCandSlice(slice); vsubshowerhandle.SetClusterID(ClusterType::kUnknown); subshowerlist->AddDaughterLink(vsubshowerhandle); CandHandle *usubshower=subshowerlist->FindDaughter(&usubshowerhandle); CandHandle *vsubshower=subshowerlist->FindDaughter(&vsubshowerhandle); //build shower TObjArray newshower; newshower.Add(usubshower); newshower.Add(vsubshower); cxx.SetDataIn(&newshower); AlgHandle showerah = af.GetAlgHandle("AlgShowerSS","default"); CandShowerSRHandle showerhandle = CandShowerSR::MakeCandidate(showerah,cxx); showerhandle.SetCandSlice(slice); showerhandle.SetCandRecord(slicelist->GetCandRecord()); Bool_t newShowerOverlapsOld=false; CandShowerHandleItr showerItr(showerlist->GetDaughterIterator()); while (CandShowerHandle *shower = dynamic_cast<CandShowerHandle*>(showerItr())){ Double_t dz = showerhandle.GetVtxZ()-shower->GetVtxZ(); Double_t du = showerhandle.GetVtxU()-shower->GetVtxU(); Double_t dv = showerhandle.GetVtxV()-shower->GetVtxV(); Double_t dt = showerhandle.GetVtxT()-shower->GetVtxT(); if( ( (du*du+dv*dv<pShwShwDtpos2 &&fabs(dz)<pShwShwDz) || (du*du<pShwShwDtpos2/4. && dv*dv<pShwShwDtpos2*4. && fabs(dz)<pShwShwDz/2.) || (du*du<pShwShwDtpos2*4 && dv*dv<pShwShwDtpos2/4. && fabs(dz)<pShwShwDz/2.) ) && fabs(dt)<pShwShwDt ) { newShowerOverlapsOld = true; MSG("EventSS",Msg::kDebug) << " new shower overlaps with previous - don't make new event " << endl; break; } } if(!newShowerOverlapsOld){ if(!showerhandle.IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) { MSG("EventSS",Msg::kDebug)<< " Building Event! " << endl; showerlist->AddDaughterLink(showerhandle); TObjArray recolist; recolist.Add(&showerhandle); cxx.SetDataIn(&recolist); AlgHandle eventah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig); CandEventHandle eventhandle = CandEvent::MakeCandidate(eventah,cxx); eventhandle.SetCandSlice(slice); eventhandle.SetPrimaryShower(pMinShwEFract,pShwShwDz); ch.AddDaughterLink(eventhandle); } } } } }

void AlgEventSSList::CreatePrimaryShower (  AlgConfig &  ac, CandContext &  cxx, CandEventHandle *  closestevent, CandShowerListHandle *  showerlist, CandSubShowerSRListHandle *  subshowerlist, CandStripHandle *  closeststrip )

Definition at line 1087 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(). { // need to create a shower -- include strips from first plane in each // view plus this previously unassociated hit to generate a 3D object MSG("EventSS",Msg::kVerbose) << "In CreatePrimaryShower" << endl; // Get singleton instance of AlgFactory const char *pEventAlgConfig = 0; ac.Get("EventAlgConfig",pEventAlgConfig); //Double_t shwshwdz = ac.GetDouble("ShwShwDz"); //Double_t minShwEFract = ac.GetDouble("MinShwEFract"); Double_t pMinShwStripPE = ac.GetDouble("MinShwStripPE"); Double_t pMaxNewShwLen = ac.GetDouble("MaxNewShwLen"); AlgFactory &af = AlgFactory::GetInstance(); CandTrackHandle *track = closestevent->GetPrimaryTrack(); if (!showerlist) { MSG("EventSS",Msg::kError) <<"CandShowerListHandle missing\n" << endl; return; } closestevent->AddDaughterLink(*closeststrip); // create arrays of strips for each planeview to create // shower from MSG("EventSS",Msg::kVerbose) << "creating shower \n"; TObjArray ustriparray; TObjArray vstriparray; Int_t uplane=track->GetBegPlane(PlaneView::kU); Int_t vplane=track->GetBegPlane(PlaneView::kV); // seed this shower with strips in most upstream track hits // in each view, along with hits above a charge threshold and // within pMaxNewShwLen of the track start TIter stripItr(track->GetDaughterIterator()); while (CandStripHandle *tstrip = dynamic_cast<CandStripHandle*>(stripItr())) { if ((tstrip->GetPlane() == uplane) || (tstrip->GetPlaneView()==PlaneView::kU && tstrip->GetCharge()>pMinShwStripPE && abs(tstrip->GetPlane()-uplane)<pMaxNewShwLen)) { ustriparray.Add(tstrip); } if ((tstrip->GetPlane() == vplane) || (tstrip->GetPlaneView()==PlaneView::kV && tstrip->GetCharge()>pMinShwStripPE && abs(tstrip->GetPlane()-vplane)<pMaxNewShwLen)) { vstriparray.Add(tstrip); } } // add unassoc strip to appropriate array switch (closeststrip->GetPlaneView()) { case PlaneView::kU: ustriparray.Add(closeststrip); break; case PlaneView::kV: vstriparray.Add(closeststrip); break; default: break; } // build shower from U/V clusters or subshowers if(ustriparray.GetLast()>=0 && vstriparray.GetLast()>=0){ if(showerlist->InheritsFrom("CandShowerSRListHandle") && subshowerlist){ MSG("EventSS",Msg::kVerbose) << "U striparray has " << ustriparray.GetLast()+1 << " entries" << endl; MSG("EventSS",Msg::kVerbose) << "V striparray has " << vstriparray.GetLast()+1 << " entries" << endl; // now construct CandSubShowers AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default"); cxx.SetDataIn(&ustriparray); CandSubShowerSRHandle usubshowerhandle = CandSubShowerSR::MakeCandidate(subshowerah,cxx); usubshowerhandle.SetCandSlice(closestevent->GetCandSlice()); usubshowerhandle.SetClusterID(ClusterType::kUnknown); subshowerlist->AddDaughterLink(usubshowerhandle); MSG("EventSS",Msg::kVerbose) << "U subshower view: " << usubshowerhandle.GetPlaneView() <<endl; cxx.SetDataIn(&vstriparray); CandSubShowerSRHandle vsubshowerhandle = CandSubShowerSR::MakeCandidate(subshowerah,cxx); vsubshowerhandle.SetCandSlice(closestevent->GetCandSlice()); vsubshowerhandle.SetClusterID(ClusterType::kUnknown); subshowerlist->AddDaughterLink(vsubshowerhandle); MSG("EventSS",Msg::kVerbose) << "V subshower view: " << vsubshowerhandle.GetPlaneView() <<endl; //add subshowers to list CandShowerSRListHandle* showerlistSR = dynamic_cast<CandShowerSRListHandle*>(showerlist); TObjArray newshower; CandSubShowerSRHandle * usubshower = dynamic_cast<CandSubShowerSRHandle*> (subshowerlist->FindDaughter(&usubshowerhandle)); CandSubShowerSRHandle * vsubshower = dynamic_cast<CandSubShowerSRHandle*> (subshowerlist->FindDaughter(&vsubshowerhandle)); newshower.Add(usubshower); newshower.Add(vsubshower); cxx.SetDataIn(&newshower); //get ShowerSS alghandle and make shower from subshowers AlgHandle showerahSS = af.GetAlgHandle("AlgShowerSS","default"); CandShowerSRHandle showerhandleSR = CandShowerSR::MakeCandidate(showerahSS,cxx); showerhandleSR.SetCandSlice(closestevent->GetCandSlice()); showerlistSR->AddDaughterLink(showerhandleSR); CandHandle *shw = showerlistSR->FindDaughter(&showerhandleSR); CandShowerHandle *shwh = dynamic_cast<CandShowerHandle*>(shw); closestevent->AddShower(shwh); // Change: cbs 9Aug06 // NEED TO FORCE CREATED SHOWER TO BE PRIMARY // due to a change in assignment of primary shower, can now create // a shower here which is *not* called primary by SetPrimaryShower() // Causes multiple CreatePrimaryShower calls in event builder and // large numbers of reco'd showers per events. // SetPrimaryShowers is called at the end of RunAlg so a correct // assignment will be made later //closestevent->SetPrimaryShower(minShwEFract,shwshwdz); closestevent->SetPrimaryShower(shwh); } } }

void AlgEventSSList::FillDist2Map (  AlgConfig &  ac, TObjArray &  unassociated, CandHandle &  ch, std::vector< std::map< const CandEventHandle *, Double_t > > &  dist2map )

Definition at line 1911 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(). { // grab all needed algorithm parameters Double_t pShwTrkDz = ac.GetDouble("ShwTrkDz"); Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0"); Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1"); Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm"); Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm"); Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm"); // construct the map, first looping over strips for (Int_t i=0; i<=unassociated.GetLast(); i++) { CandStripHandle *strip = dynamic_cast<CandStripHandle*>(unassociated.At(i)); TIter eventItr(ch.GetDaughterIterator()); // inner loop is over events while (CandEventHandle *event = dynamic_cast<CandEventHandle*>(eventItr())) { // calculate separation in space/time Double_t dtime = (strip->GetBegTime()-event->GetVtxT()); Double_t bestdist2=-1.; Double_t bestdplane=-1; Double_t bestdtpos=-1; Double_t bestdtime=-1; dist2map[i][event] = 999999999.; MSG("EventSS",Msg::kVerbose) << i << " " << strip->GetPlane() << " " << strip->GetStrip() << " stp beg time " << strip->GetBegTime() << " evt vtx time " << event->GetVtxT() << " dtime " << dtime*1e9 << "/" << pHitAssocDTime0*1e9 << "/" << pHitAssocDTime1*1e9 << endl; // if there is a rough time match proceed if (dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1) { Float_t vertexsep = 0; // If event has both shower // and track calc separation between the two. If this // separation is greater // than pSHwTrkDz, we won't bother adding the strip to the // primary shower, and will instead consider whether // to create a new shower at track vertex. if (event->GetPrimaryShower() && event->GetPrimaryTrack()) vertexsep = fabs(event->GetPrimaryShower()->GetVtxZ() - event->GetPrimaryTrack()->GetVtxZ()); // the event has a shower, and either no track, or track and // shower vertices are in agreement. In this case, we // consider adding this strip to this shower, so calc separation. if (event->GetPrimaryShower() && (!event->GetPrimaryTrack() || vertexsep<pShwTrkDz)) { // create iterator of shower daughter strips in same view // as unassoc strip CandStripHandleItr shwstripItr(event->GetPrimaryShower()-> GetDaughterIterator()); CandStripHandleKeyFunc *shwstripKf = shwstripItr.CreateKeyFunc(); shwstripKf->SetFun(CandStripHandle::KeyFromView); shwstripItr.GetSet()->AdoptSortKeyFunc(shwstripKf); shwstripKf = 0; shwstripItr.GetSet()->Slice(); shwstripItr.GetSet()->Slice(strip->GetPlaneView(), strip->GetPlaneView()); // iterate over the daughter list, and find minimum distane // between daughter and unassoc strip while (CandStripHandle *shwstrip = dynamic_cast<CandStripHandle*>(shwstripItr())) { Double_t dplane = (Double_t)(strip->GetPlane()-shwstrip->GetPlane()); Double_t dtpos = strip->GetTPos()-shwstrip->GetTPos(); Double_t dist2 = (pHitAssocZParm*dplane*dplane) + (pHitAssocPParm*dtpos*dtpos) + (pHitAssocTParm*dtime*dtime); if (bestdist2<0. || dist2<bestdist2) { bestdist2 = dist2; bestdplane=dplane; bestdtpos=dtpos; bestdtime=dtime; } } MSG("EventSS",Msg::kVerbose) << "primary shower." << " dplane:" << bestdplane <<" dtpos:" << bestdtpos <<" dtime:" << bestdtime <<" dist2:" << bestdist2<< "\n"; dist2map[i][event] = bestdist2; } // end if event has shower, etc. // else if criteria above was not satisfied, consider whether // strip is near track vertex, in which case we start new shower else if (event->GetPrimaryTrack()) { // calc distance to vertex Double_t dplane = (Double_t)(strip->GetPlane() - event->GetPrimaryTrack()-> GetBegPlane()); Double_t dtpos=0; switch (strip->GetPlaneView()) { case PlaneView::kU: dtpos=strip->GetTPos()-event->GetPrimaryTrack()-> GetU(event->GetPrimaryTrack()->GetBegPlane()); break; case PlaneView::kV: dtpos=strip->GetTPos()-event->GetPrimaryTrack()-> GetV(event->GetPrimaryTrack()->GetBegPlane()); break; default: break; } Double_t dist2 = ( (pHitAssocZParm*dplane*dplane) + (pHitAssocPParm*dtpos*dtpos) + (pHitAssocTParm*dtime*dtime) ); dist2map[i][event] = dist2; MSG("EventSS",Msg::kVerbose) << "primary track. begin:" << " dplane:" << dplane << " dtpos:" << dtpos <<" dtime:" << dtime <<" dist2:" << dist2 << "\n"; } // end if primary track } // end if time match OK } // end loop over events } // end loop over unassoc strips }

void AlgEventSSList::FillRecoList (  CandSliceHandle *  , CandShowerHandleItr *  , CandTrackHandleItr *  , TObjArray &  )

Definition at line 1791 of file AlgEventSSList.cxx. References CandRecoHandle::GetCandSlice(), CandHandle::GetNDaughters(), CandHandle::IsCloneOf(), and MSG. Referenced by RunAlg(). { if (showerItr) { showerItr->Reset(); while (CandShowerHandle *shower = dynamic_cast<CandShowerHandle*>((*showerItr)())) { if (shower->GetCandSlice()) { if (slice->IsCloneOf(*(shower->GetCandSlice()))) { recolist.Add(shower); } } else { MSG("EventSS", Msg::kError) << "Shower doesn't contain a slice. Not added to recolist." << endl; } } } if (trackItr) { trackItr->Reset(); while (CandTrackHandle *track = dynamic_cast<CandTrackHandle*>((*trackItr)())) { if (track->GetCandSlice()) { if (slice->IsCloneOf(*track->GetCandSlice())) { if(track->GetNDaughters()==0 && track->InheritsFrom("CandFitTrackHandle") && dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack()){ track = dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack(); MSG("EventSS", Msg::kDebug)<< "Finder Track being used " << endl; } recolist.Add(track); } } else { MSG("EventSS", Msg::kError) << "Track doesn't contain a slice. Not added to recolist." << endl; } } } }

void AlgEventSSList::FindUnassociated (  AlgConfig &  ac, CandSliceHandleItr &  sliceItr, CandEventHandleItr *  , TObjArray &  unassociated )

Definition at line 1858 of file AlgEventSSList.cxx. References CandHandle::FindDaughter(), CandEventHandle::GetCandSlice(), CandStripHandle::GetCharge(), CandHandle::GetDaughterIterator(), Registry::GetDouble(), and CandHandle::IsCloneOf(). Referenced by RunAlg(). { Double_t pMinUnassocStripPE = ac.GetDouble("MinUnassocStripPE"); while (CandSliceHandle *slice = dynamic_cast<CandSliceHandle*>(sliceItr())) { TIter stripItr(slice->GetDaughterIterator()); while (CandStripHandle *strip = dynamic_cast<CandStripHandle*>(stripItr())) { // Ignore strips below the required minimum PE if (strip->GetCharge(CalDigitType::kPE)<pMinUnassocStripPE) continue; Bool_t found=false; // check for duplicates for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) { CandStripHandle *strip0 = dynamic_cast<CandStripHandle*>(unassociated.At(i)); if (*strip == *strip0) { found = true; break; } } // if this strip not already in unnassociated list, then // loop over objects in this slice, and check against all daughter strips if (!found) { if (eventItr) { eventItr->Reset(); while (CandEventHandle *event = dynamic_cast<CandEventHandle*>((*eventItr)())) { if (event->GetCandSlice()) { if (slice->IsCloneOf(*(event->GetCandSlice()))) { if (event->FindDaughter(strip)) { found = true; break; } } } } } } // this strip is not in any event, add to unassociated list if (!found) { unassociated.Add(strip); } } } }

void AlgEventSSList::ReConstructShowers (  AlgConfig &  ac, CandHandle &  ch, CandContext &  cxx )

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 1245 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(). { MSG("EventSS",Msg::kVerbose) << "In ReConstructShowers: " << endl; // Retrieve persistant modifiable CandShowerList CandShowerListHandle *showerlist = const_cast<CandShowerListHandle*> (dynamic_cast<const CandShowerListHandle*>(cxx.GetDataIn())); if(!showerlist->InheritsFrom("CandShowerSRListHandle")) { MSG("EventSS",Msg::kWarning) << "ListHandle is not a CandShowerSRListHandle." << " ShowerList not processed in ReConstructShowers." << endl; return; } // Get singleton instance of AlgFactory AlgFactory &af = AlgFactory::GetInstance(); // Get an AlgHandle to AlgShowerSR with default AlgConfig const char *pEventAlgConfig = 0; ac.Get("EventAlgConfig",pEventAlgConfig); AlgHandle ahss = af.GetAlgHandle("AlgShowerSS","default"); AlgConfig acss = ahss.GetAlgConfig(); Int_t eventCounter = 0; TIter evItr(ch.GetDaughterIterator()); while (CandEventHandle *ev = dynamic_cast<CandEventHandle*>(evItr())) { //event loop MSG("EventSS",Msg::kVerbose) << "event = " << eventCounter << endl; for (Int_t ishower=0; ishower<ev->GetLastShower()+1; ishower++) { //shower loop MSG("EventSS",Msg::kVerbose) << "shower = " << ishower << endl; CandShowerHandle *shower = dynamic_cast<CandShowerHandle*>(ev->GetShowerWritable(ishower)); CandShowerSRHandle *showerSR = 0; if(shower->InheritsFrom("CandShowerSRHandle")) { showerSR = dynamic_cast<CandShowerSRHandle*>(shower); } MSG("EventSS",Msg::kVerbose) << "Shower UID = " << showerSR->GetUidInt() << endl; //assume subshower chain if(showerSR) { TObjArray newshower; // fill array of subshowers for this shower for (Int_t isubshower=0; isubshower<showerSR->GetLastSubShower()+1; isubshower++) { MSG("EventSS",Msg::kVerbose) << "subshower = " << isubshower << endl; const CandSubShowerSRHandle *ssh = showerSR->GetSubShower(isubshower); if (ssh) { // New handle on heap CandSubShowerSRHandle *sshd = ssh->DupHandle(); // Add new heap-based CandHandle to TObjArray newshower.Add(sshd); MSG("EventSS",Msg::kDebug) << "Uid's (old,new) " << ssh->GetUidInt() << " " << sshd->GetUidInt() << endl; } } // re-run algorithm for this shower if (newshower.GetEntriesFast() > 0) { cxx.SetDataIn(&newshower); ahss.RunAlg(*showerSR,cxx); } else { MSG("EventSS", Msg::kError) << "Attempt to pass empty TObjArray newshower to AlgShowerSS" << " in CandContext. AlgShowerSS::RunAlg() call ignored." << endl; } // Delete heap-based CandHandles in TObjArray newshower.Delete(); } else { MSG("EventSS",Msg::kWarning) << "Handle is not a CandShowerSRHandle." << " Shower not processed in ReConstructShowers." << endl; continue; } MSG("EventSS",Msg::kVerbose) << "Modified shower UID = " << showerSR->GetUidInt() << endl; // now loop through strips in this shower, and determine whether primary // track in this event intercepts this strip. If so, subtract 1 mip // from shower energy. // if event has primary track, extrapolate through shower // and subtract energy loss to recalibrate CandTrackHandle *primaryTrack = ev->GetPrimaryTrack(); if (primaryTrack) { shower->CalibrateEnergy(primaryTrack,acss); } // Update modified showers in persistent modifiable CandShowerList. // CandHandle could potentially supply a ReplaceDaughter() function. // Interface of ReplaceDaughter() would be debatable. if (showerlist) { TIter shiter(showerlist->GetDaughterIterator()); CandShowerHandle *target; Bool_t found = kFALSE; while (!found && (target = dynamic_cast<CandShowerHandle*>(shiter()))) { if (shower->IsCloneOf(*target)) { // Tests clone or Cand == found = kTRUE; // Replace target only if shower is a modified version of target if (!(shower->IsEqual(target))) { // CandBase inequality shower->SetCandSlice(ev->GetCandSlice()); if (!showerlist->RemoveDaughter(target)) { // Failure MSG("EventSS",Msg::kError) << endl << "Failure of ShowerList daughter removal " << endl << "attempted during replacement by modified Shower." << endl << "Will result in double counted Shower." << endl; } showerlist->AddDaughterLink(*shower); } } } // Add new (unfound) shower to persistent modifiable CandShowerList. if (!found){ shower->SetCandSlice(ev->GetCandSlice()); showerlist->AddDaughterLink(*shower); } } } eventCounter++; } }

void AlgEventSSList::ReFillDist2Map (  AlgConfig &  ac, CandStripHandle *  closeststrip, CandEventHandle *  closestevent, TObjArray &  unassociated, std::vector< std::map< const CandEventHandle *, Double_t > > &  dist2map )

Definition at line 1217 of file AlgEventSSList.cxx. References CandStripHandle::GetBegTime(), Registry::GetDouble(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandStripHandle::GetTPos(), and CandRecoHandle::GetVtxT(). Referenced by RunAlg(). { Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm"); Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm"); Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm"); for (Int_t i=0; i<=unassociated.GetLast(); i++) { CandStripHandle *strip = dynamic_cast<CandStripHandle*>(unassociated.At(i)); if(strip->GetPlaneView()==closeststrip->GetPlaneView()){ Double_t dtime = strip->GetBegTime()-closestevent->GetVtxT(); Double_t dplane = (Double_t)(strip->GetPlane()-closeststrip->GetPlane()); Double_t dtpos = strip->GetTPos()-closeststrip->GetTPos(); Double_t dist2 = (pHitAssocZParm*dplane*dplane) + (pHitAssocPParm*dtpos*dtpos) + (pHitAssocTParm*dtime*dtime); if (dist2<dist2map[i][closestevent]) { dist2map[i][closestevent] = dist2; } } } // loop over unassoc hits }

void AlgEventSSList::RunAlg (  AlgConfig &  ac, CandHandle &  ch, CandContext &  cx )  [virtual]

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(). { assert(cx.GetDataIn()); if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) return; CandContext cxx(this,cx.GetMom()); AlgFactory &af = AlgFactory::GetInstance(); const char *pEventAlgConfig = 0; ac.Get("EventAlgConfig",pEventAlgConfig); AlgHandle ah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig); // grab all needed algorithm parameters //generic parameters: Double_t pTrkTrkDtpos2 = ac.GetDouble("TrkTrkDtpos2"); Double_t pTrkTrkDz = ac.GetDouble("TrkTrkDz"); Double_t pTrkTrkDt = ac.GetDouble("TrkTrkDt"); Double_t pShwTrkDtpos2 = ac.GetDouble("ShwTrkDtpos2"); Double_t pShwTrkDz = ac.GetDouble("ShwTrkDz"); Double_t pShwTrkDt = ac.GetDouble("ShwTrkDt"); Double_t pShwShwDtpos2 = ac.GetDouble("ShwShwDtpos2"); Double_t pShwShwDz = ac.GetDouble("ShwShwDz"); Double_t pShwShwDt = ac.GetDouble("ShwShwDt"); Double_t pMaxNewShwLen = ac.GetDouble("MaxNewShwLen"); Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0"); Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1"); Double_t pHitAssocMaxDist2 = ac.GetDouble("HitAssocMaxDist2"); Double_t pMinTrackIsolation = ac.GetDouble("MinTrackIsolation"); //SS specific parameters Double_t pHardBuriedFrac = ac.GetDouble("HardBuriedFrac"); Double_t pSoftBuriedFrac = ac.GetDouble("SoftBuriedFrac"); Double_t pDCosUVCut = ac.GetDouble("DCosUVCut"); Double_t pTrkGapFracCut = ac.GetDouble("TrkGapFracCut"); Double_t pTrkAsymCut = ac.GetDouble("TrkAsymCut"); Double_t pTrkXTalkFracCut = ac.GetDouble("TrkXTalkFracCut"); Double_t pTrkXTalkDef = ac.GetDouble("TrkXTalkDef"); Double_t pShwXTalkFracCut = ac.GetDouble("ShwXTalkFracCut"); Double_t pShwXTalkDef = ac.GetDouble("ShwXTalkDef"); Int_t useChopper = ac.GetInt("UseChopper"); MSG("EventSS",Msg::kDebug) << pHardBuriedFrac << " " << pSoftBuriedFrac << " " << pDCosUVCut << " " << pTrkGapFracCut << " " << pTrkAsymCut << " " << pTrkXTalkFracCut << " " << pTrkXTalkDef << " " << pShwXTalkFracCut << " " << pShwXTalkDef << endl; // grab all needed input lists const CandRecord *candrec = cx.GetCandRecord(); assert(candrec); const VldContext *vldcptr = candrec->GetVldContext(); assert(vldcptr); VldContext vldc = *vldcptr; const CandSliceListHandle *slicelist = 0; CandShowerListHandle *showerlist = 0; const CandTrackListHandle *tracklist = 0; CandSubShowerSRListHandle *subshowerlist = 0; const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn()); for (Int_t i=0; i<=cxin->GetLast(); i++) { TObject *tobj = cxin->At(i); if (tobj->InheritsFrom("CandSliceListHandle")) { slicelist = dynamic_cast<CandSliceListHandle*>(tobj); } if (tobj->InheritsFrom("CandShowerListHandle")) { showerlist = dynamic_cast<CandShowerListHandle*>(tobj); } if (tobj->InheritsFrom("CandTrackListHandle")) { tracklist = dynamic_cast<CandTrackListHandle*>(tobj); } if (tobj->InheritsFrom("CandSubShowerSRListHandle")) { subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj); } } if (!slicelist) { MSG("EventSS",Msg::kError) << "CandSliceListHandle missing\n"; return; } CandShowerHandleItr * showerItr=0; // create iterators for track, shower, and cluster lists, keying on slice if (showerlist && showerlist->GetNDaughters()>0) { showerItr = new CandShowerHandleItr(showerlist->GetDaughterIterator()); CandShowerHandleKeyFunc *showerKf = showerItr->CreateKeyFunc(); showerKf->SetFun(CandShowerHandle::KeyFromSlice); showerItr->GetSet()->AdoptSortKeyFunc(showerKf); showerKf = 0; } CandTrackHandleItr * trackItr = 0; if (tracklist && tracklist->GetNDaughters()>0) { trackItr = new CandTrackHandleItr(tracklist->GetDaughterIterator()); CandTrackHandleKeyFunc *trackKf = trackItr->CreateKeyFunc(); trackKf->SetFun(CandTrackHandle::KeyFromSlice); trackItr->GetSet()->AdoptSortKeyFunc(trackKf); trackKf = 0; } CandSubShowerSRHandleItr * subshowerItr = 0; if (subshowerlist && subshowerlist->GetNDaughters()>0) { subshowerItr = new CandSubShowerSRHandleItr(subshowerlist-> GetDaughterIterator()); CandSubShowerSRHandleKeyFunc *subshowerKf = subshowerItr->CreateKeyFunc(); subshowerKf->SetFun(CandSubShowerSRHandle::KeyFromSlice); subshowerItr->GetSet()->AdoptSortKeyFunc(subshowerKf); subshowerKf = 0; } Int_t nslice=0; // loop over slices CandSliceHandleItr sliceItr(slicelist->GetDaughterIterator()); while (CandSliceHandle *slice = dynamic_cast<CandSliceHandle*>(sliceItr())) { //loop over slices ++nslice; MSG("EventSS", Msg::kDebug) << " ****** Slice " << nslice << " *******" << endl; //TObjArray newevent; TObjArray recolist; TObjArray eventlist; //construct list of all showers and tracks in this slice FillRecoList(slice,showerItr,trackItr,recolist); //first try to remove tracks buried in showers std::vector<Bool_t> removeMe(recolist.GetLast()+1); for (Int_t i=0;i<=recolist.GetLast();i++) { removeMe[i] = false; if(recolist.At(i)->InheritsFrom("CandTrackHandle")) { CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(i)); Bool_t removeTrack = false; for (Int_t j=0;j<=recolist.GetLast();j++) { if (recolist.At(j)->InheritsFrom("CandShowerSRHandle")) { CandShowerSRHandle * shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(j)); if(shower->BuriedTrack(track,pMinTrackIsolation)) { removeTrack = true; break; } } } if(removeTrack){ MSG("EventSS",Msg::kDebug) << " REMOVING TRACK" << endl; removeMe[i] = true; } } } for (Int_t i=0;i<=recolist.GetLast();i++) { if(removeMe[i]) recolist.RemoveAt(i); } recolist.Compress(); MSG("EventSS",Msg::kDebug) << "Number of reco objects in slice: " << recolist.GetLast()+1 << endl; //now make matrix to hold match results between all objects. //this also includes buried tracks that have not been removed //NOTE: matching requires both views match; // burial just counts up matching planes in either view; // One view being significantly buried may indicate a bad track, which // should be associated with the shower, not used to make a new event //gmi Change by G. Irwin 28Feb2009 - //gmi For GCC 4.4, change old statement below to model (following that) of //gmi http://www.codepedia.com/1/Cpp2dVector (see "One-line construction") //gmi Original line (doesn't compile in GCC 4.4): //gmi std::vector<std::vector<Int_t> > objectMatch(recolist.GetLast()+1,0); //gmi Model for replacement line: //gmi std::vector<std::vector<int> > vec(maxx,std::vector<int>(maxy,0)); //gmi Replacement line itself: std::vector<std::vector<Int_t> > objectMatch(recolist.GetLast()+1, std::vector<Int_t>(recolist.GetLast()+1,0)); //initialize: for(int i=0;i<=recolist.GetLast();i++){ for(int j=0;j<=recolist.GetLast();j++) { //gmi Replace this line: "objectMatch[i].push_back(0);" by: objectMatch[i][j] = 0; //gmi May be redundant with constructor init } } //use this vector to keep track of objects already assigned std::vector<Bool_t> objectUsed(recolist.GetLast()+1,0); //loop over all objects for(Int_t i=0;i<=recolist.GetLast();i++) { //loop over i objectUsed[i] = false; //initialization CandShowerSRHandle *shower=0; CandTrackHandle *track=0; if (recolist.At(i)->InheritsFrom("CandShowerSRHandle")) shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i)); if (recolist.At(i)->InheritsFrom("CandTrackHandle")) track = dynamic_cast<CandTrackHandle*>(recolist.At(i)); //loop over all later objects only (since matching matrix is symmetric) for(Int_t j=i+1;j<=recolist.GetLast();j++) { //loop over j MSG("EventSS",Msg::kDebug) << "Checking objectMatch between " << i << " and " << j << endl; CandShowerSRHandle *shower2=0; CandTrackHandle *track2=0; if (recolist.At(j)->InheritsFrom("CandShowerSRHandle")) shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(j)); if (recolist.At(j)->InheritsFrom("CandTrackHandle")) track2 = dynamic_cast<CandTrackHandle*>(recolist.At(j)); if(shower) { if(shower2) { if(!shower->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef) && !shower2->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) { if(shower->BelongsWithShower(shower2,ac,vldcptr,pShwShwDtpos2, pShwShwDz,pShwShwDt)) objectMatch[i][j] = 1; } objectMatch[j][i] = objectMatch[i][j]; } else if(track2) { if(shower->BelongsWithTrack(track2,ac,vldcptr,pShwTrkDtpos2, pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1; if(objectMatch[i][j]!=1) { //check track isn't buried: Int_t *stats = new Int_t[8]; shower->BuriedTrack(track2,1000,stats); //1000 forces all calculations if(stats[3]>0) { //at least one buried track planes Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1; MSG("EventSS",Msg::kDebug) << "BuriedTrack Stats: " << "\n" << "Number of track planes = " << nplanes << "\n" << "nSharedPlanes = " << stats[0] << "\n" << "nMissingIsoPlanes = " << stats[1] << "\n" << "nMissingSharedPlanes = " << stats[2] << "\n" << "nBuriedPlanes = " << stats[3] << "\n" << "nPhysBuriedPlanes = " << stats[4] << "\n" << "nTrkLike = " << stats[5] << "\n" << "nTrkLikeTopol = " << stats[6] << "\n" << "nXTalk = " << stats[7] << endl; nplanes -= (stats[1] + stats[2]); //subtract missing track planes //if more than half of the track planes are in the shower core: match if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac) objectMatch[i][j] = -1; //else if more than half the track planes are anywhere in //the shower or are x-talk like and are not "trklike": match else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > pHardBuriedFrac) objectMatch[i][j] = -1; //check whether track passes through shower at high angle in either view //expect significant number of buried planes else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && ( TMath::Abs(track2->GetDirCosU())>pDCosUVCut || TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) { objectMatch[i][j] = -1; } } delete [] stats; } objectMatch[j][i] = objectMatch[i][j]; } } if(track) { if(shower2) { if(shower2->BelongsWithTrack(track,ac,vldcptr,pShwTrkDtpos2, pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1; if(objectMatch[i][j]!=1) { //check track isn't buried: Int_t *stats = new Int_t[8]; shower2->BuriedTrack(track,1000,stats); //1000 forces all calculations if(stats[3]>0) { //at least one buried track planes Int_t nplanes = track->GetEndPlane() - track->GetBegPlane() + 1; MSG("EventSS",Msg::kDebug) << "BuriedTrack Stats: " << "\n" << "Number of track planes = " << nplanes << "\n" << "nSharedPlanes = " << stats[0] << "\n" << "nMissingIsoPlanes = " << stats[1] << "\n" << "nMissingSharedPlanes = " << stats[2] << "\n" << "nBuriedPlanes = " << stats[3] << "\n" << "nPhysBuriedPlanes = " << stats[4] << "\n" << "nTrkLike = " << stats[5] << "\n" << "nTrkLikeTopol = " << stats[6] << "\n" << "nXTalk = " << stats[7] << endl; nplanes -= (stats[1] + stats[2]); //subtract missing track planes //if more than half of the track planes are in the shower core: match if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac) objectMatch[i][j] = -1; //else if more than half the track planes are anywhere in //the shower or are x-talk like and are not "trklike": match else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > pHardBuriedFrac) objectMatch[i][j] = -1; //check whether track passes through shower at high angle in either view //expect significant number of buried planes else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) { objectMatch[i][j] = -1; } } delete [] stats; } objectMatch[j][i] = objectMatch[i][j]; } else if(track2) { if(track->BelongsWithTrack(track2,ac,vldcptr,pTrkTrkDtpos2, pTrkTrkDz,pTrkTrkDt)) objectMatch[i][j] = 1; objectMatch[j][i] = objectMatch[i][j]; } } MSG("EventSS",Msg::kDebug) << "Match for objects " << i << " and " << j << " : " << objectMatch[i][j] << endl; } } MSG("EventSS",Msg::kDebug) << "Starting to form events" << endl; //start to make events: Int_t nRecoObjects = recolist.GetLast()+1; //get total charge from each object std::vector<Float_t> ph(nRecoObjects,0); for(int i=0;i<nRecoObjects;i++){ if (recolist.At(i)->InheritsFrom("CandShowerSRHandle")) { CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i)); ph[i] = shower->GetCharge(CalStripType::kMIP); } if (recolist.At(i)->InheritsFrom("CandTrackHandle")) { CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(i)); //ph[i] = track->GetCharge(CalStripType::kMIP)+1e7; //force tracks to come first //force tracks to come first, order tracks by number of strips ph[i] = track->GetNStrip()+1e7; } } Int_t nEvents = 0; Int_t nUnusedObjects = nRecoObjects; while(nUnusedObjects>0) { //loop until all reco objects are in events //find largest object: Float_t maxPH = 0; Int_t maxPHIndex = -1; for(int i=0;i<nRecoObjects;i++){ if(ph[i]>maxPH) { maxPH = ph[i]; maxPHIndex = i; } } //add largest objects first: TObjArray *event = new TObjArray; event->Add(recolist.At(maxPHIndex)); ph[maxPHIndex] = 0; objectUsed[maxPHIndex] = true; nUnusedObjects-=1; MSG("EventSS",Msg::kDebug) << "Added object " << maxPHIndex << " to event " << nEvents << endl; if(recolist.At(maxPHIndex)->InheritsFrom("CandTrackHandle")) { MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a track" << endl; //maxPH object is a track CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(maxPHIndex)); //now add objects clearly associated with largest object: for(int i=0;i<nRecoObjects;i++){ if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) { if(recolist.At(i)->InheritsFrom("CandTrackHandle")) { MSG("EventSS",Msg::kDebug) << "Track match. Adding object " << i << " to event " << nEvents << endl; //if this match is a track, add it to event event->Add(recolist.At(i)); ph[i] = 0; objectUsed[i] = true; nUnusedObjects-=1; } // end if for CandTrackHandle else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) { MSG("EventSS",Msg::kDebug) << "Shower match. Testing whether to add object " << i << " to event " << nEvents << endl; CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i)); //if this match is a shower, check where shower is relative to track vertex //if(TMath::Abs(shower2->GetVtxZ()-track->GetVtxZ())<pShwTrkDz/2.){ if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) - track->GetBegPlane(PlaneView::kU)) < pShwTrkDz/(PITCHINMETRES*2.) || TMath::Abs(shower2->GetBegPlane(PlaneView::kV) - track->GetBegPlane(PlaneView::kV)) < pShwTrkDz/(PITCHINMETRES*2.)){ MSG("EventSS",Msg::kDebug) << "Shower " << i << " is close to track vertex. " << "Adding to event " << nEvents << endl; event->Add(recolist.At(i)); ph[i] = 0; objectUsed[i] = true; nUnusedObjects-=1; //look for objects associated with this vertex shower for(int j=0;j<nRecoObjects;j++){ if(j!=maxPHIndex && objectMatch[i][j]!=0 && !objectUsed[j]) { //got a match if(objectMatch[maxPHIndex][j]==0) { // j does not match maxPHIndex if(recolist.At(j)->InheritsFrom("CandTrackHandle")) { MSG("EventSS",Msg::kDebug) << "Shower " << i << " has another track match. " << "Test how buried track " << j << " is in shower" << endl; if(objectMatch[i][j]==-1){ event->Add(recolist.At(j)); ph[j] = 0; objectUsed[j] = true; nUnusedObjects-=1; } else { CandTrackHandle *track2 = dynamic_cast<CandTrackHandle*>(recolist.At(j)); //check if track2 is buried Int_t *stats = new Int_t[8]; shower2->BuriedTrack(track2,1000,stats); if(stats[3]>0){ Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1; MSG("EventSS",Msg::kDebug) << "BuriedTrack Stats: " << "\n" << "Number of track planes = " << nplanes << "\n" << "nSharedPlanes = " << stats[0] << "\n" << "nMissingIsoPlanes = " << stats[1] << "\n" << "nMissingSharedPlanes = " << stats[2] << "\n" << "nBuriedPlanes = " << stats[3] << "\n" << "nPhysBuriedPlanes = " << stats[4] << "\n" << "nTrkLike = " << stats[5] << "\n" << "nTrkLikeTopol = " << stats[6] << "\n" << "nXTalk = " << stats[7] << endl; nplanes -= (stats[1] + stats[2]); if( (Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac) || (Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > pHardBuriedFrac) || (stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) { MSG("EventSS",Msg::kDebug) << "Adding track " << j << " to event " << nEvents << endl; event->Add(recolist.At(j)); ph[j] = 0; objectUsed[j] = true; nUnusedObjects-=1; } } delete [] stats; } } else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) { MSG("EventSS",Msg::kDebug) << "Shower " << i << " has another shower match. " << "Test whether shower " << j << " has a better association" << endl; Bool_t isAssocTrack = false; //is this shower assoc with a track for(int k=0;k<nRecoObjects;k++) { //check for match if(k!=maxPHIndex && k!=i && objectMatch[j][k]!=0){ //check for track if(recolist.At(k)->InheritsFrom("CandTrackHandle")) //check it does not match with shower2 if(objectMatch[i][k]==0) isAssocTrack = true; } } if(!isAssocTrack) { //no track associated MSG("EventSS",Msg::kDebug) << "Shower " << j << " has no track association. " << "Add to event " << nEvents << endl; event->Add(recolist.At(j)); ph[j] = 0; objectUsed[j] = true; nUnusedObjects-=1; } } } else { //j matches maxPHIndex, so add it to event MSG("EventSS",Msg::kDebug) << "Shower " << j << " matches orginal track " << maxPHIndex << " Add to event " << nEvents << endl; event->Add(recolist.At(j)); ph[j] = 0; objectUsed[j] = true; nUnusedObjects-=1; } } } } //end if for vertex shower else { //if not a vertex shower, check other associations MSG("EventSS",Msg::kDebug) << "Testing whether shower " << i << " has better associations. " << endl; Bool_t otherTrackVertexMatch = false; Bool_t otherTrackMatch = false; Bool_t otherShowerMatch = false; for(int j=0;j<nRecoObjects;j++){ if(j!=maxPHIndex && objectMatch[i][j]!=0 && !objectUsed[j] && objectMatch[maxPHIndex][j]==0) { if(recolist.At(j)->InheritsFrom("CandTrackHandle")) { CandTrackHandle *track2 = dynamic_cast<CandTrackHandle*>(recolist.At(j)); if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) - track->GetBegPlane(PlaneView::kU)) < pShwTrkDz/(PITCHINMETRES*2.) || TMath::Abs(shower2->GetBegPlane(PlaneView::kV) - track->GetBegPlane(PlaneView::kV)) < pShwTrkDz/(PITCHINMETRES*2.) ) otherTrackVertexMatch = true; else { //first test whether track2 is buried //if so then can still associate shower to track if(objectMatch[i][j]!=-1) { //object match is not flagged as buried //perform check: Bool_t track2Buried = false; Int_t *stats2 = new Int_t[8]; shower2->BuriedTrack(track2,1000,stats2); if(stats2[3]>0){ Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1; MSG("EventSS",Msg::kDebug) << "BuriedTrack Stats: " << "\n" << "Number of track planes = " << nplanes << "\n" << "nSharedPlanes = " << stats2[0] << "\n" << "nMissingIsoPlanes = " << stats2[1] << "\n" << "nMissingSharedPlanes = " << stats2[2] << "\n" << "nBuriedPlanes = " << stats2[3] << "\n" << "nPhysBuriedPlanes = " << stats2[4] << "\n" << "nTrkLike = " << stats2[5] << "\n" << "nTrkLikeTopol = " << stats2[6] << "\n" << "nXTalk = " << stats2[7] << endl; nplanes -= (stats2[1] + stats2[2]); if( (Float_t(stats2[4])/Float_t(nplanes) > pHardBuriedFrac) || (Float_t(stats2[3]-stats2[5]-stats2[6])/Float_t(nplanes) > pHardBuriedFrac) || (stats2[3]>=Int_t(nplanes*pSoftBuriedFrac) && (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) { track2Buried = true; objectMatch[i][j]=-1; objectMatch[j][i]=-1; } } if(!track2Buried){ //test whether track2 matches better than track Int_t *stats = new Int_t[8]; shower2->BuriedTrack(track,1000,stats); if( ( stats2[3] - stats2[5] - stats2[6] ) >= ( stats[3] - stats[5] - stats[6] ) && ( stats2[4] >= stats[4] ) ) { otherTrackMatch = true; } delete [] stats; } delete [] stats2; } } } else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) { Int_t *stats = new Int_t[8]; shower2->BuriedTrack(track,1000,stats); if(stats[4]==0) { //how convincing is first track/shower match? otherShowerMatch = true; } delete [] stats; } } } //if there are no other viable matches: if(otherTrackVertexMatch || otherShowerMatch || otherTrackMatch) { MSG("EventSS",Msg::kDebug) << "Shower " << i << " has better association: " << " otherTrackVertexMatch = " << otherTrackVertexMatch << " otherShowerMatch = " << otherShowerMatch << " otherTrackMatch = " << otherTrackMatch << endl; continue; } //then add this shower to this track MSG("EventSS",Msg::kDebug) << "Shower " << i << " has no better association. " << "Add to event " << nEvents << endl; event->Add(recolist.At(i)); ph[i] = 0; objectUsed[i] = true; nUnusedObjects-=1; } //end if for non-vertex shower } //end if for CandShowerSRHandle } //end if for valid association } //end for loop over all objects } //end if maxPH object is a track else if(recolist.At(maxPHIndex)->InheritsFrom("CandShowerSRHandle")) { MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a shower" << endl; MSG("EventSS",Msg::kDebug) << "Check whether shower " << maxPHIndex << " has other associations" << endl; //maxPH object is a shower //CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(maxPHIndex)); //now add objects clearly associated with largest object: for(int i=0;i<nRecoObjects;i++){ if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) { if(recolist.At(i)->InheritsFrom("CandTrackHandle")) { MSG("EventSS",Msg::kWarning) << "Arrgh - Logic flaw in eventSR. " << "Seeing tracks in event building, " << "when all tracks should already be " << "assigned!" << endl; } //end if CandTrackHandle else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) { //check for two more levels of associated showers //CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i)); MSG("EventSS",Msg::kDebug) << "Adding shower " << i << " to event " << nEvents << endl; MSG("EventSS",Msg::kDebug) << "Check whether shower " << i << " has other associations" << endl; event->Add(recolist.At(i)); ph[i] = 0; objectUsed[i] = true; nUnusedObjects-=1; for(int j=0;j<nRecoObjects;j++){ if(j!=maxPHIndex && objectMatch[i][j]!=0 && !objectUsed[j] && objectMatch[maxPHIndex][j]==0) { if(recolist.At(j)->InheritsFrom("CandTrackHandle")) { MSG("EventSS",Msg::kWarning) << "Arrgh again - Logic flaw in eventSR. " << "Seeing tracks in event building, " << "when all tracks should already be " << "assigned!" << endl; } else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) { MSG("EventSS",Msg::kDebug) << "Adding shower " << j << " to event " << nEvents << endl; event->Add(recolist.At(j)); ph[j] = 0; objectUsed[j] = true; nUnusedObjects-=1; } } } } //end if CandShowerSRHandle } //end if valid association } //end loop over all objects } //end if maxPH object is a shower //add event to eventlist eventlist.Add(event); nEvents+=1; } // remove all empty events from event list for (Int_t i=0; i<=eventlist.GetLast(); i++) { TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i)); if (objectlist->GetLast()<0) eventlist.Remove(objectlist); else if(objectlist->GetLast()==0) { //consider single shower/track events if(objectlist->At(0)->InheritsFrom("CandShowerSRHandle")) { MSG("EventSS",Msg::kDebug) << "Have a single candshower event " << endl; CandShowerSRHandle *showerSR = dynamic_cast<CandShowerSRHandle*> (objectlist->At(0)); //if shower seems unphysical, remove it if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) { eventlist.Remove(objectlist); MSG("EventSS",Msg::kDebug) << "REMOVING SHOWER-ONLY EVENT" << endl; continue; } } else if(objectlist->At(0)->InheritsFrom("CandTrackHandle")) { CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(0)); if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) { //if track seems unphysical, remove it eventlist.Remove(objectlist); MSG("EventSS",Msg::kDebug) << "REMOVING TRACK-ONLY EVENT" << endl; continue; } } } else if(objectlist->GetLast()==1) { //multiple shower/track events for(int ind_loop = 0;ind_loop<2;ind_loop++){ Int_t ind_loop1 = ind_loop; Int_t ind_loop2 = ind_loop+1; if(ind_loop2==2) ind_loop2 = 0; if(objectlist->At(ind_loop1)->InheritsFrom("CandShowerSRHandle")){ CandShowerSRHandle *showerSR = dynamic_cast<CandShowerSRHandle*> (objectlist->At(ind_loop1)); if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) { if(objectlist->At(ind_loop2)->InheritsFrom("CandTrackHandle")){ CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(ind_loop2)); if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) { eventlist.Remove(objectlist); MSG("EventSS",Msg::kDebug) << "REMOVING UNPHYSICAL SHOWER+TRACK EVENT" << endl; continue; } } } } } } if(useChopper==1) { std::vector<CandHandle> evtVec; for(int i=0;i<=objectlist->GetLast();i++){ CandHandle *event = dynamic_cast<CandHandle*>(objectlist->At(i)); evtVec.push_back(*event); } ChopHelper chopper(cx.GetMom()); ChopHelp *chop = chopper.GetChopHelp(evtVec); if(chop->nchop==1) continue; else { chop->Print(); //do something with chop info //for(int i=0;i<ncandidates;i++){ //} } } } eventlist.Compress(); // now make CandEvents for each event we have found Bool_t buildEvent=false; for (Int_t i=0; i<=eventlist.GetLast(); i++) { TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i)); MSG("EventSS",Msg::kDebug) << "making event of " << objectlist->GetLast()+1 << " objects \n"; cxx.SetDataIn(objectlist); buildEvent=true; CandEventHandle eventhandle = CandEvent::MakeCandidate(ah,cxx); MSG("EventSS",Msg::kDebug) << " # of strips: " << eventhandle.GetNDaughters() << "\n"; eventhandle.SetCandSlice(slice); ch.AddDaughterLink(eventhandle); delete eventlist.At(i); } if (!buildEvent) BuildEventFromUnassoc(ac,ch,cx,slice); } // end loop over slice list delete trackItr; delete showerItr; delete subshowerItr; MSG("EventSS",Msg::kDebug) << "starting unassociated hits assignment \n"; // find unassociated hits, and place them in tobjarray sliceItr.Reset(); TObjArray unassociated; CandEventListHandle &eventlist = dynamic_cast<CandEventListHandle &>(ch); CandEventHandleItr * eventItr2 = 0; if (eventlist.GetNDaughters()>0) { eventItr2 = new CandEventHandleItr(eventlist.GetDaughterIterator()); } FindUnassociated(ac, sliceItr,eventItr2,unassociated); delete eventItr2; SetPrimaryShowers(ch,ac); // now loop over unnassociated hits, and check for adjacency with all // showers, and track vertices (in the case that no shower exists). // If a match is found in the former case, add strip to shower. // In the latter case, create a new shower. // this map contains the separation between each event/unassociated hit vector<map<const CandEventHandle*, Double_t> > dist2map(unassociated.GetLast()+1); FillDist2Map(ac,unassociated,ch,dist2map); // initialize a vector which holds a Bool_t // indicating whether a hit has already been used vector<Bool_t> alreadyUsed(unassociated.GetLast()+1); for (Int_t i=0; i<=unassociated.GetLast(); i++) alreadyUsed[i] = false; // loop over unnassociated hits, adding them to events. // After each loop the affected distance map entries are recalculated, // and the loop is repeated. This continues until a loop adds no strips. Bool_t found=true; Int_t n_found=0; while (found) { //found indicates whether a strip was added to an event in this loop found=false; CandStripHandle *closeststrip= 0; CandEventHandle *closestevent= 0; Double_t closestdist2 = 9999999.; Int_t closesti=-1; // loop over unassoc strips for (Int_t i=0; i<=unassociated.GetLast(); i++) { if(alreadyUsed[i]) continue; CandStripHandle *strip = dynamic_cast<CandStripHandle*>(unassociated.At(i)); CandEventHandle *bestevent = 0; Double_t bestdist2 = 9999999.; // loop over events, and find event which is closest to strip // and in time TIter eventItr(ch.GetDaughterIterator()); while (CandEventHandle *event = dynamic_cast<CandEventHandle*>(eventItr())) { Double_t dtime = strip->GetBegTime()-event->GetVtxT(); if ( dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1 && (!bestevent || dist2map[i][event]<bestdist2) ) { bestevent = event; bestdist2 = dist2map[i][event]; } } if (bestevent && (!closeststrip || bestdist2<closestdist2)) { closeststrip = strip; closesti=i; closestdist2 = bestdist2; closestevent = bestevent; } } MSG("EventSS",Msg::kVerbose)<< "closest strip " << closestdist2 << "/"<<pHitAssocMaxDist2<< endl; // if distance to closest event less than maximum allowable add to event if(closeststrip && closestevent && closestdist2<pHitAssocMaxDist2){ MSG("EventSS",Msg::kVerbose) << "closest strip " << closeststrip->GetPlane() << " " << closeststrip->GetTPos() << " " << closestdist2 << "\n"; alreadyUsed[closesti] = true; found = true; ++n_found; // now either add this strip to existing primary shower, or // create a new shower at the track vertex if one does not exist // Get PrimaryShower CandHandle from fShowerList for modification. CandShowerHandle *shower = 0; if (CandShowerHandle *primsh = closestevent->GetPrimaryShower()) { for (Int_t ishower=0; !shower && ishower<closestevent->GetLastShower()+1; ishower++) { const CandShowerHandle *target = closestevent->GetShower(ishower); if (primsh->IsCloneOf(*target)) // Tests clone or Cand == shower = closestevent->GetShowerWritable(ishower); } } Float_t minShwPlane=0; Float_t maxShwPlane=0; if(shower){ minShwPlane=min(shower->GetVtxPlane(),shower->GetEndPlane()); maxShwPlane=max(shower->GetVtxPlane(),shower->GetEndPlane()); } Float_t vertexsep = 0; CandTrackHandle *track = closestevent->GetPrimaryTrack(); if (shower && track) { vertexsep = abs(shower->GetBegPlane() - track->GetBegPlane()); vertexsep *= PITCHINMETRES; MSG("EventSS",Msg::kVerbose) << "Found track and primary shower. \n" << "Vertex separation (limit) = " << vertexsep << "(" << pShwTrkDz << ")" << " \n"; } // if event has shower and either no track, or shower is close // to track vertex, add this hit to the shower which was // compared in position to this strip earler. if (shower && (!track || vertexsep<pShwTrkDz)) { // NOTE: CandShowerHandle *shower, used to modify CandShower, is owned // by fShowerList. // add new strip to this shower, and the appropriate cluster MSG("EventSS",Msg::kVerbose) << "Shower min/max plane (limit): " << minShwPlane << "/" << maxShwPlane << " (" << pMaxNewShwLen << ")" << "\n"; if(closeststrip->GetPlane()>minShwPlane-pMaxNewShwLen && closeststrip->GetPlane()<maxShwPlane+pMaxNewShwLen ){ MSG("EventSS",Msg::kVerbose) << "adding strip to shower \n"; AddStripToEvent(closestevent,shower,subshowerlist,closeststrip); } } // no shower consistent with track vertex, and this strip // is close to vertex, so start new shower. else if (track) CreatePrimaryShower(ac,cxx,closestevent,showerlist, subshowerlist,closeststrip); // need to recalculate distances bewteen remaining // unassoc hits and this modified event, if the added strip // results in a smaller separation, replace value in dist2map ReFillDist2Map(ac,closeststrip,closestevent,unassociated,dist2map); } // if adding strip to event in form of new shower } // if found (if an unassoc hit was matched to event in last loop over them // because we have added strips to previously constructed showers, their // runalg method should be run on them again. In the future, an // optimization would be to do this only on showers which have been // touched. MSG("EventSS",Msg::kVerbose) << "added " << n_found << " of " << unassociated.GetLast()+1 << " unassociated hits " << endl; cxx.SetDataIn(showerlist); // Pass along persistable CandShowerList ReConstructShowers(ac,ch,cxx); //Final step: use the chopper to combine events across slices: if(useChopper==1) { std::vector<CandHandle> evtVec; CandEventHandleItr evitr(ch.GetDaughterIterator()); while(CandEventHandle *event = dynamic_cast<CandEventHandle*>(evitr())){ evtVec.push_back(*event); } ChopHelper chopper(cx.GetMom()); ChopHelp *chop = chopper.GetChopHelp(evtVec); chop->Print(); } // For each event, set primary shower based on // distance from vertex and energy SetPrimaryShowers(ch,ac); }

void AlgEventSSList::SetPrimaryShowers (  CandHandle &  ch, AlgConfig &  ac )

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 1380 of file AlgEventSSList.cxx. References CandHandle::GetDaughterIterator(), Registry::GetDouble(), and CandEventHandle::SetPrimaryShower(). Referenced by RunAlg(). { // This method loops over events in the eventlist, finding the largest // shower and setting that to the primary shower in the event. Double_t pShwShwDz = ac.GetDouble("ShwShwDz"); Double_t pMinShwEFract = ac.GetDouble("MinShwEFract"); TIter evItr(ch.GetDaughterIterator()); while (CandEventHandle *ev = dynamic_cast<CandEventHandle*>(evItr())) { ev->SetPrimaryShower(pMinShwEFract,pShwShwDz); } }

void AlgEventSSList::Trace (  const char *  c  )  const [virtual]

Reimplemented from AlgBase. Definition at line 2033 of file AlgEventSSList.cxx. { }

---

The documentation for this class was generated from the following files:

• AlgEventSSList.h
• AlgEventSSList.cxx

---