CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Attributes
GSFTrackImporter Class Reference
Inheritance diagram for GSFTrackImporter:
BlockElementImporterBase

Public Member Functions

 GSFTrackImporter (const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)
 
void importToBlock (const edm::Event &, ElementList &) const override
 
- Public Member Functions inherited from BlockElementImporterBase
 BlockElementImporterBase (const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)
 
 BlockElementImporterBase (const BlockElementImporterBase &)=delete
 
const std::string & name () const
 
BlockElementImporterBaseoperator= (const BlockElementImporterBase &)=delete
 
virtual void updateEventSetup (const edm::EventSetup &)
 
virtual ~BlockElementImporterBase ()=default
 

Private Attributes

const bool _isSecondary
 
edm::EDGetTokenT< reco::GsfPFRecTrackCollection_src
 
const bool _superClustersArePF
 

Additional Inherited Members

- Public Types inherited from BlockElementImporterBase
typedef std::vector< std::unique_ptr< reco::PFBlockElement > > ElementList
 

Detailed Description

Definition at line 14 of file GSFTrackImporter.cc.

Constructor & Destructor Documentation

GSFTrackImporter::GSFTrackImporter ( const edm::ParameterSet conf,
edm::ConsumesCollector sumes 
)
inline

Definition at line 16 of file GSFTrackImporter.cc.

References importToBlock().

17  : BlockElementImporterBase(conf, sumes),
19  _isSecondary(conf.getParameter<bool>("gsfsAreSecondary")),
20  _superClustersArePF(conf.getParameter<bool>("superClustersArePF")) {}
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
T getParameter(std::string const &) const
const bool _isSecondary
std::vector< GsfPFRecTrack > GsfPFRecTrackCollection
collection of GsfPFRecTrack objects
edm::EDGetTokenT< reco::GsfPFRecTrackCollection > _src
BlockElementImporterBase(const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)
const bool _superClustersArePF

Member Function Documentation

void GSFTrackImporter::importToBlock ( const edm::Event e,
BlockElementImporterBase::ElementList elems 
) const
overridevirtual

Implements BlockElementImporterBase.

Definition at line 31 of file GSFTrackImporter.cc.

References _isSecondary, _src, _superClustersArePF, a, HLT_2018_cff::distance, DetId::Forward, edm::Event::getHandle(), DetId::HGCalEE, DetId::HGCalHSc, DetId::HGCalHSi, edm::Ref< C, T, F >::isAvailable(), edm::Ref< C, T, F >::isNonnull(), reco::PFBlockElement::SC, reco::PFBlockElementSuperCluster::setFromGsfElectron(), reco::PFBlockElementSuperCluster::setFromPFSuperCluster(), reco::PFBlockElementGsfTrack::setTrackType(), reco::PFBlockElement::T_FROM_GAMMACONV, and groupFilesInBlocks::temp.

Referenced by GSFTrackImporter().

