AlgEventSSList Class Reference

#include <AlgEventSSList.h>

Inheritance diagram for AlgEventSSList:

AlgBase 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)

Detailed Description

Definition at line 31 of file AlgEventSSList.h.


Constructor & Destructor Documentation

AlgEventSSList::AlgEventSSList (  ) 

Definition at line 62 of file AlgEventSSList.cxx.

00063 {
00064 }

AlgEventSSList::~AlgEventSSList (  )  [virtual]

Definition at line 67 of file AlgEventSSList.cxx.

00068 {
00069 }


Member Function Documentation

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

Definition at line 1836 of file AlgEventSSList.cxx.

References Msg::kVerbose, and MSG.

01837                                                                            {
01838   if (!merge) {
01839     // if this is the first match we have found for this object, 
01840     //add object to list
01841     // and point prevlist to this object list
01842     objectlist->Add(reco);
01843     MSG("EventSS",Msg::kVerbose) << "    adding to list\n";
01844   }
01845   else {                  
01846     MSG("EventSS",Msg::kVerbose) << "    combining lists\n";
01847     // this object matches more than one event - merge their member
01848     // objects to list prevlist, and add the new object
01849     for (Int_t ireco2=0; ireco2<=objectlist->GetLast();ireco2++)
01850       prevlist->Add(objectlist->At(ireco2));
01851     prevlist->Add(reco);
01852     objectlist->Clear();
01853   }
01854 }

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(), CandRecoHandle::GetCharge(), OscFit::GetCharge(), CandHandle::GetDaughterIterator(), CandSubShowerSRHandle::GetEnergy(), CandShowerSRHandle::GetLastSubShower(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandShowerSRHandle::GetSubShower(), CandStripHandle::GetTPos(), CandHandle::GetUidInt(), CandStripHandle::GetZPos(), Calibrator::Instance(), CalStripType::kSigCorr, CalDigitType::kSigCorr, PlaneView::kU, PlaneView::kV, Msg::kVerbose, Msg::kWarning, PlaneView::kX, PlaneView::kY, MSG, CandHandle::RemoveDaughter(), CandShowerSRHandle::RemoveSubShower(), CandRecoHandle::SetCandSlice(), and CandSubShowerSRHandle::SetEnergy().

Referenced by RunAlg().

00980 {
00981   MSG("EventSS",Msg::kVerbose) << "In AddStripToEvent:" << endl;
00982   closestevent->AddDaughterLink(*closeststrip);
00983   shower->AddDaughterLink(*closeststrip);
00984   if(subshowerlist && shower->InheritsFrom("CandShowerSRHandle")) {
00985     Calibrator& calibrator=Calibrator::Instance();
00986     Int_t foundsubshower=0;
00987     Int_t closestSubShower=-1;
00988     Double_t BestStripDistance = 99999;
00989     Double_t StripDistance = 99999;
00990     CandShowerSRHandle *showerSR = 
00991       dynamic_cast<CandShowerSRHandle*> (shower);
00992     MSG("EventSS",Msg::kVerbose) << "Shower UID = " 
00993                                  << shower->GetUidInt() << endl;
00994 
00995     for (Int_t isubshower=0;
00996          isubshower<showerSR->GetLastSubShower()+1;
00997          isubshower++) {
00998       CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
00999         (showerSR->GetSubShower(isubshower));
01000       //add strip to cluster of the same view
01001       if (ssl->GetPlaneView()==closeststrip->GetPlaneView()) {
01002         //decide which is most appropriate subshower to add strip to
01003         //based on distance to/spread of subshower:
01004         Double_t tposVtx = 0;
01005         if(ssl->GetPlaneView()==PlaneView::kU || 
01006            ssl->GetPlaneView()==PlaneView::kX) 
01007           tposVtx = ssl->GetVtxU();
01008         else if(ssl->GetPlaneView()==PlaneView::kV || 
01009                 ssl->GetPlaneView()==PlaneView::kY)
01010           tposVtx = ssl->GetVtxV();
01011         StripDistance = 99999;
01012         if(ssl->GetAvgDev()>0) {
01013           StripDistance = 
01014             (TMath::Abs(closeststrip->GetTPos() - 
01015                         (tposVtx + ssl->GetSlope() * 
01016                          (closeststrip->GetZPos() - 
01017                           ssl->GetVtxZ()))) * 
01018              TMath::Cos(TMath::ATan(ssl->GetSlope()))) * 
01019             calibrator.GetMIP(closeststrip->
01020                               GetCharge(CalDigitType::kSigCorr)) / 
01021             ssl->GetAvgDev();
01022         }
01023         if(StripDistance>9999){ //bad slope
01024           StripDistance = 9999 + 
01025             TMath::Sqrt(TMath::Power(TMath::Abs(closeststrip->GetZPos() - 
01026                                                 ssl->GetVtxZ()),2) +
01027                         TMath::Power(TMath::Abs(closeststrip->GetTPos() - 
01028                                                 tposVtx),2));
01029         }
01030         if(StripDistance<BestStripDistance || foundsubshower==0){
01031           foundsubshower=1;
01032           BestStripDistance = StripDistance;
01033           closestSubShower = isubshower;
01034         }
01035       }
01036     }
01037     if(foundsubshower) {
01038       MSG("EventSS",Msg::kVerbose) << "adding strip to subshower "
01039                                    << closeststrip->GetPlane() 
01040                                    << endl;
01041       CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
01042         (showerSR->GetSubShower(closestSubShower));
01043       //extract mip -> gev conversion:
01044       Double_t GeVperMIP = 0;
01045       if(ssl->GetEnergy()!=0) GeVperMIP = ssl->GetEnergy()/
01046         Double_t(calibrator.GetMIP(ssl->GetCharge(CalStripType::kSigCorr)));
01047       TIter csshItr(subshowerlist->GetDaughterIterator());
01048       while (CandSubShowerSRHandle *cssh = 
01049              dynamic_cast<CandSubShowerSRHandle*>(csshItr())) { 
01050         if(cssh->IsEquivalent(ssl) && ssl->GetUidInt()==cssh->GetUidInt()) {
01051           MSG("EventSS",Msg::kVerbose) << "SubShower (old,new) UID " 
01052                                        << ssl->GetUidInt()
01053                                        << " " << cssh->GetUidInt() << endl; 
01054           showerSR->RemoveSubShower(ssl);
01055           CandSubShowerSRHandle *newss = cssh->DupHandle();
01056           newss->SetCandSlice(cssh->GetCandSlice());
01057           newss->AddDaughterLink(*closeststrip);
01058           newss->SetEnergy(GeVperMIP*Double_t(calibrator.GetMIP(closeststrip->
01059                  GetCharge(CalDigitType::kSigCorr))) + newss->GetEnergy());
01060           subshowerlist->RemoveDaughter(cssh);
01061           subshowerlist->AddDaughterLink(*newss);
01062           showerSR->AddSubShower(newss);
01063           delete newss;
01064           break;
01065         }
01066       }
01067     }
01068     else {
01069       MSG("EventSS",Msg::kWarning) << "Strip added to shower but not to subshower!"
01070                                    << " StripDistance = " << StripDistance
01071                                    << endl;
01072       MSG("EventSS",Msg::kWarning) << "Number of subshowers in shower = " 
01073                                    << showerSR->GetLastSubShower()+1
01074                                    << endl;
01075       for(int i=0;i<showerSR->GetLastSubShower()+1;i++){
01076         CandSubShowerSRHandle *ssl = const_cast<CandSubShowerSRHandle*>
01077           (showerSR->GetSubShower(i));
01078         MSG("EventSS",Msg::kWarning) << "SS " << i << " UID = " 
01079                                      << ssl->GetUidInt() << endl;
01080       }
01081     }
01082   }
01083 }

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 CandHandle::AddDaughterLink(), CandHandle::FindDaughter(), Registry::Get(), CandHandle::GetCandRecord(), CandContext::GetCandRecord(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), VldContext::GetDetector(), Registry::GetDouble(), AlgFactory::GetInstance(), Registry::GetInt(), CandContext::GetMom(), RecMinos::GetVldContext(), CandHandle::GetVldContext(), CandRecoHandle::GetVtxT(), CandRecoHandle::GetVtxU(), CandRecoHandle::GetVtxV(), CandRecoHandle::GetVtxZ(), CandShowerSRHandle::IsUnphysical(), Msg::kDebug, Detector::kNear, CalDigitType::kPE, PlaneView::kU, ClusterType::kUnknown, PlaneView::kV, Msg::kVerbose, CandEvent::MakeCandidate(), CandShowerSR::MakeCandidate(), CandSubShowerSR::MakeCandidate(), MSG, CandHandle::SetCandRecord(), CandEventHandle::SetCandSlice(), CandRecoHandle::SetCandSlice(), CandSubShowerSRHandle::SetClusterID(), and CandEventHandle::SetPrimaryShower().

Referenced by RunAlg().

01403 {
01404   assert(cx.GetDataIn());
01405   MSG("EventSS",Msg::kDebug) <<" Building Event - no showers or tracks " 
01406                              << endl;
01407   if (!(cx.GetDataIn()->InheritsFrom("TObjArray"))) {
01408     MSG("EventSS",Msg::kDebug)<<" dataIn inherits from TObjArray " << endl;
01409     return;
01410   }
01411 
01412   Double_t pShwXTalkFracCut = ac.GetDouble("ShwXTalkFracCut");
01413   Double_t pShwXTalkDef = ac.GetDouble("ShwXTalkDef");
01414 
01415   const CandSliceListHandle *slicelist = 0;
01416   CandShowerListHandle *showerlist = 0;
01417   CandSubShowerSRListHandle *subshowerlist = 0;
01418   const TObjArray *cxin = dynamic_cast<const TObjArray *>(cx.GetDataIn());
01419   for (Int_t i=0; i<=cxin->GetLast(); i++) {
01420     TObject *tobj = cxin->At(i);
01421     if (tobj->InheritsFrom("CandSliceListHandle")) {
01422       slicelist = dynamic_cast<CandSliceListHandle*>(tobj);
01423     }
01424     // Retrieve persistant modifiable CandShowerList
01425     if (tobj->InheritsFrom("CandShowerListHandle")) {
01426       showerlist = dynamic_cast<CandShowerListHandle*>(tobj);
01427     }
01428     if (tobj->InheritsFrom("CandSubShowerSRListHandle")) {
01429       subshowerlist = dynamic_cast<CandSubShowerSRListHandle*>(tobj);
01430     }
01431   }
01432   
01433   if (!showerlist || !slicelist || !subshowerlist) {
01434     MSG("EventSS",Msg::kDebug)
01435       << "Missing either slicelist, showerlist or subshowerlist" << endl;
01436     return;
01437   }
01438   Int_t detectorType = slicelist->GetVldContext()->GetDetector();
01439   
01440   CandContext cxx(this,cx.GetMom());
01441   AlgFactory &af = AlgFactory::GetInstance();
01442   
01443   const char *pEventAlgConfig = 0;
01444   ac.Get("EventAlgConfig",pEventAlgConfig);
01445   AlgHandle ah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig);
01446   
01447   Double_t pHitAssocDPlane = ac.GetInt("HitAssocDPlane");
01448   Double_t pHitAssocDTime0 = ac.GetDouble("HitAssocDTime0");
01449   Double_t pHitAssocDTime1 = ac.GetDouble("HitAssocDTime1");
01450   Double_t pHitAssocMaxDist2 = ac.GetDouble("HitAssocMaxDist2");
01451   Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm");
01452   Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm");
01453   Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm");
01454   Double_t pShwShwDtpos2 = ac.GetDouble("ShwShwDtpos2");
01455   Double_t pShwShwDz = ac.GetDouble("ShwShwDz");
01456   Double_t pShwShwDt = ac.GetDouble("ShwShwDt");
01457   Double_t pMinShwEFract = ac.GetDouble("MinShwEFract");
01458   Double_t pMinUnassocStripPE = ac.GetDouble("MinUnassocStripPE");
01459 
01460 
01461   const CandRecord *candrec = cx.GetCandRecord();
01462   assert(candrec);
01463   const VldContext *vldcptr = candrec->GetVldContext();
01464   assert(vldcptr);
01465   VldContext vldc = *vldcptr;
01466   UgliGeomHandle ugh(vldc);
01467 
01468   // Iterate over strips , finding sets which are 'in time', and within
01469   // a given distance of the nearest neighbor in the set, in both views.
01470   // Each of these sets is then used to construct a shower, which becomes
01471   // the primary shower of a new event.  We continue until an iteration
01472   // in which no new event is created.
01473   
01474   // fill array of unique unassociated strips
01475   TObjArray unassociated;
01476   TIter stripItr(slice->GetDaughterIterator());
01477   while (CandStripHandle *strip =
01478          dynamic_cast<CandStripHandle*>(stripItr())) {
01479 
01480     // Ignore strips below the required minimum PE
01481     if (strip->GetCharge(CalDigitType::kPE)<pMinUnassocStripPE) continue;
01482 
01483     Bool_t found=false;
01484     for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) {
01485       CandStripHandle *strip0 =
01486         dynamic_cast<CandStripHandle*>(unassociated.At(i));
01487       if (*strip == *strip0) {
01488         found = true;
01489       }
01490     }
01491     if (!found) {
01492       if(!(detectorType==Detector::kNear && strip->GetPlane()>120)){
01493         unassociated.Add(strip);
01494       }
01495     }
01496   }
01497   
01498   // loop until no new matches are found in the loop over unassoc hits.
01499   Bool_t found=false;
01500   while (unassociated.GetLast()>0 && !found) {
01501     MSG("EventSS",Msg::kDebug)<< " # unassoc " 
01502                               << unassociated.GetLast()+1 << endl; 
01503     Bool_t firstu=true;
01504     Bool_t firstv=true;
01505     TObjArray neweventu;
01506     TObjArray neweventv;
01507     Bool_t addedstrip = true;
01508     Float_t totcharge=0;
01509     while (addedstrip) {
01510       addedstrip=false;      
01511       // loop over unassoc strips.
01512       for (Int_t i=0; i<=unassociated.GetLast(); i++) {
01513         CandStripHandle *strip =
01514           dynamic_cast<CandStripHandle*>(unassociated.At(i));
01515         switch (strip->GetPlaneView()) {
01516         case PlaneView::kU:
01517           // strip is in U plane
01518           // if this is the first strip, add to newevent strip list
01519           if (firstu) {
01520             // if we don't have v strips yet, just add this strip to u list
01521             // seed the event
01522             if (firstv) {
01523               neweventu.Add(strip);
01524               MSG("EventSS",Msg::kDebug)<< " add  first  " << endl;
01525               totcharge +=strip->GetCharge();
01526               addedstrip = true;
01527               firstu=false;
01528             }
01529             // if we already have one or more v strips check that they are
01530             // in time, and  agree in Z position.
01531             else {
01532               Bool_t found2=false;
01533               for (Int_t j=0; j<=neweventv.GetLast() && !found2;
01534                    j++) {
01535                 CandStripHandle *estrip =
01536                   dynamic_cast<CandStripHandle*>(neweventv.At(j));
01537                 Double_t dtime =
01538                   strip->GetBegTime()-estrip->GetBegTime(); 
01539                 MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 
01540                                              << " " <<pHitAssocDTime0*1e9 
01541                                              << "/" << pHitAssocDTime1*1e9 
01542                                              << " plane diff " 
01543                                              <<  abs(strip->GetPlane() - 
01544                                                      estrip->GetPlane()) 
01545                                              << "/" << pHitAssocDPlane << endl;
01546                 if (dtime>pHitAssocDTime0 &&
01547                     dtime<pHitAssocDTime1 &&
01548                     abs(strip->GetPlane()-estrip->GetPlane()) <=
01549                     pHitAssocDPlane) {
01550                   MSG("EventSS",Msg::kVerbose)<< " add  " << endl;
01551                   // match is found, add this strip to array
01552                   neweventu.Add(strip);
01553                   totcharge +=strip->GetCharge();
01554                   addedstrip = true;
01555                   firstu=0;
01556                   found2=true;
01557                 } // if match found
01558               }  // for loop over V strips in new event
01559             } // end not first V strip
01560           } // end first U strip
01561           else {  
01562             // if we already have u strips, make sure that the one we
01563             // add is adjaccent, using distance criteria based
01564             // on plane, time, and transverse postion separation
01565             Bool_t found2=false;
01566             for (Int_t j=0; j<=neweventu.GetLast() && !found2; j++) {
01567               CandStripHandle *estrip =
01568                 dynamic_cast<CandStripHandle*>(neweventu.At(j));
01569               Double_t dtime =
01570                 strip->GetBegTime()-estrip->GetBegTime(); 
01571               MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 
01572                                            << " " <<pHitAssocDTime0*1e9 
01573                                            << "/" << pHitAssocDTime1*1e9 
01574                                            << " plane diff " 
01575                                            << abs(strip->GetPlane() - 
01576                                                   estrip->GetPlane())
01577                                            << " dtpos " << (strip->GetTPos() - 
01578                                                             estrip->GetTPos())
01579                                            << endl;
01580               if (dtime>pHitAssocDTime0 &&
01581                   dtime<pHitAssocDTime1) {
01582                 Double_t dplane =
01583                   (Double_t)(strip->GetPlane()-estrip->GetPlane());
01584                 Double_t dtpos = strip->GetTPos()-estrip->GetTPos();
01585                 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01586                   (pHitAssocPParm*dtpos*dtpos) +
01587                   (pHitAssocTParm*dtime*dtime);
01588                 // if dist2 2 parameter meets match criteria
01589                 // add strip to array
01590                 if (dist2<pHitAssocMaxDist2*4) {
01591                   MSG("EventSS",Msg::kVerbose)<< " add " << endl;
01592                   neweventu.Add(strip);
01593                   totcharge +=strip->GetCharge();
01594                   addedstrip = true;
01595                   firstu=0;
01596                   found2=true;
01597                 } // if positions match
01598               } // if in time
01599             } // end loop over event u strips
01600           } // end if already have u strip
01601           break;
01602           
01603         case PlaneView::kV:
01604           
01605           // if this is the first strip, add to newevent strip list
01606           if (firstv) {
01607             // strip is in V plane
01608             // if this is the first strip, add to newevent strip list
01609             if (firstu) {
01610               neweventv.Add(strip);
01611               totcharge +=strip->GetCharge();
01612               MSG("EventSS",Msg::kDebug)<< " add  first" << endl;
01613               addedstrip = true;
01614               firstv = 0;
01615             }
01616             else {
01617               // if we already have one or more u strips check that they are
01618               // in time, and require agreement in Z position.
01619               Bool_t found2=false;
01620               for (Int_t j=0; j<=neweventu.GetLast() && !found2;
01621                    j++) {
01622                 CandStripHandle *estrip =
01623                   dynamic_cast<CandStripHandle*>(neweventu.At(j));
01624                 Double_t dtime =
01625                   strip->GetBegTime()-estrip->GetBegTime();
01626                 MSG("EventSS",Msg::kVerbose) << " dtime " << dtime*1e9 
01627                                              << " " << pHitAssocDTime0*1e9 
01628                                              << "/" << pHitAssocDTime1*1e9 
01629                                              << " plane diff " 
01630                                              <<  abs(strip->GetPlane() - 
01631                                                      estrip->GetPlane()) 
01632                                              << "/" << pHitAssocDPlane << endl;
01633                 if (dtime>pHitAssocDTime0 &&
01634                     dtime<pHitAssocDTime1 &&
01635                     abs(strip->GetPlane()-estrip->GetPlane()) <=
01636                     pHitAssocDPlane) {
01637                   // if dist2 2 parameter meets match criteria, 
01638                   //add strip to array
01639                   MSG("EventSS",Msg::kVerbose) << " add  " << endl;
01640                   neweventv.Add(strip);
01641                   totcharge +=strip->GetCharge();
01642                   addedstrip = true;
01643                   firstv = false;
01644                   found2 = true;
01645                 }
01646               } // end loop over u strips in event
01647             } // end if already have u strips
01648           } // end if this is the first v strip
01649           else {
01650             // if we already have v strips, make sure that the one we
01651             // add is adjaccent, using same distance criteria base
01652             // on plane, time, and transverse postion separation
01653             Bool_t found2=false;
01654             for (Int_t j=0; j<=neweventv.GetLast() && !found2; j++) {
01655               CandStripHandle *estrip =
01656                 dynamic_cast<CandStripHandle*>(neweventv.At(j));
01657               Double_t dtime = strip->GetBegTime() -
01658                 estrip->GetBegTime();              
01659               MSG("EventSS",Msg::kVerbose)<< " dtime " << dtime*1e9 << " " 
01660                                           <<  pHitAssocDTime0*1e9 << "/" 
01661                                           << pHitAssocDTime1*1e9 
01662                                           << " plane diff " 
01663                                           <<  abs(strip->GetPlane() - 
01664                                                   estrip->GetPlane())
01665                                           << " dtpos " << (strip->GetTPos() - 
01666                                                            estrip->GetTPos())
01667                                           << endl;
01668               if (dtime>pHitAssocDTime0 &&
01669                   dtime<pHitAssocDTime1) {
01670                 Double_t dplane =
01671                   (Double_t)(strip->GetPlane()-estrip->GetPlane());
01672                 Double_t dtpos = strip->GetTPos()-estrip->GetTPos();
01673                 Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01674                   (pHitAssocPParm*dtpos*dtpos) +
01675                   (pHitAssocTParm*dtime*dtime);
01676                 // if dist2 2 parameter meets match criteria, 
01677                 // add strip to array
01678                 if (dist2<pHitAssocMaxDist2*4) {
01679                   MSG("EventSS",Msg::kVerbose) << " add  " << endl;
01680                   neweventv.Add(strip);
01681                   totcharge +=strip->GetCharge();
01682                   addedstrip = true;
01683                   firstv = false;
01684                   found2 = true;
01685                 } // end if add strip
01686               } //end if in time
01687             }  // end loop over v strips
01688           }
01689           break;
01690           
01691         default:
01692           break;
01693         }
01694       }
01695       
01696       // end of loop over hits - remove new 'associations' and compress
01697       for (Int_t i=0; i<=neweventu.GetLast(); i++) {
01698         CandStripHandle *estrip =
01699           dynamic_cast<CandStripHandle*>(neweventu.At(i));
01700         unassociated.Remove(estrip);
01701       }
01702       for (Int_t i=0; i<=neweventv.GetLast(); i++) {
01703         CandStripHandle *estrip =
01704           dynamic_cast<CandStripHandle*>(neweventv.At(i));
01705         unassociated.Remove(estrip);
01706       }
01707       unassociated.Compress();
01708     }
01709     
01710     // if we have hits in both views, build an event from them
01711     // if pulse height>10 pes
01712     MSG("EventSS",Msg::kDebug)<< " considering event formation  u/v: " 
01713                               << !firstu  << "/" << !firstv << " charge " 
01714                               << totcharge << endl;
01715     if (!firstu && !firstv && totcharge>20.0) {      
01716       found=true;
01717 
01718       //construct CandSubShowerSRs
01719       cxx.SetDataIn(&neweventu);
01720       AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default");
01721       CandSubShowerSRHandle usubshowerhandle =
01722         CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01723       usubshowerhandle.SetCandSlice(slice);
01724       usubshowerhandle.SetClusterID(ClusterType::kUnknown);
01725       subshowerlist->AddDaughterLink(usubshowerhandle);
01726       cxx.SetDataIn(&neweventv);
01727       CandSubShowerSRHandle vsubshowerhandle =
01728         CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01729       vsubshowerhandle.SetCandSlice(slice);
01730       vsubshowerhandle.SetClusterID(ClusterType::kUnknown);
01731       subshowerlist->AddDaughterLink(vsubshowerhandle);
01732       CandHandle *usubshower=subshowerlist->FindDaughter(&usubshowerhandle);
01733       CandHandle *vsubshower=subshowerlist->FindDaughter(&vsubshowerhandle);
01734       
01735       //build shower
01736       TObjArray newshower;
01737       newshower.Add(usubshower);
01738       newshower.Add(vsubshower);
01739       cxx.SetDataIn(&newshower);
01740       AlgHandle showerah = af.GetAlgHandle("AlgShowerSS","default");
01741       CandShowerSRHandle showerhandle =
01742         CandShowerSR::MakeCandidate(showerah,cxx);
01743       showerhandle.SetCandSlice(slice);
01744       showerhandle.SetCandRecord(slicelist->GetCandRecord());
01745 
01746       Bool_t newShowerOverlapsOld=false;
01747       CandShowerHandleItr showerItr(showerlist->GetDaughterIterator());      
01748       while (CandShowerHandle *shower =
01749              dynamic_cast<CandShowerHandle*>(showerItr())){
01750 
01751         Double_t dz = showerhandle.GetVtxZ()-shower->GetVtxZ();
01752         Double_t du = showerhandle.GetVtxU()-shower->GetVtxU();
01753         Double_t dv = showerhandle.GetVtxV()-shower->GetVtxV();
01754         Double_t dt = showerhandle.GetVtxT()-shower->GetVtxT();
01755 
01756         if( ( (du*du+dv*dv<pShwShwDtpos2 &&fabs(dz)<pShwShwDz) ||
01757               (du*du<pShwShwDtpos2/4. && dv*dv<pShwShwDtpos2*4. && 
01758                fabs(dz)<pShwShwDz/2.) ||
01759               (du*du<pShwShwDtpos2*4 && dv*dv<pShwShwDtpos2/4. && 
01760                fabs(dz)<pShwShwDz/2.) ) &&  
01761             fabs(dt)<pShwShwDt ) {
01762           newShowerOverlapsOld = true;
01763           MSG("EventSS",Msg::kDebug) 
01764             << " new shower overlaps with previous - don't make new event " 
01765             << endl;
01766           break;
01767         }
01768       }
01769         
01770       if(!newShowerOverlapsOld){
01771         if(!showerhandle.IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
01772           MSG("EventSS",Msg::kDebug)<< " Building Event! " << endl;
01773           showerlist->AddDaughterLink(showerhandle);
01774           TObjArray recolist;
01775           recolist.Add(&showerhandle);
01776           cxx.SetDataIn(&recolist);
01777           AlgHandle eventah = af.GetAlgHandle("AlgEventSR",pEventAlgConfig);
01778           CandEventHandle eventhandle =
01779             CandEvent::MakeCandidate(eventah,cxx);
01780           eventhandle.SetCandSlice(slice);
01781           eventhandle.SetPrimaryShower(pMinShwEFract,pShwShwDz);
01782           ch.AddDaughterLink(eventhandle);      
01783         }
01784       }
01785     }
01786   }
01787 }

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

