RecoUtils/src/PFCand_AssoMapAlgos.cc

0001
0002 #include "CommonTools/RecoUtils/interface/PFCand_AssoMapAlgos.h"
0003
0004 #include "TrackingTools/IPTools/interface/IPTools.h"
0005
0006 using namespace edm;
0007 using namespace std;
0008 using namespace reco;
0009
0010 /*************************************************************************************/
0011 /* dedicated constructor for the algorithms                                          */
0012 /*************************************************************************************/
0013
0014 PFCand_AssoMapAlgos::PFCand_AssoMapAlgos(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& iC)
0015     : PF_PU_AssoMapAlgos(iConfig, iC), token_bField_(iC.esConsumes()), token_TrackingGeometry_(iC.esConsumes()) {
0016   input_MaxNumAssociations_ = iConfig.getParameter<int>("MaxNumberOfAssociations");
0017
0018   token_VertexCollection_ = iC.consumes<VertexCollection>(iConfig.getParameter<InputTag>("VertexCollection"));
0019
0020   token_BeamSpot_ = iC.consumes<BeamSpot>(iConfig.getParameter<InputTag>("BeamSpot"));
0021 }
0022
0023 /*************************************************************************************/
0024 /* get all needed collections at the beginning                                       */
0025 /*************************************************************************************/
0026
0027 void PFCand_AssoMapAlgos::GetInputCollections(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0028   PF_PU_AssoMapAlgos::GetInputCollections(iEvent, iSetup);
0029
0030   //get the offline beam spot
0031   iEvent.getByToken(token_BeamSpot_, beamspotH);
0032
0033   //get the input vertex collection
0034   iEvent.getByToken(token_VertexCollection_, vtxcollH);
0035
0036   bFieldH = iSetup.getHandle(token_bField_);
0037   trackingGeometryH = iSetup.getHandle(token_TrackingGeometry_);
0038 }
0039
0040 /*************************************************************************************/
0041 /* create the pf candidate to vertex association and the inverse map                 */
0042 /*************************************************************************************/
0043 std::pair<std::unique_ptr<PFCandToVertexAssMap>, std::unique_ptr<VertexToPFCandAssMap>>
0044 PFCand_AssoMapAlgos::createMappings(edm::Handle<reco::PFCandidateCollection> pfCandH) {
0045   unique_ptr<PFCandToVertexAssMap> pfcand2vertex(new PFCandToVertexAssMap(vtxcollH, pfCandH));
0046   unique_ptr<VertexToPFCandAssMap> vertex2pfcand(new VertexToPFCandAssMap(pfCandH, vtxcollH));
0047
0048   int num_vertices = vtxcollH->size();
0049   if (num_vertices < input_MaxNumAssociations_)
0050     input_MaxNumAssociations_ = num_vertices;
0051   vector<VertexRef> vtxColl_help;
0052   if (input_MaxNumAssociations_ == 1)
0053     vtxColl_help = CreateVertexVector(vtxcollH);
0054
0055   for (unsigned i = 0; i < pfCandH->size(); i++) {
0056     PFCandidateRef candref(pfCandH, i);
0057
0058     if (input_MaxNumAssociations_ > 1)
0059       vtxColl_help = CreateVertexVector(vtxcollH);
0060
0061     VertexPfcQuality VtxPfcQual;
0062
0063     TrackRef PFCtrackref = candref->trackRef();
0064
0065     if (PFCtrackref.isNull()) {
0066       for (int assoc_ite = 0; assoc_ite < input_MaxNumAssociations_; ++assoc_ite) {
0067         int quality = -1 - assoc_ite;
0068
0069         // Insert the best vertex and the pair of track and the quality of this association in the map
0070         pfcand2vertex->insert(vtxColl_help.at(0), make_pair(candref, quality));
0071         vertex2pfcand->insert(candref, make_pair(vtxColl_help.at(0), quality));
0072
0073         //cleanup only if multiple iterations are made
0074         if (input_MaxNumAssociations_ > 1)
0075           PF_PU_AssoMapAlgos::EraseVertex(vtxColl_help, vtxColl_help.at(0));
0076       }
0077
0078     } else {
0079       TransientTrack transtrk(PFCtrackref, &(*bFieldH));
0080       transtrk.setBeamSpot(*beamspotH);
0081       transtrk.setTrackingGeometry(trackingGeometryH);
0082
0083       for (int assoc_ite = 0; assoc_ite < input_MaxNumAssociations_; ++assoc_ite) {
0084         VertexStepPair assocVtx =
0085             FindAssociation(PFCtrackref, vtxColl_help, bFieldH, trackingGeometryH, beamspotH, assoc_ite);
0086         int step = assocVtx.second;
0087         double distance = (IPTools::absoluteImpactParameter3D(transtrk, *(assocVtx.first))).second.value();
0088
0089         int quality = DefineQuality(assoc_ite, step, distance);
0090
0091         // Insert the best vertex and the pair of track and the quality of this association in the map
0092         pfcand2vertex->insert(assocVtx.first, make_pair(candref, quality));
0093         vertex2pfcand->insert(candref, make_pair(assocVtx.first, quality));
0094
0095         //cleanup only if multiple iterations are made
0096         if (input_MaxNumAssociations_ > 2)
0097           PF_PU_AssoMapAlgos::EraseVertex(vtxColl_help, assocVtx.first);
0098       }
0099
0100     }  //check PFCtrackref.isNull
0101
0102   }  //i on pfCandH
0103
0104   return {std::move(pfcand2vertex), std::move(vertex2pfcand)};
0105 }
0106
0107 /*************************************************************************************/
0108 /* create the pf candidate to vertex association map                                 */
0109 /*************************************************************************************/
0110
0111 std::unique_ptr<PFCandToVertexAssMap> PFCand_AssoMapAlgos::CreatePFCandToVertexMap(
0112     edm::Handle<reco::PFCandidateCollection> pfCandH) {
0113   return createMappings(pfCandH).first;
0114 }
0115
0116 /*************************************************************************************/
0117 /* create the vertex to pf candidate association map                                 */
0118 /*************************************************************************************/
0119
0120 std::unique_ptr<VertexToPFCandAssMap> PFCand_AssoMapAlgos::CreateVertexToPFCandMap(
0121     edm::Handle<reco::PFCandidateCollection> pfCandH) {
0122   return createMappings(pfCandH).second;
0123 }
0124
0125 /*************************************************************************************/
0126 /* create the vertex to pf candidate association map                                 */
0127 /*************************************************************************************/
0128
0129 std::unique_ptr<PFCandToVertexAssMap> PFCand_AssoMapAlgos::SortPFCandAssociationMap(
0130     PFCandToVertexAssMap* pfcvertexassInput, edm::EDProductGetter const* getter) {
0131   //create a new PFCandVertexAssMap for the Output which will be sorted
0132   unique_ptr<PFCandToVertexAssMap> pfcvertexassOutput(new PFCandToVertexAssMap(getter));
0133
0134   //Create and fill a vector of pairs of vertex and the summed (pT)**2 of the pfcandidates associated to the vertex
0135   VertexPtsumVector vertexptsumvector;
0136
0137   //loop over all vertices in the association map
0138   for (PFCandToVertexAssMap::const_iterator assomap_ite = pfcvertexassInput->begin();
0139        assomap_ite != pfcvertexassInput->end();
0140        assomap_ite++) {
0141     const VertexRef assomap_vertexref = assomap_ite->key;
0142     const PFCandQualityPairVector pfccoll = assomap_ite->val;
0143
0144     float ptsum = 0;
0145
0146     PFCandidateRef pfcandref;
0147
0148     //get the pfcandidates associated to the vertex and calculate the pT**2
0149     for (unsigned int pfccoll_ite = 0; pfccoll_ite < pfccoll.size(); pfccoll_ite++) {
0150       pfcandref = pfccoll[pfccoll_ite].first;
0151       int quality = pfccoll[pfccoll_ite].second;
0152
0153       if ((quality <= 2) && (quality != -1))
0154         continue;
0155
0156       double man_pT = pfcandref->pt();
0157       if (man_pT > 0.)
0158         ptsum += man_pT * man_pT;
0159     }
0160
0161     vertexptsumvector.push_back(make_pair(assomap_vertexref, ptsum));
0162   }
0163
0164   while (!vertexptsumvector.empty()) {
0165     VertexRef vertexref_highestpT;
0166     float highestpT = 0.;
0167     int highestpT_index = 0;
0168
0169     for (unsigned int vtxptsumvec_ite = 0; vtxptsumvec_ite < vertexptsumvector.size(); vtxptsumvec_ite++) {
0170       if (vertexptsumvector[vtxptsumvec_ite].second > highestpT) {
0171         vertexref_highestpT = vertexptsumvector[vtxptsumvec_ite].first;
0172         highestpT = vertexptsumvector[vtxptsumvec_ite].second;
0173         highestpT_index = vtxptsumvec_ite;
0174       }
0175     }
0176
0177     //loop over all vertices in the association map
0178     for (PFCandToVertexAssMap::const_iterator assomap_ite = pfcvertexassInput->begin();
0179          assomap_ite != pfcvertexassInput->end();
0180          assomap_ite++) {
0181       const VertexRef assomap_vertexref = assomap_ite->key;
0182       const PFCandQualityPairVector pfccoll = assomap_ite->val;
0183
0184       //if the vertex from the association map the vertex with the highest pT
0185       //insert all associated pfcandidates in the output Association Map
0186       if (assomap_vertexref == vertexref_highestpT)
0187         for (unsigned int pfccoll_ite = 0; pfccoll_ite < pfccoll.size(); pfccoll_ite++)
0188           pfcvertexassOutput->insert(assomap_vertexref, pfccoll[pfccoll_ite]);
0189     }
0190
0191     vertexptsumvector.erase(vertexptsumvector.begin() + highestpT_index);
0192   }
0193
0194   return pfcvertexassOutput;
0195 }