31  {
32  auto gsftracks = e.getHandle(_src);
33  elems.reserve(elems.size() + gsftracks->size());
34  // setup our elements so that all the SCs are grouped together
35  auto SCs_end =
36  std::partition(elems.begin(), elems.end(), [](const auto& a) { return a->type() == reco::PFBlockElement::SC; });
37  // insert gsf tracks and SCs, binding pre-existing SCs to ECAL-Driven GSF
38  auto bgsf = gsftracks->cbegin();
39  auto egsf = gsftracks->cend();
40  for (auto gsftrack = bgsf; gsftrack != egsf; ++gsftrack) {
41  reco::GsfPFRecTrackRef gsfref(gsftracks, std::distance(bgsf, gsftrack));
42  const reco::GsfTrackRef& basegsfref = gsftrack->gsfTrackRef();
43  const auto& gsfextraref = basegsfref->extra();
44  // import associated super clusters
45  if (gsfextraref.isAvailable() && gsfextraref->seedRef().isAvailable()) {
46  reco::ElectronSeedRef seedref = gsfextraref->seedRef().castTo<reco::ElectronSeedRef>();
47  if (seedref.isAvailable() && seedref->isEcalDriven()) {
48  reco::SuperClusterRef scref = seedref->caloCluster().castTo<reco::SuperClusterRef>();
49  if (scref.isNonnull()) {
50  // explicitly veto HGCal super clusters
51  if (scref->seed()->seed().det() == DetId::HGCalEE || scref->seed()->seed().det() == DetId::HGCalHSi ||
52  scref->seed()->seed().det() == DetId::HGCalHSc || scref->seed()->seed().det() == DetId::Forward)
53  continue;
54  PFBlockElementSCEqual myEqual(scref);
55  auto sc_elem = std::find_if(elems.begin(), SCs_end, myEqual);
56  if (sc_elem != SCs_end) {
57  reco::PFBlockElementSuperCluster* scbe = static_cast<reco::PFBlockElementSuperCluster*>(sc_elem->get());
58  scbe->setFromGsfElectron(true);
59  } else {
61  scbe->setFromGsfElectron(true);
63  SCs_end = elems.emplace(SCs_end, scbe);
64  ++SCs_end; // point to element *after* the new one
65  }
66  }
67  }
68  } // gsf extra ref?
69  // cache the SC_end offset
70  size_t SCs_end_position = std::distance(elems.begin(), SCs_end);
71  // get track momentum information
72  const std::vector<reco::PFTrajectoryPoint>& PfGsfPoint = gsftrack->trajectoryPoints();
73  unsigned int c_gsf = 0;
74  bool PassTracker = false;
75  unsigned int IndexPout = 0;
76  for (const auto& pfGsfPoint : PfGsfPoint) {
77  if (pfGsfPoint.isValid()) {
78  int layGsfP = pfGsfPoint.layer();
79  if (layGsfP == -1)
80  PassTracker = true;
81  else if (PassTracker && layGsfP > 0) {
82  IndexPout = c_gsf - 1;
83  break;
84  }
85  ++c_gsf;
86  }
87  }
88  const math::XYZTLorentzVector& pin = PfGsfPoint[0].momentum();
89  const math::XYZTLorentzVector& pout = PfGsfPoint[IndexPout].momentum();
91  if (_isSecondary) {
93  }
94  elems.emplace_back(temp);
95  // import brems from this primary gsf
96  for (const auto& brem : gsfref->PFRecBrem()) {
97  const unsigned TrajP = brem.indTrajPoint();
98  if (TrajP != 99) {
99  const double DP = brem.DeltaP();
100  const double sDP = brem.SigmaDeltaP();
101  elems.emplace_back(new reco::PFBlockElementBrem(gsfref, DP, sDP, TrajP));
102  }
103  }
104  // protect against reallocations, create a fresh iterator
105  SCs_end = elems.begin() + SCs_end_position;
106  } // loop on gsf tracks
107  elems.shrink_to_fit();
108 }
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
const bool _isSecondary
void setFromGsfElectron(bool val)
set provenance
bool isAvailable() const
Definition: Ref.h:537
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
Handle< PROD > getHandle(EDGetTokenT< PROD > token) const
Definition: Event.h:547
void setTrackType(TrackType trType, bool value) override
the trackType
double a
Definition: hdecay.h:119
edm::EDGetTokenT< reco::GsfPFRecTrackCollection > _src
const bool _superClustersArePF

Member Data Documentation

const bool GSFTrackImporter::_isSecondary
private

Definition at line 26 of file GSFTrackImporter.cc.

Referenced by importToBlock().

edm::EDGetTokenT<reco::GsfPFRecTrackCollection> GSFTrackImporter::_src
private

Definition at line 25 of file GSFTrackImporter.cc.

Referenced by importToBlock().

const bool GSFTrackImporter::_superClustersArePF
private

Definition at line 26 of file GSFTrackImporter.cc.

Referenced by importToBlock().