Definition at line 1087 of file AlgEventSSList.cxx.

References CandHandle::AddDaughterLink(), CandEventHandle::AddShower(), CandHandle::FindDaughter(), Registry::Get(), AlgFactory::GetAlgHandle(), CandRecoHandle::GetBegPlane(), CandEventHandle::GetCandSlice(), CandHandle::GetDaughterIterator(), Registry::GetDouble(), AlgFactory::GetInstance(), CandSubShowerSRHandle::GetPlaneView(), CandStripHandle::GetPlaneView(), CandEventHandle::GetPrimaryTrack(), Msg::kError, PlaneView::kU, ClusterType::kUnknown, PlaneView::kV, Msg::kVerbose, CandShowerSR::MakeCandidate(), CandSubShowerSR::MakeCandidate(), MSG, CandRecoHandle::SetCandSlice(), CandSubShowerSRHandle::SetClusterID(), CandContext::SetDataIn(), and CandEventHandle::SetPrimaryShower().

Referenced by RunAlg().

01092 {
01093   // need to create a shower -- include strips from first plane in each
01094   // view plus this previously unassociated hit to generate a 3D object
01095   MSG("EventSS",Msg::kVerbose) << "In CreatePrimaryShower" << endl;
01096   // Get singleton instance of AlgFactory
01097   const char *pEventAlgConfig = 0;
01098   ac.Get("EventAlgConfig",pEventAlgConfig);
01099 
01100   //Double_t shwshwdz = ac.GetDouble("ShwShwDz");
01101   //Double_t minShwEFract = ac.GetDouble("MinShwEFract");
01102   Double_t pMinShwStripPE = ac.GetDouble("MinShwStripPE");
01103   Double_t pMaxNewShwLen = ac.GetDouble("MaxNewShwLen");
01104 
01105   AlgFactory &af = AlgFactory::GetInstance();
01106   
01107   CandTrackHandle *track = closestevent->GetPrimaryTrack();
01108   if (!showerlist) {
01109     MSG("EventSS",Msg::kError) <<"CandShowerListHandle missing\n" << endl;
01110     return;
01111   }
01112   
01113   closestevent->AddDaughterLink(*closeststrip);
01114   
01115   // create arrays of strips for each planeview to create
01116   // shower from
01117   
01118   MSG("EventSS",Msg::kVerbose) <<  "creating shower \n";
01119   TObjArray ustriparray;
01120   TObjArray vstriparray;
01121   Int_t uplane=track->GetBegPlane(PlaneView::kU);
01122   Int_t vplane=track->GetBegPlane(PlaneView::kV);
01123   // seed this shower with strips in most upstream track hits 
01124   // in each view, along with hits above a charge threshold and
01125   // within pMaxNewShwLen of  the track start
01126   TIter stripItr(track->GetDaughterIterator());
01127   while (CandStripHandle *tstrip =
01128          dynamic_cast<CandStripHandle*>(stripItr())) {
01129     if ((tstrip->GetPlane() == uplane) || 
01130         (tstrip->GetPlaneView()==PlaneView::kU &&
01131             tstrip->GetCharge()>pMinShwStripPE && 
01132             abs(tstrip->GetPlane()-uplane)<pMaxNewShwLen)) {
01133       ustriparray.Add(tstrip);
01134     }
01135     if ((tstrip->GetPlane() == vplane) || 
01136         (tstrip->GetPlaneView()==PlaneView::kV &&
01137             tstrip->GetCharge()>pMinShwStripPE && 
01138             abs(tstrip->GetPlane()-vplane)<pMaxNewShwLen)) {
01139       vstriparray.Add(tstrip);
01140     }
01141   }
01142   // add unassoc strip to appropriate array
01143   switch (closeststrip->GetPlaneView()) {
01144   case PlaneView::kU:
01145     ustriparray.Add(closeststrip);
01146     break;
01147   case PlaneView::kV:
01148     vstriparray.Add(closeststrip);
01149     break;
01150   default:
01151     break;
01152   }
01153   // build shower from U/V clusters or subshowers
01154   if(ustriparray.GetLast()>=0 && vstriparray.GetLast()>=0){
01155     if(showerlist->InheritsFrom("CandShowerSRListHandle") && subshowerlist){
01156       MSG("EventSS",Msg::kVerbose) << "U striparray has " 
01157                                    << ustriparray.GetLast()+1 << " entries" << endl;
01158       MSG("EventSS",Msg::kVerbose) << "V striparray has " 
01159                                    << vstriparray.GetLast()+1 << " entries" << endl;  
01160 
01161       // now construct CandSubShowers
01162       AlgHandle subshowerah = af.GetAlgHandle("AlgSubShowerSR","default");      
01163       cxx.SetDataIn(&ustriparray);
01164       CandSubShowerSRHandle usubshowerhandle =
01165         CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01166       usubshowerhandle.SetCandSlice(closestevent->GetCandSlice());
01167       usubshowerhandle.SetClusterID(ClusterType::kUnknown);
01168       subshowerlist->AddDaughterLink(usubshowerhandle);
01169       MSG("EventSS",Msg::kVerbose) << "U subshower view: " 
01170                                    << usubshowerhandle.GetPlaneView() <<endl;
01171       
01172       cxx.SetDataIn(&vstriparray);
01173       CandSubShowerSRHandle vsubshowerhandle =
01174         CandSubShowerSR::MakeCandidate(subshowerah,cxx);
01175       vsubshowerhandle.SetCandSlice(closestevent->GetCandSlice());
01176       vsubshowerhandle.SetClusterID(ClusterType::kUnknown);
01177       subshowerlist->AddDaughterLink(vsubshowerhandle);
01178       MSG("EventSS",Msg::kVerbose) << "V subshower view: " 
01179                                    << vsubshowerhandle.GetPlaneView() <<endl;
01180       
01181       //add subshowers to list
01182       CandShowerSRListHandle* showerlistSR = 
01183         dynamic_cast<CandShowerSRListHandle*>(showerlist);
01184       TObjArray newshower;
01185       CandSubShowerSRHandle * usubshower = dynamic_cast<CandSubShowerSRHandle*>
01186         (subshowerlist->FindDaughter(&usubshowerhandle));
01187       
01188       CandSubShowerSRHandle * vsubshower = dynamic_cast<CandSubShowerSRHandle*>
01189         (subshowerlist->FindDaughter(&vsubshowerhandle));           
01190       newshower.Add(usubshower);
01191       newshower.Add(vsubshower);
01192       cxx.SetDataIn(&newshower);
01193       
01194       //get ShowerSS alghandle and make shower from subshowers
01195       AlgHandle showerahSS = af.GetAlgHandle("AlgShowerSS","default");
01196       CandShowerSRHandle showerhandleSR = CandShowerSR::MakeCandidate(showerahSS,cxx);
01197       showerhandleSR.SetCandSlice(closestevent->GetCandSlice());
01198       showerlistSR->AddDaughterLink(showerhandleSR);
01199       CandHandle *shw = showerlistSR->FindDaughter(&showerhandleSR);
01200       CandShowerHandle *shwh = dynamic_cast<CandShowerHandle*>(shw);
01201       closestevent->AddShower(shwh);
01202       // Change: cbs 9Aug06
01203       // NEED TO FORCE CREATED SHOWER TO BE PRIMARY
01204       // due to a change in assignment of primary shower, can now create
01205       // a shower here which is *not* called primary by SetPrimaryShower()
01206       // Causes multiple CreatePrimaryShower calls in event builder and 
01207       // large numbers of reco'd showers per events. 
01208       // SetPrimaryShowers is called at the end of RunAlg so a correct
01209       // assignment will be made later
01210       //closestevent->SetPrimaryShower(minShwEFract,shwshwdz);
01211       closestevent->SetPrimaryShower(shwh);
01212     }
01213   }
01214 }

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

Referenced by RunAlg().

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

Definition at line 1791 of file AlgEventSSList.cxx.

References CandHandle::IsCloneOf(), Msg::kDebug, Msg::kError, and MSG.

Referenced by RunAlg().

01792                                                                                      {
01793 
01794   if (showerItr) {
01795     showerItr->Reset();
01796     while (CandShowerHandle *shower =
01797            dynamic_cast<CandShowerHandle*>((*showerItr)())) {
01798       if (shower->GetCandSlice()) {
01799         if (slice->IsCloneOf(*(shower->GetCandSlice()))) {
01800           recolist.Add(shower);
01801         }
01802       }
01803       else {
01804         MSG("EventSS", Msg::kError)
01805           << "Shower doesn't contain a slice.  Not added to recolist."
01806           << endl;
01807       }
01808     }
01809   }
01810   
01811   if (trackItr) {
01812     trackItr->Reset();
01813     while (CandTrackHandle *track =
01814            dynamic_cast<CandTrackHandle*>((*trackItr)())) {
01815       if (track->GetCandSlice()) {
01816         if (slice->IsCloneOf(*track->GetCandSlice())) {
01817           if(track->GetNDaughters()==0 && 
01818              track->InheritsFrom("CandFitTrackHandle") &&  
01819              dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack()){
01820             track = dynamic_cast<CandFitTrackHandle*>(track)->GetFinderTrack();    
01821             MSG("EventSS", Msg::kDebug)<< "Finder Track being used " << endl;
01822           }
01823           recolist.Add(track);
01824         }
01825       }
01826       else {
01827         MSG("EventSS", Msg::kError)
01828           << "Track doesn't contain a slice.  Not added to recolist."
01829           << endl;
01830       }
01831     }
01832   }
01833 }

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

Definition at line 1858 of file AlgEventSSList.cxx.

References Registry::GetDouble(), and CalDigitType::kPE.

Referenced by RunAlg().

01860                                                                {
01861 
01862   Double_t pMinUnassocStripPE = ac.GetDouble("MinUnassocStripPE");
01863 
01864   while (CandSliceHandle *slice =
01865          dynamic_cast<CandSliceHandle*>(sliceItr())) {
01866     TIter stripItr(slice->GetDaughterIterator());
01867     while (CandStripHandle *strip =
01868            dynamic_cast<CandStripHandle*>(stripItr())) {
01869 
01870       // Ignore strips below the required minimum PE
01871       if (strip->GetCharge(CalDigitType::kPE)<pMinUnassocStripPE) continue;
01872 
01873       Bool_t found=false;
01874       // check for duplicates      
01875       for (Int_t i=0; i<=unassociated.GetLast() && !found; i++) {
01876         CandStripHandle *strip0 =
01877           dynamic_cast<CandStripHandle*>(unassociated.At(i));
01878         if (*strip == *strip0) {
01879           found = true;
01880           break;
01881         }
01882       }
01883       // if this strip not already in unnassociated list, then
01884       // loop over objects in this slice, and check against all daughter strips
01885       if (!found) {
01886         if (eventItr) {
01887           eventItr->Reset();
01888           while (CandEventHandle *event =
01889                  dynamic_cast<CandEventHandle*>((*eventItr)())) {
01890             if (event->GetCandSlice()) {
01891               if (slice->IsCloneOf(*(event->GetCandSlice()))) { 
01892                 if (event->FindDaughter(strip)) {
01893                   found = true;
01894                   break;
01895                 }
01896               }
01897             }
01898           }
01899         }
01900       }
01901       // this strip is not in any event, add to unassociated list
01902       if (!found) {
01903         unassociated.Add(strip);
01904       }
01905     }
01906   }
01907 }

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(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), AlgFactory::GetInstance(), CandHandle::GetUidInt(), CandHandle::IsCloneOf(), CandHandle::IsEqual(), Msg::kDebug, Msg::kError, Msg::kVerbose, Msg::kWarning, MSG, CandHandle::RemoveDaughter(), AlgHandle::RunAlg(), CandRecoHandle::SetCandSlice(), and CandContext::SetDataIn().

Referenced by RunAlg().

01248 {
01249   MSG("EventSS",Msg::kVerbose) << "In ReConstructShowers: " << endl;
01250   // Retrieve persistant modifiable CandShowerList
01251   CandShowerListHandle *showerlist = const_cast<CandShowerListHandle*>
01252     (dynamic_cast<const CandShowerListHandle*>(cxx.GetDataIn()));
01253 
01254   if(!showerlist->InheritsFrom("CandShowerSRListHandle")) {
01255     MSG("EventSS",Msg::kWarning)
01256       << "ListHandle is not a CandShowerSRListHandle."
01257       << "  ShowerList not processed in ReConstructShowers." << endl;
01258     return;
01259   }
01260 
01261   // Get singleton instance of AlgFactory
01262   AlgFactory &af = AlgFactory::GetInstance();
01263   
01264   // Get an AlgHandle to AlgShowerSR with default AlgConfig
01265   const char *pEventAlgConfig = 0;
01266   ac.Get("EventAlgConfig",pEventAlgConfig);
01267 
01268   AlgHandle ahss = af.GetAlgHandle("AlgShowerSS","default");
01269   AlgConfig acss = ahss.GetAlgConfig();
01270   Int_t eventCounter = 0;
01271   TIter evItr(ch.GetDaughterIterator());
01272   while (CandEventHandle *ev =
01273          dynamic_cast<CandEventHandle*>(evItr())) {  //event loop
01274     MSG("EventSS",Msg::kVerbose) << "event = " << eventCounter << endl;
01275     for (Int_t ishower=0; ishower<ev->GetLastShower()+1; ishower++) { //shower loop
01276       MSG("EventSS",Msg::kVerbose) << "shower = " << ishower << endl;
01277       CandShowerHandle *shower =
01278         dynamic_cast<CandShowerHandle*>(ev->GetShowerWritable(ishower));
01279       CandShowerSRHandle *showerSR = 0;
01280       if(shower->InheritsFrom("CandShowerSRHandle")) {
01281         showerSR = dynamic_cast<CandShowerSRHandle*>(shower);
01282       }
01283       MSG("EventSS",Msg::kVerbose) << "Shower UID = " << showerSR->GetUidInt() << endl;
01284       //assume subshower chain
01285       if(showerSR) {  
01286         TObjArray newshower;
01287         // fill array of subshowers for this shower
01288         for (Int_t isubshower=0;
01289              isubshower<showerSR->GetLastSubShower()+1;
01290              isubshower++) {
01291           MSG("EventSS",Msg::kVerbose) << "subshower = " << isubshower << endl;
01292           const CandSubShowerSRHandle *ssh = 
01293             showerSR->GetSubShower(isubshower);
01294           if (ssh) {
01295             // New handle on heap
01296             CandSubShowerSRHandle *sshd = ssh->DupHandle();
01297             // Add new heap-based CandHandle to TObjArray
01298             newshower.Add(sshd);
01299             MSG("EventSS",Msg::kDebug) << "Uid's (old,new) " 
01300                                        << ssh->GetUidInt() << " " 
01301                                        << sshd->GetUidInt() << endl;
01302           }
01303         }
01304         // re-run algorithm for this shower
01305         if (newshower.GetEntriesFast() > 0) {
01306           cxx.SetDataIn(&newshower);
01307           ahss.RunAlg(*showerSR,cxx);
01308         }
01309         else {
01310           MSG("EventSS", Msg::kError)
01311             << "Attempt to pass empty TObjArray newshower to AlgShowerSS"
01312             << " in CandContext.  AlgShowerSS::RunAlg() call ignored."
01313             << endl;
01314         }
01315         // Delete heap-based CandHandles in TObjArray
01316         newshower.Delete();
01317       }
01318       else {
01319         MSG("EventSS",Msg::kWarning)
01320           << "Handle is not a CandShowerSRHandle."
01321           << "  Shower not processed in ReConstructShowers."
01322           << endl;
01323         continue;
01324       }
01325       MSG("EventSS",Msg::kVerbose) << "Modified shower UID = " 
01326                                    << showerSR->GetUidInt() << endl; 
01327       // now loop through strips in this shower, and determine whether primary
01328       // track in this event intercepts this strip.  If so, subtract 1 mip
01329       // from shower energy.
01330       
01331       // if event has primary track, extrapolate through shower 
01332       // and subtract energy loss to recalibrate 
01333       CandTrackHandle *primaryTrack = ev->GetPrimaryTrack();
01334       if (primaryTrack) {
01335         shower->CalibrateEnergy(primaryTrack,acss);
01336       }
01337       
01338       // Update modified showers in persistent modifiable CandShowerList.
01339       // CandHandle could potentially supply a ReplaceDaughter() function.
01340       // Interface of ReplaceDaughter() would be debatable.
01341       if (showerlist) {
01342         TIter shiter(showerlist->GetDaughterIterator());
01343         CandShowerHandle *target;
01344         Bool_t found = kFALSE;
01345         while (!found &&
01346                (target = dynamic_cast<CandShowerHandle*>(shiter()))) {
01347           if (shower->IsCloneOf(*target)) {  // Tests clone or Cand ==
01348             found = kTRUE;
01349             
01350             // Replace target only if shower is a modified version of target
01351             if (!(shower->IsEqual(target))) {   // CandBase inequality
01352               shower->SetCandSlice(ev->GetCandSlice());
01353               if (!showerlist->RemoveDaughter(target)) {    // Failure
01354                 MSG("EventSS",Msg::kError) 
01355                   << endl
01356                   << "Failure of ShowerList daughter removal " << endl
01357                   << "attempted during replacement by modified Shower."
01358                   << endl << "Will result in double counted Shower."
01359                   << endl;
01360               }
01361               showerlist->AddDaughterLink(*shower);
01362             }
01363           }
01364         }
01365         
01366         // Add new (unfound) shower to persistent modifiable CandShowerList.
01367         if (!found){
01368           shower->SetCandSlice(ev->GetCandSlice());
01369           showerlist->AddDaughterLink(*shower);
01370         }
01371       }
01372     }
01373     eventCounter++;
01374   }
01375 }

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 Registry::GetDouble(), CandStripHandle::GetPlane(), CandStripHandle::GetPlaneView(), CandStripHandle::GetTPos(), and CandRecoHandle::GetVtxT().

Referenced by RunAlg().

01218 {
01219   Double_t pHitAssocZParm = ac.GetDouble("HitAssocZParm");
01220   Double_t pHitAssocPParm = ac.GetDouble("HitAssocPParm");
01221   Double_t pHitAssocTParm = ac.GetDouble("HitAssocTParm");
01222 
01223   for (Int_t i=0; i<=unassociated.GetLast(); i++) {
01224     CandStripHandle *strip =
01225       dynamic_cast<CandStripHandle*>(unassociated.At(i));
01226     if(strip->GetPlaneView()==closeststrip->GetPlaneView()){
01227       Double_t dtime = strip->GetBegTime()-closestevent->GetVtxT();
01228       Double_t dplane =
01229         (Double_t)(strip->GetPlane()-closeststrip->GetPlane());
01230       Double_t dtpos = strip->GetTPos()-closeststrip->GetTPos();
01231       Double_t dist2 = (pHitAssocZParm*dplane*dplane) +
01232         (pHitAssocPParm*dtpos*dtpos) +
01233         (pHitAssocTParm*dtime*dtime);
01234       if (dist2<dist2map[i][closestevent]) {
01235         dist2map[i][closestevent] = dist2;
01236       }
01237     }
01238   } // loop over unassoc hits
01239 }

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

Implements AlgBase.

Definition at line 72 of file AlgEventSSList.cxx.

References CandHandle::AddDaughterLink(), AddStripToEvent(), CandTrackHandle::BelongsWithTrack(), BuildEventFromUnassoc(), CandShowerSRHandle::BuriedTrack(), CreatePrimaryShower(), FillDist2Map(), FillRecoList(), FindUnassociated(), Registry::Get(), CandRecoHandle::GetBegPlane(), CandStripHandle::GetBegTime(), CandContext::GetCandRecord(), CandContext::GetDataIn(), CandHandle::GetDaughterIterator(), CandRecoHandle::GetDirCosU(), CandRecoHandle::GetDirCosV(), Registry::GetDouble(), CandRecoHandle::GetEndPlane(), AlgFactory::GetInstance(), Registry::GetInt(), CandEventHandle::GetLastShower(), CandContext::GetMom(), CandHandle::GetNDaughters(), CandStripHandle::GetPlane(), CandEventHandle::GetPrimaryShower(), CandEventHandle::GetPrimaryTrack(), CandEventHandle::GetShower(), CandEventHandle::GetShowerWritable(), CandStripHandle::GetTPos(), RecMinos::GetVldContext(), CandRecoHandle::GetVtxPlane(), CandShowerSRHandle::IsUnphysical(), Msg::kDebug, Msg::kError, CandSubShowerSRHandle::KeyFromSlice(), CandTrackHandle::KeyFromSlice(), CandShowerHandle::KeyFromSlice(), CalStripType::kMIP, PlaneView::kU, PlaneView::kV, Msg::kVerbose, Msg::kWarning, CandEvent::MakeCandidate(), max, min, MSG, PITCHINMETRES, 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 //gmi Change by G. Irwin 28Feb2009 -
00237 //gmi For GCC 4.4, change old statement below to model (following that) of
00238 //gmi http://www.codepedia.com/1/Cpp2dVector (see "One-line construction")
00239 //gmi Original line (doesn't compile in GCC 4.4):
00240 //gmi std::vector<std::vector<Int_t> > objectMatch(recolist.GetLast()+1,0);
00241 //gmi Model for replacement line:
00242 //gmi std::vector<std::vector<int> > vec(maxx,std::vector<int>(maxy,0));
00243 //gmi Replacement line itself:
00244     std::vector<std::vector<Int_t> > objectMatch(recolist.GetLast()+1,
00245                 std::vector<Int_t>(recolist.GetLast()+1,0));
00246     //initialize:
00247     for(int i=0;i<=recolist.GetLast();i++){
00248       for(int j=0;j<=recolist.GetLast();j++) {
00249 //gmi Replace this line:  "objectMatch[i].push_back(0);"  by:
00250         objectMatch[i][j] = 0;  //gmi May be redundant with constructor init
00251       }
00252     }
00253 
00254     //use this vector to keep track of objects already assigned
00255     std::vector<Bool_t> objectUsed(recolist.GetLast()+1,0);
00256     
00257     //loop over all objects
00258     for(Int_t i=0;i<=recolist.GetLast();i++) { //loop over i
00259       objectUsed[i] = false;  //initialization
00260       
00261       CandShowerSRHandle *shower=0;
00262       CandTrackHandle *track=0;
00263       if (recolist.At(i)->InheritsFrom("CandShowerSRHandle"))
00264         shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00265       if (recolist.At(i)->InheritsFrom("CandTrackHandle"))
00266         track = dynamic_cast<CandTrackHandle*>(recolist.At(i));
00267       
00268       //loop over all later objects only (since matching matrix is symmetric)
00269       for(Int_t j=i+1;j<=recolist.GetLast();j++) { //loop over j
00270 
00271         MSG("EventSS",Msg::kDebug) << "Checking objectMatch between " 
00272                                    << i << " and " << j << endl;
00273 
00274         CandShowerSRHandle *shower2=0;
00275         CandTrackHandle *track2=0;
00276         if (recolist.At(j)->InheritsFrom("CandShowerSRHandle"))
00277           shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(j));
00278         if (recolist.At(j)->InheritsFrom("CandTrackHandle"))
00279           track2 = dynamic_cast<CandTrackHandle*>(recolist.At(j));
00280         
00281         if(shower) {
00282           if(shower2) {
00283             if(!shower->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef) && 
00284                !shower2->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00285               if(shower->BelongsWithShower(shower2,ac,vldcptr,pShwShwDtpos2,
00286                                            pShwShwDz,pShwShwDt)) objectMatch[i][j] = 1;
00287             }
00288             objectMatch[j][i] = objectMatch[i][j];
00289           }
00290           else if(track2) {
00291             if(shower->BelongsWithTrack(track2,ac,vldcptr,pShwTrkDtpos2,
00292                                         pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1;
00293 
00294             if(objectMatch[i][j]!=1) {
00295               //check track isn't buried:
00296               Int_t *stats = new Int_t[8];
00297               shower->BuriedTrack(track2,1000,stats); //1000 forces all calculations
00298               if(stats[3]>0) { //at least one buried track planes
00299                 Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00300                 MSG("EventSS",Msg::kDebug) 
00301                   << "BuriedTrack Stats: " << "\n"
00302                   << "Number of track planes = " << nplanes << "\n"
00303                   << "nSharedPlanes = " << stats[0] << "\n" 
00304                   << "nMissingIsoPlanes = " << stats[1] << "\n" 
00305                   << "nMissingSharedPlanes = " << stats[2] << "\n" 
00306                   << "nBuriedPlanes = " << stats[3] << "\n" 
00307                   << "nPhysBuriedPlanes = " << stats[4] << "\n" 
00308                   << "nTrkLike = " << stats[5] << "\n" 
00309                   << "nTrkLikeTopol = " << stats[6] << "\n" 
00310                   << "nXTalk = " << stats[7] << endl;
00311                 nplanes -= (stats[1] + stats[2]); //subtract missing track planes
00312                 //if more than half of the track planes are in the shower core: match
00313                 if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac)
00314                   objectMatch[i][j] = -1;
00315                 //else if more than half the track planes are anywhere in 
00316                 //the shower or are x-talk like and are not "trklike": match
00317                 else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > 
00318                         pHardBuriedFrac) objectMatch[i][j] = -1;
00319                 //check whether track passes through shower at high angle in either view
00320                 //expect significant number of buried planes
00321                 else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && 
00322                         ( TMath::Abs(track2->GetDirCosU())>pDCosUVCut || 
00323                           TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) {
00324                   objectMatch[i][j] = -1;
00325                 }
00326               }
00327               delete [] stats;
00328             }
00329             objectMatch[j][i] = objectMatch[i][j];
00330           }
00331         }
00332         
00333         if(track) {
00334           if(shower2) {
00335             if(shower2->BelongsWithTrack(track,ac,vldcptr,pShwTrkDtpos2,
00336                                          pShwTrkDz,pShwTrkDt)) objectMatch[i][j] = 1;
00337             if(objectMatch[i][j]!=1) {
00338               //check track isn't buried:
00339               Int_t *stats = new Int_t[8];
00340               shower2->BuriedTrack(track,1000,stats); //1000 forces all calculations
00341               if(stats[3]>0) { //at least one buried track planes
00342                 Int_t nplanes = track->GetEndPlane() - track->GetBegPlane() + 1;
00343                 MSG("EventSS",Msg::kDebug) 
00344                   << "BuriedTrack Stats: " << "\n"
00345                   << "Number of track planes = " << nplanes << "\n"
00346                   << "nSharedPlanes = " << stats[0] << "\n" 
00347                   << "nMissingIsoPlanes = " << stats[1] << "\n" 
00348                   << "nMissingSharedPlanes = " << stats[2] << "\n" 
00349                   << "nBuriedPlanes = " << stats[3] << "\n" 
00350                   << "nPhysBuriedPlanes = " << stats[4] << "\n" 
00351                   << "nTrkLike = " << stats[5] << "\n" 
00352                   << "nTrkLikeTopol = " << stats[6] << "\n" 
00353                   << "nXTalk = " << stats[7] << endl;
00354                 nplanes -= (stats[1] + stats[2]); //subtract missing track planes
00355                 //if more than half of the track planes are in the shower core: match
00356                 if(Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac)
00357                   objectMatch[i][j] = -1;
00358                 //else if more than half the track planes are anywhere in 
00359                 //the shower or are x-talk like and are not "trklike": match
00360                 else if(Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > 
00361                         pHardBuriedFrac) objectMatch[i][j] = -1;
00362                 //check whether track passes through shower at high angle in either view
00363                 //expect significant number of buried planes
00364                 else if(stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && 
00365                         (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || 
00366                          TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) {
00367                   objectMatch[i][j] = -1;
00368                 }
00369               }
00370               delete [] stats;
00371             }
00372             objectMatch[j][i] = objectMatch[i][j];
00373           }
00374           else if(track2) {
00375             if(track->BelongsWithTrack(track2,ac,vldcptr,pTrkTrkDtpos2,
00376                                        pTrkTrkDz,pTrkTrkDt)) objectMatch[i][j] = 1;
00377             objectMatch[j][i] = objectMatch[i][j];
00378           }
00379         }
00380         MSG("EventSS",Msg::kDebug) << "Match for objects " << i << " and " 
00381                                    << j << " : " << objectMatch[i][j] << endl;
00382       }
00383     }
00384 
00385     MSG("EventSS",Msg::kDebug) << "Starting to form events" << endl;
00386 
00387     //start to make events:
00388     Int_t nRecoObjects = recolist.GetLast()+1;
00389     //get total charge from each object
00390     std::vector<Float_t> ph(nRecoObjects,0);
00391     for(int i=0;i<nRecoObjects;i++){
00392       if (recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {
00393         CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00394         ph[i] = shower->GetCharge(CalStripType::kMIP);
00395       }
00396       if (recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00397         CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(i));
00398         //ph[i] = track->GetCharge(CalStripType::kMIP)+1e7; //force tracks to come first
00399         //force tracks to come first, order tracks by number of strips
00400         ph[i] = track->GetNStrip()+1e7; 
00401       }
00402     }
00403 
00404     Int_t nEvents = 0;
00405     Int_t nUnusedObjects = nRecoObjects;
00406     while(nUnusedObjects>0) { //loop until all reco objects are in events
00407       
00408       //find largest object:
00409       Float_t maxPH = 0;
00410       Int_t maxPHIndex = -1;
00411       for(int i=0;i<nRecoObjects;i++){
00412         if(ph[i]>maxPH) {
00413           maxPH = ph[i];
00414           maxPHIndex = i;
00415         }
00416       }
00417       
00418       //add largest objects first:
00419       TObjArray *event = new TObjArray;
00420       event->Add(recolist.At(maxPHIndex));     
00421       ph[maxPHIndex] = 0;
00422       objectUsed[maxPHIndex] = true;
00423       nUnusedObjects-=1;
00424 
00425       MSG("EventSS",Msg::kDebug) << "Added object " << maxPHIndex 
00426                                  << " to event " << nEvents << endl;
00427       
00428       if(recolist.At(maxPHIndex)->InheritsFrom("CandTrackHandle")) {
00429         MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a track" << endl;
00430 
00431         //maxPH object is a track
00432         CandTrackHandle *track = dynamic_cast<CandTrackHandle*>(recolist.At(maxPHIndex));
00433         
00434         //now add objects clearly associated with largest object:
00435         for(int i=0;i<nRecoObjects;i++){
00436           if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) {       
00437             if(recolist.At(i)->InheritsFrom("CandTrackHandle")) { 
00438               MSG("EventSS",Msg::kDebug) << "Track match. Adding object " << i 
00439                                          << " to event " << nEvents << endl;
00440               //if this match is a track, add it to event
00441               event->Add(recolist.At(i));
00442               ph[i] = 0;
00443               objectUsed[i] = true;
00444               nUnusedObjects-=1;
00445             } // end if for CandTrackHandle
00446             else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {           
00447               MSG("EventSS",Msg::kDebug) << "Shower match. Testing whether to add object " 
00448                                          << i << " to event " << nEvents << endl;
00449               CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00450               //if this match is a shower, check where shower is relative to track vertex
00451               //if(TMath::Abs(shower2->GetVtxZ()-track->GetVtxZ())<pShwTrkDz/2.){
00452               if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) - 
00453                             track->GetBegPlane(PlaneView::kU)) < 
00454                  pShwTrkDz/(PITCHINMETRES*2.) || 
00455                  TMath::Abs(shower2->GetBegPlane(PlaneView::kV) - 
00456                             track->GetBegPlane(PlaneView::kV)) < 
00457                  pShwTrkDz/(PITCHINMETRES*2.)){
00458                 MSG("EventSS",Msg::kDebug) << "Shower " << i << " is close to track vertex. "
00459                                            << "Adding to event " << nEvents << endl;
00460                 event->Add(recolist.At(i));
00461                 ph[i] = 0;
00462                 objectUsed[i] = true;
00463                 nUnusedObjects-=1;
00464                 //look for objects associated with this vertex shower
00465                 for(int j=0;j<nRecoObjects;j++){
00466                   if(j!=maxPHIndex && objectMatch[i][j]!=0 && !objectUsed[j]) { //got a match
00467                     if(objectMatch[maxPHIndex][j]==0) { // j does not match maxPHIndex
00468                       if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00469                         MSG("EventSS",Msg::kDebug) << "Shower " << i 
00470                                                    << " has another track match. " 
00471                                                    << "Test how buried track " << j
00472                                                    << " is in shower" << endl; 
00473                         if(objectMatch[i][j]==-1){
00474                           event->Add(recolist.At(j));
00475                           ph[j] = 0;
00476                           objectUsed[j] = true;
00477                           nUnusedObjects-=1;
00478                         }
00479                         else {
00480                           CandTrackHandle *track2 = 
00481                             dynamic_cast<CandTrackHandle*>(recolist.At(j));
00482                           //check if track2 is buried
00483                           Int_t *stats = new Int_t[8];
00484                           shower2->BuriedTrack(track2,1000,stats);
00485                           if(stats[3]>0){
00486                             Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00487                             MSG("EventSS",Msg::kDebug) 
00488                               << "BuriedTrack Stats: " << "\n"
00489                               << "Number of track planes = " << nplanes << "\n"
00490                               << "nSharedPlanes = " << stats[0] << "\n" 
00491                               << "nMissingIsoPlanes = " << stats[1] << "\n" 
00492                               << "nMissingSharedPlanes = " << stats[2] << "\n" 
00493                               << "nBuriedPlanes = " << stats[3] << "\n" 
00494                               << "nPhysBuriedPlanes = " << stats[4] << "\n" 
00495                               << "nTrkLike = " << stats[5] << "\n" 
00496                               << "nTrkLikeTopol = " << stats[6] << "\n" 
00497                               << "nXTalk = " << stats[7] << endl;
00498                             nplanes -= (stats[1] + stats[2]);
00499                             if( (Float_t(stats[4])/Float_t(nplanes) > pHardBuriedFrac) ||
00500                                 (Float_t(stats[3]-stats[5]-stats[6])/Float_t(nplanes) > 
00501                                  pHardBuriedFrac) ||
00502                                 (stats[3]>=Int_t(nplanes*pSoftBuriedFrac) && 
00503                                  (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || 
00504                                   TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) {
00505                               MSG("EventSS",Msg::kDebug) << "Adding track " << j 
00506                                                          << " to event " << nEvents << endl;
00507                               event->Add(recolist.At(j));
00508                               ph[j] = 0;
00509                               objectUsed[j] = true;
00510                               nUnusedObjects-=1;
00511                             }
00512                           }
00513                           delete [] stats;
00514                         }
00515                       }
00516                       else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00517                         MSG("EventSS",Msg::kDebug) << "Shower " << i 
00518                                                    << " has another shower match. " 
00519                                                    << "Test whether shower " << j 
00520                                                    << " has a better association" << endl; 
00521                         Bool_t isAssocTrack = false; //is this shower assoc with a track
00522                         for(int k=0;k<nRecoObjects;k++) { 
00523                           //check for match
00524                           if(k!=maxPHIndex && k!=i && objectMatch[j][k]!=0){
00525                             //check for track
00526                             if(recolist.At(k)->InheritsFrom("CandTrackHandle"))
00527                               //check it does not match with shower2
00528                               if(objectMatch[i][k]==0) isAssocTrack = true;
00529                           }
00530                         }
00531                         if(!isAssocTrack) { //no track associated
00532                           MSG("EventSS",Msg::kDebug) << "Shower " << j 
00533                                                      << " has no track association. " 
00534                                                      << "Add to event " << nEvents << endl; 
00535                           event->Add(recolist.At(j));
00536                           ph[j] = 0;
00537                           objectUsed[j] = true;
00538                           nUnusedObjects-=1;
00539                         }
00540                       }
00541                     }
00542                     else { //j matches maxPHIndex, so add it to event
00543                       MSG("EventSS",Msg::kDebug) << "Shower " << j 
00544                                                  << " matches orginal track " << maxPHIndex 
00545                                                  << " Add to event " << nEvents << endl; 
00546                       event->Add(recolist.At(j));
00547                       ph[j] = 0;
00548                       objectUsed[j] = true;
00549                       nUnusedObjects-=1;
00550                     }
00551                   }
00552                 }
00553               } //end if for vertex shower
00554               else { //if not a vertex shower, check other associations
00555                 MSG("EventSS",Msg::kDebug) << "Testing whether shower " << i 
00556                                            << " has better associations. " << endl;
00557                 Bool_t otherTrackVertexMatch = false;
00558                 Bool_t otherTrackMatch = false;
00559                 Bool_t otherShowerMatch = false;
00560                 for(int j=0;j<nRecoObjects;j++){
00561                   if(j!=maxPHIndex && objectMatch[i][j]!=0 && 
00562                      !objectUsed[j] && objectMatch[maxPHIndex][j]==0) {
00563                     if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00564                       CandTrackHandle *track2 = 
00565                         dynamic_cast<CandTrackHandle*>(recolist.At(j));
00566                       if(TMath::Abs(shower2->GetBegPlane(PlaneView::kU) - 
00567                                     track->GetBegPlane(PlaneView::kU)) < 
00568                          pShwTrkDz/(PITCHINMETRES*2.) || 
00569                          TMath::Abs(shower2->GetBegPlane(PlaneView::kV) -
00570                                     track->GetBegPlane(PlaneView::kV)) < 
00571                          pShwTrkDz/(PITCHINMETRES*2.) )
00572                         otherTrackVertexMatch = true;
00573                       else {
00574                         //first test whether track2 is buried
00575                         //if so then can still associate shower to track
00576                         if(objectMatch[i][j]!=-1) {
00577                           //object match is not flagged as buried
00578                           //perform check:
00579                           Bool_t track2Buried = false;
00580                           Int_t *stats2 = new Int_t[8];
00581                           shower2->BuriedTrack(track2,1000,stats2);                       
00582                           if(stats2[3]>0){
00583                             Int_t nplanes = track2->GetEndPlane() - track2->GetBegPlane() + 1;
00584                             MSG("EventSS",Msg::kDebug) 
00585                               << "BuriedTrack Stats: " << "\n"
00586                               << "Number of track planes = " << nplanes << "\n"
00587                               << "nSharedPlanes = " << stats2[0] << "\n" 
00588                               << "nMissingIsoPlanes = " << stats2[1] << "\n" 
00589                               << "nMissingSharedPlanes = " << stats2[2] << "\n" 
00590                               << "nBuriedPlanes = " << stats2[3] << "\n" 
00591                               << "nPhysBuriedPlanes = " << stats2[4] << "\n" 
00592                               << "nTrkLike = " << stats2[5] << "\n" 
00593                               << "nTrkLikeTopol = " << stats2[6] << "\n" 
00594                               << "nXTalk = " << stats2[7] << endl;
00595                             nplanes -= (stats2[1] + stats2[2]);
00596                             if( (Float_t(stats2[4])/Float_t(nplanes) > pHardBuriedFrac) ||
00597                                 (Float_t(stats2[3]-stats2[5]-stats2[6])/Float_t(nplanes) > 
00598                                  pHardBuriedFrac) ||
00599                                 (stats2[3]>=Int_t(nplanes*pSoftBuriedFrac) && 
00600                                  (TMath::Abs(track2->GetDirCosU())>pDCosUVCut || 
00601                                   TMath::Abs(track2->GetDirCosV())>pDCosUVCut) ) ) {
00602                                track2Buried = true;
00603                               objectMatch[i][j]=-1;
00604                               objectMatch[j][i]=-1;
00605                             }
00606                           }
00607                           if(!track2Buried){
00608                             //test whether track2 matches better than track
00609                             Int_t *stats = new Int_t[8];
00610                             shower2->BuriedTrack(track,1000,stats);
00611                             if( ( stats2[3] - stats2[5] - stats2[6] ) >=
00612                                 ( stats[3]  - stats[5]  - stats[6]  ) && 
00613                                 ( stats2[4] >= stats[4] ) ) {
00614                               otherTrackMatch = true;
00615                             }
00616                             delete [] stats;
00617                           }
00618                           delete [] stats2;
00619                         }
00620                       }
00621                     }
00622                     else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00623                       Int_t *stats = new Int_t[8];
00624                       shower2->BuriedTrack(track,1000,stats);
00625                       if(stats[4]==0) { //how convincing is first track/shower match?
00626                         otherShowerMatch = true;
00627                       }
00628                       delete [] stats;
00629                     }
00630                   }
00631                 }
00632                 
00633                 //if there are no other viable matches:
00634                 if(otherTrackVertexMatch || otherShowerMatch || 
00635                    otherTrackMatch) {
00636                   MSG("EventSS",Msg::kDebug) << "Shower " << i << " has better association: "
00637                                              << " otherTrackVertexMatch = " 
00638                                              << otherTrackVertexMatch 
00639                                              << " otherShowerMatch = " 
00640                                              << otherShowerMatch 
00641                                              << " otherTrackMatch = "
00642                                              << otherTrackMatch << endl;
00643                   continue;
00644                 }
00645                 //then add this shower to this track
00646                 MSG("EventSS",Msg::kDebug) << "Shower " << i << " has no better association. "
00647                                            << "Add to event " << nEvents << endl;
00648                 event->Add(recolist.At(i));
00649                 ph[i] = 0;
00650                 objectUsed[i] = true;
00651                 nUnusedObjects-=1;                
00652               } //end if for non-vertex shower
00653             } //end if for CandShowerSRHandle
00654           } //end if for valid association
00655         } //end for loop over all objects
00656       } //end if maxPH object is a track
00657       
00658       else if(recolist.At(maxPHIndex)->InheritsFrom("CandShowerSRHandle")) {
00659         MSG("EventSS",Msg::kDebug) << "Object " << maxPHIndex << " is a shower" << endl;
00660         MSG("EventSS",Msg::kDebug) << "Check whether shower " << maxPHIndex 
00661                                    << " has other associations" << endl;
00662         //maxPH object is a shower
00663         //CandShowerSRHandle *shower = dynamic_cast<CandShowerSRHandle*>(recolist.At(maxPHIndex));
00664         //now add objects clearly associated with largest object:
00665         for(int i=0;i<nRecoObjects;i++){
00666           if(objectMatch[maxPHIndex][i]!=0 && !objectUsed[i]) {   
00667             if(recolist.At(i)->InheritsFrom("CandTrackHandle")) {
00668               MSG("EventSS",Msg::kWarning) << "Arrgh - Logic flaw in eventSR. "
00669                                            << "Seeing tracks in event building, "
00670                                            << "when all tracks should already be "
00671                                            << "assigned!" << endl;
00672             } //end if CandTrackHandle
00673             else if(recolist.At(i)->InheritsFrom("CandShowerSRHandle")) {
00674               //check for two more levels of associated showers
00675               //CandShowerSRHandle *shower2 = dynamic_cast<CandShowerSRHandle*>(recolist.At(i));
00676               MSG("EventSS",Msg::kDebug) << "Adding shower " << i << " to event " 
00677                                          << nEvents << endl;
00678               MSG("EventSS",Msg::kDebug) << "Check whether shower " << i
00679                                          << " has other associations" << endl;
00680               event->Add(recolist.At(i));
00681               ph[i] = 0;
00682               objectUsed[i] = true;
00683               nUnusedObjects-=1;
00684               for(int j=0;j<nRecoObjects;j++){
00685                 if(j!=maxPHIndex && objectMatch[i][j]!=0 && 
00686                    !objectUsed[j] && objectMatch[maxPHIndex][j]==0) {
00687                   if(recolist.At(j)->InheritsFrom("CandTrackHandle")) {
00688                     MSG("EventSS",Msg::kWarning) << "Arrgh again - Logic flaw in eventSR. "
00689                                                  << "Seeing tracks in event building, "
00690                                                  << "when all tracks should already be "
00691                                                  << "assigned!" << endl;
00692                   }
00693                   else if(recolist.At(j)->InheritsFrom("CandShowerSRHandle")) {
00694                     MSG("EventSS",Msg::kDebug) << "Adding shower " << j << " to event " 
00695                                                << nEvents << endl;
00696                     event->Add(recolist.At(j));
00697                     ph[j] = 0;
00698                     objectUsed[j] = true;
00699                     nUnusedObjects-=1;
00700                   }
00701                 }
00702               }
00703             } //end if CandShowerSRHandle
00704           } //end if valid association
00705         } //end loop over all objects
00706       } //end if maxPH object is a shower
00707         //add event to eventlist
00708       eventlist.Add(event);
00709       nEvents+=1;
00710     }
00711     
00712     // remove all empty events from event list
00713 
00714     for (Int_t i=0; i<=eventlist.GetLast(); i++) {
00715       TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i));
00716       if (objectlist->GetLast()<0) eventlist.Remove(objectlist);
00717       else if(objectlist->GetLast()==0) { //consider single shower/track events
00718         if(objectlist->At(0)->InheritsFrom("CandShowerSRHandle")) {
00719           MSG("EventSS",Msg::kDebug) << "Have a single candshower event " << endl;
00720           CandShowerSRHandle *showerSR = 
00721             dynamic_cast<CandShowerSRHandle*> (objectlist->At(0));
00722           //if shower seems unphysical, remove it
00723           if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00724             eventlist.Remove(objectlist);
00725             MSG("EventSS",Msg::kDebug) << "REMOVING SHOWER-ONLY EVENT" << endl;
00726             continue;
00727           }
00728         }
00729         else if(objectlist->At(0)->InheritsFrom("CandTrackHandle")) {     
00730           CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(0));
00731           if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) {
00732             //if track seems unphysical, remove it
00733            
00734             eventlist.Remove(objectlist);
00735             MSG("EventSS",Msg::kDebug) << "REMOVING TRACK-ONLY EVENT" << endl;
00736             continue;
00737           }
00738         }
00739       }
00740       else if(objectlist->GetLast()==1) { //multiple shower/track events
00741         for(int ind_loop = 0;ind_loop<2;ind_loop++){
00742           Int_t ind_loop1 = ind_loop;
00743           Int_t ind_loop2 = ind_loop+1; if(ind_loop2==2) ind_loop2 = 0;
00744           if(objectlist->At(ind_loop1)->InheritsFrom("CandShowerSRHandle")){
00745             CandShowerSRHandle *showerSR = 
00746               dynamic_cast<CandShowerSRHandle*> (objectlist->At(ind_loop1));
00747             if(showerSR->IsUnphysical(pShwXTalkFracCut,pShwXTalkDef)) {
00748               if(objectlist->At(ind_loop2)->InheritsFrom("CandTrackHandle")){
00749                 CandTrackHandle *track = dynamic_cast <CandTrackHandle*> (objectlist->At(ind_loop2));
00750                 if(track->IsUnphysical(pTrkGapFracCut,pTrkAsymCut,pTrkXTalkFracCut,pTrkXTalkDef)) {
00751                   eventlist.Remove(objectlist);
00752                   MSG("EventSS",Msg::kDebug) << "REMOVING UNPHYSICAL SHOWER+TRACK EVENT" << endl;
00753                   continue;
00754                 }
00755               }
00756             }
00757           }
00758         }
00759       }
00760       if(useChopper==1) {
00761         std::vector<CandHandle> evtVec;
00762         for(int i=0;i<=objectlist->GetLast();i++){
00763           CandHandle *event = dynamic_cast<CandHandle*>(objectlist->At(i));
00764           evtVec.push_back(*event);       
00765         }
00766         ChopHelper chopper(cx.GetMom());
00767         ChopHelp *chop = chopper.GetChopHelp(evtVec);           
00768         if(chop->nchop==1) continue;
00769         else {
00770           chop->Print();
00771           //do something with chop info
00772           //for(int i=0;i<ncandidates;i++){
00773           //}
00774         }
00775       }
00776     }
00777     eventlist.Compress();
00778 
00779 
00780     // now make CandEvents for each event we have found
00781     Bool_t buildEvent=false;
00782     for (Int_t i=0; i<=eventlist.GetLast(); i++) {
00783       TObjArray *objectlist = dynamic_cast<TObjArray*>(eventlist.At(i));
00784       MSG("EventSS",Msg::kDebug)  << "making event of " 
00785                                     << objectlist->GetLast()+1 
00786                                     << " objects \n";
00787       cxx.SetDataIn(objectlist);
00788       buildEvent=true;
00789       CandEventHandle eventhandle = CandEvent::MakeCandidate(ah,cxx);
00790       MSG("EventSS",Msg::kDebug) << " # of strips:  " 
00791                                    << eventhandle.GetNDaughters() << "\n";
00792       eventhandle.SetCandSlice(slice);
00793       ch.AddDaughterLink(eventhandle);
00794       delete eventlist.At(i);
00795     }
00796     if (!buildEvent) BuildEventFromUnassoc(ac,ch,cx,slice);
00797   } // end loop over slice list
00798 
00799   delete trackItr;
00800   delete showerItr;
00801   delete subshowerItr;
00802 
00803   MSG("EventSS",Msg::kDebug) << "starting unassociated hits assignment \n";
00804 
00805   // find unassociated hits, and place them in tobjarray
00806   sliceItr.Reset();
00807   TObjArray unassociated;
00808   CandEventListHandle &eventlist = dynamic_cast<CandEventListHandle &>(ch);
00809   CandEventHandleItr * eventItr2 = 0;
00810   if (eventlist.GetNDaughters()>0) {    
00811     eventItr2 = new CandEventHandleItr(eventlist.GetDaughterIterator());
00812   }
00813   FindUnassociated(ac, sliceItr,eventItr2,unassociated);
00814   delete eventItr2;
00815   SetPrimaryShowers(ch,ac);
00816   
00817   // now loop over unnassociated hits, and check for adjacency with all
00818   // showers, and track vertices (in the case that no shower exists).  
00819   // If a match is found in the former case, add strip to shower. 
00820   // In the latter case, create a new shower.
00821  
00822   // this map contains the separation between each event/unassociated hit
00823   vector<map<const CandEventHandle*, Double_t> > dist2map(unassociated.GetLast()+1);
00824   FillDist2Map(ac,unassociated,ch,dist2map);
00825   
00826   // initialize a vector which holds a Bool_t
00827   // indicating whether a hit has already been used
00828   vector<Bool_t> alreadyUsed(unassociated.GetLast()+1);
00829   for (Int_t i=0; i<=unassociated.GetLast(); i++) alreadyUsed[i] = false;  
00830 
00831   // loop over unnassociated hits, adding them to events.
00832   // After each loop the affected distance map entries are recalculated, 
00833   // and the  loop is repeated.  This continues until a loop adds no strips.
00834   Bool_t found=true;
00835   Int_t n_found=0;
00836   while (found) {
00837     
00838     //found indicates whether a strip was added to an event in this loop
00839     found=false;
00840     
00841     CandStripHandle *closeststrip= 0;
00842     CandEventHandle *closestevent= 0;
00843     Double_t closestdist2 = 9999999.;
00844     Int_t closesti=-1;
00845     
00846     // loop over unassoc strips
00847     for (Int_t i=0; i<=unassociated.GetLast(); i++) {
00848       if(alreadyUsed[i]) continue;
00849 
00850       CandStripHandle *strip =
00851         dynamic_cast<CandStripHandle*>(unassociated.At(i));
00852       
00853       CandEventHandle *bestevent = 0;
00854       Double_t bestdist2 = 9999999.;
00855       // loop over events, and find event which is closest to strip
00856       // and in time
00857       TIter eventItr(ch.GetDaughterIterator());
00858       while (CandEventHandle *event =
00859              dynamic_cast<CandEventHandle*>(eventItr())) {
00860         Double_t dtime = strip->GetBegTime()-event->GetVtxT();
00861         if ( dtime>pHitAssocDTime0 && dtime<pHitAssocDTime1 &&
00862              (!bestevent || dist2map[i][event]<bestdist2) ) {
00863           bestevent = event;
00864           bestdist2 = dist2map[i][event];
00865         }
00866       }
00867       if (bestevent && (!closeststrip || bestdist2<closestdist2)) {
00868         closeststrip = strip;
00869         closesti=i;
00870         closestdist2 = bestdist2;
00871         closestevent = bestevent;
00872       }
00873     }
00874     MSG("EventSS",Msg::kVerbose)<< "closest strip " << closestdist2
00875                                 << "/"<<pHitAssocMaxDist2<< endl;
00876 
00877     // if distance to closest event less than maximum allowable add to event   
00878     if(closeststrip && closestevent && closestdist2<pHitAssocMaxDist2){
00879       MSG("EventSS",Msg::kVerbose) << "closest strip "
00880                                    << closeststrip->GetPlane() << " "
00881                                    << closeststrip->GetTPos()
00882                                    << " " << closestdist2 << "\n";       
00883       alreadyUsed[closesti] = true;
00884       found = true;
00885       ++n_found;
00886       
00887       // now either add this strip to existing primary shower, or
00888       // create a new shower at the track vertex if one does not exist
00889       // Get PrimaryShower CandHandle from fShowerList for modification.
00890       CandShowerHandle *shower = 0;
00891       if (CandShowerHandle *primsh = closestevent->GetPrimaryShower()) {
00892         for (Int_t ishower=0; 
00893              !shower && ishower<closestevent->GetLastShower()+1; 
00894              ishower++) {
00895           const CandShowerHandle *target = closestevent->GetShower(ishower);
00896           if (primsh->IsCloneOf(*target))     // Tests clone or Cand ==
00897             shower = closestevent->GetShowerWritable(ishower);
00898         }
00899       }
00900 
00901       Float_t minShwPlane=0;
00902       Float_t maxShwPlane=0;
00903       if(shower){
00904         minShwPlane=min(shower->GetVtxPlane(),shower->GetEndPlane());
00905         maxShwPlane=max(shower->GetVtxPlane(),shower->GetEndPlane());
00906       }
00907       Float_t vertexsep = 0;
00908       CandTrackHandle *track = closestevent->GetPrimaryTrack();
00909       if (shower && track) {
00910         vertexsep = abs(shower->GetBegPlane() - track->GetBegPlane());
00911         vertexsep *= PITCHINMETRES;
00912         MSG("EventSS",Msg::kVerbose) << "Found track and primary shower. \n" 
00913                                      << "Vertex separation (limit) = " 
00914                                      << vertexsep << "(" 
00915                                      << pShwTrkDz << ")" << " \n";
00916       }
00917       // if event has shower and either no track, or shower is close
00918       // to track vertex, add this hit to the shower which was
00919       // compared in position to this strip earler.
00920       if (shower && (!track || vertexsep<pShwTrkDz)) {
00921         // NOTE: CandShowerHandle *shower, used to modify CandShower, is owned
00922         // by fShowerList.
00923         // add new strip to this shower, and the appropriate cluster
00924         MSG("EventSS",Msg::kVerbose) << "Shower min/max plane (limit): "
00925                                      << minShwPlane << "/" << maxShwPlane
00926                                      << " (" << pMaxNewShwLen << ")" << "\n";
00927         if(closeststrip->GetPlane()>minShwPlane-pMaxNewShwLen &&
00928            closeststrip->GetPlane()<maxShwPlane+pMaxNewShwLen ){
00929           MSG("EventSS",Msg::kVerbose) << "adding strip to shower \n"; 
00930           AddStripToEvent(closestevent,shower,subshowerlist,closeststrip);
00931         }
00932       }
00933       // no shower consistent with track vertex, and this strip
00934       // is close to vertex, so start new shower.    
00935       else if (track) CreatePrimaryShower(ac,cxx,closestevent,showerlist,
00936                                           subshowerlist,closeststrip);      
00937       
00938       // need to recalculate distances bewteen remaining
00939       // unassoc hits and this modified event, if the added strip
00940       // results in a smaller separation, replace value in dist2map 
00941       ReFillDist2Map(ac,closeststrip,closestevent,unassociated,dist2map);      
00942     } // if adding strip to event in form of new shower
00943   } // if found (if an unassoc hit was matched to event in last loop over them
00944 
00945   // because we have added strips to previously constructed showers, their
00946   // runalg method should be run on them again.  In the future, an
00947   // optimization would be to do this only on showers which have been
00948   // touched.
00949   MSG("EventSS",Msg::kVerbose) << "added " << n_found << " of "
00950                                << unassociated.GetLast()+1 
00951                                << " unassociated hits " << endl;
00952   
00953   cxx.SetDataIn(showerlist);  // Pass along persistable CandShowerList
00954   ReConstructShowers(ac,ch,cxx);
00955 
00956   //Final step: use the chopper to combine events across slices:
00957   if(useChopper==1) {
00958     std::vector<CandHandle> evtVec;
00959     CandEventHandleItr evitr(ch.GetDaughterIterator());
00960     while(CandEventHandle *event = dynamic_cast<CandEventHandle*>(evitr())){
00961       evtVec.push_back(*event);
00962     }
00963     ChopHelper chopper(cx.GetMom());
00964     ChopHelp *chop = chopper.GetChopHelp(evtVec);
00965     chop->Print();
00966   }
00967 
00968   // For each event, set primary shower based on 
00969   // distance from vertex and energy
00970   SetPrimaryShowers(ch,ac);
00971 }

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(), and Registry::GetDouble().

Referenced by RunAlg().

01381 {
01382   
01383   // This method loops over events in the eventlist, finding the largest
01384   // shower and setting that to the primary shower in the event.
01385   
01386   Double_t pShwShwDz = ac.GetDouble("ShwShwDz");
01387   Double_t pMinShwEFract = ac.GetDouble("MinShwEFract");
01388   
01389   TIter evItr(ch.GetDaughterIterator());
01390   while (CandEventHandle *ev =
01391                               dynamic_cast<CandEventHandle*>(evItr())) {
01392     ev->SetPrimaryShower(pMinShwEFract,pShwShwDz);
01393   }
01394 }

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

Reimplemented from AlgBase.

Definition at line 2033 of file AlgEventSSList.cxx.

02034 {
02035 }


The documentation for this class was generated from the following files:
Generated on Wed Dec 10 22:49:03 2014 for loon by  doxygen 1.4.7