File indexing completed on 2021-07-15 08:48:08
0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "FWCore/Framework/interface/ESHandle.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "CommonTools/ParticleFlow/interface/PFClusterWidthAlgo.h"
0006 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0007 #include "DataFormats/EgammaReco/interface/PreshowerCluster.h"
0008 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0009 #include "DataFormats/EgammaCandidates/interface/PhotonCore.h"
0010 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0011 #include "DataFormats/VertexReco/interface/Vertex.h"
0012 #include "DataFormats/ParticleFlowReco/interface/PFBlockElement.h"
0013 #include "DataFormats/Math/interface/Vector3D.h"
0014 #include "DataFormats/Math/interface/LorentzVector.h"
0015 #include "RecoEcal/EgammaCoreTools/interface/Mustache.h"
0016 #include "FWCore/Framework/interface/stream/EDProducer.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018 #include "DataFormats/Common/interface/ValueMap.h"
0019 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0020 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0021 #include "DataFormats/Common/interface/Handle.h"
0022 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidatePhotonExtra.h"
0023 #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
0024 #include "CondFormats/EcalObjects/interface/EcalMustacheSCParameters.h"
0025 #include "CondFormats/DataRecord/interface/EcalMustacheSCParametersRcd.h"
0026
0027 class CaloGeometry;
0028 class CaloTopology;
0029 class DetId;
0030 namespace edm {
0031 class EventSetup;
0032 }
0033
0034 class PFPhotonTranslator : public edm::stream::EDProducer<> {
0035 public:
0036 explicit PFPhotonTranslator(const edm::ParameterSet &);
0037 ~PFPhotonTranslator() override;
0038
0039 void produce(edm::Event &, const edm::EventSetup &) override;
0040
0041 typedef std::vector<edm::Handle<edm::ValueMap<double> > > IsolationValueMaps;
0042
0043 private:
0044
0045 bool fetchCandidateCollection(edm::Handle<reco::PFCandidateCollection> &c,
0046 const edm::InputTag &tag,
0047 const edm::Event &iEvent) const;
0048
0049
0050
0051 void createBasicCluster(const reco::PFBlockElement &,
0052 reco::BasicClusterCollection &basicClusters,
0053 std::vector<const reco::PFCluster *> &,
0054 const reco::PFCandidate &coCandidate) const;
0055
0056 void createPreshowerCluster(const reco::PFBlockElement &PFBE,
0057 reco::PreshowerClusterCollection &preshowerClusters,
0058 unsigned plane) const;
0059
0060
0061 void createBasicClusterPtrs(const edm::OrphanHandle<reco::BasicClusterCollection> &basicClustersHandle);
0062
0063
0064 void createPreshowerClusterPtrs(const edm::OrphanHandle<reco::PreshowerClusterCollection> &preshowerClustersHandle);
0065
0066
0067 void createSuperClusters(const reco::PFCandidateCollection &, reco::SuperClusterCollection &superClusters) const;
0068
0069 void createOneLegConversions(const edm::OrphanHandle<reco::SuperClusterCollection> &superClustersHandle,
0070 reco::ConversionCollection &oneLegConversions);
0071
0072
0073 void createPhotonCores(const edm::OrphanHandle<reco::SuperClusterCollection> &superClustersHandle,
0074 const edm::OrphanHandle<reco::ConversionCollection> &oneLegConversionHandle,
0075 reco::PhotonCoreCollection &photonCores);
0076
0077 void createPhotons(reco::VertexCollection &vertexCollection,
0078 edm::Handle<reco::PhotonCollection> &egPhotons,
0079 const edm::OrphanHandle<reco::PhotonCoreCollection> &photonCoresHandle,
0080 const IsolationValueMaps &isolationValues,
0081 reco::PhotonCollection &photons);
0082
0083 const reco::PFCandidate &correspondingDaughterCandidate(const reco::PFCandidate &cand,
0084 const reco::PFBlockElement &pfbe) const;
0085
0086 edm::InputTag inputTagPFCandidates_;
0087 std::vector<edm::InputTag> inputTagIsoVals_;
0088 std::string PFBasicClusterCollection_;
0089 std::string PFPreshowerClusterCollection_;
0090 std::string PFSuperClusterCollection_;
0091 std::string PFPhotonCoreCollection_;
0092 std::string PFPhotonCollection_;
0093 std::string PFConversionCollection_;
0094 std::string EGPhotonCollection_;
0095 std::string vertexProducer_;
0096 edm::InputTag barrelEcalHits_;
0097 edm::InputTag endcapEcalHits_;
0098 edm::InputTag hcalTowers_;
0099 double hOverEConeSize_;
0100
0101
0102 std::vector<reco::BasicClusterCollection> basicClusters_;
0103
0104 std::vector<std::vector<const reco::PFCluster *> > pfClusters_;
0105
0106 std::vector<reco::PreshowerClusterCollection> preshowerClusters_;
0107
0108 std::vector<reco::SuperClusterCollection> superClusters_;
0109
0110 std::vector<reco::CaloClusterPtrVector> basicClusterPtr_;
0111
0112 std::vector<reco::CaloClusterPtrVector> preshowerClusterPtr_;
0113
0114 std::vector<int> photPFCandidateIndex_;
0115
0116 std::vector<reco::CandidatePtr> CandidatePtr_;
0117
0118 std::vector<reco::SuperClusterRef> egSCRef_;
0119
0120 std::vector<reco::PhotonRef> egPhotonRef_;
0121
0122 std::vector<float> pfPhotonMva_;
0123 std::vector<float> energyRegression_;
0124 std::vector<float> energyRegressionError_;
0125
0126
0127 std::vector<reco::ConversionRefVector> pfConv_;
0128 std::vector<std::vector<reco::TrackRef> > pfSingleLegConv_;
0129 std::vector<std::vector<float> > pfSingleLegConvMva_;
0130
0131 std::vector<int> conv1legPFCandidateIndex_;
0132 std::vector<int> conv2legPFCandidateIndex_;
0133
0134 edm::ESHandle<CaloTopology> theCaloTopo_;
0135 edm::ESHandle<CaloGeometry> theCaloGeom_;
0136
0137
0138 edm::ESGetToken<EcalMustacheSCParameters, EcalMustacheSCParametersRcd> ecalMustacheSCParametersToken_;
0139 const EcalMustacheSCParameters *mustacheSCParams_;
0140
0141 bool emptyIsOk_;
0142 };
0143
0144 DEFINE_FWK_MODULE(PFPhotonTranslator);
0145
0146 using namespace edm;
0147 using namespace std;
0148 using namespace reco;
0149
0150 typedef math::XYZTLorentzVector LorentzVector;
0151 typedef math::XYZPoint Point;
0152 typedef math::XYZVector Vector;
0153
0154 PFPhotonTranslator::PFPhotonTranslator(const edm::ParameterSet &iConfig) {
0155
0156
0157 inputTagPFCandidates_ = iConfig.getParameter<edm::InputTag>("PFCandidate");
0158
0159 edm::ParameterSet isoVals = iConfig.getParameter<edm::ParameterSet>("isolationValues");
0160 inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfChargedHadrons"));
0161 inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfPhotons"));
0162 inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfNeutralHadrons"));
0163
0164 PFBasicClusterCollection_ = iConfig.getParameter<std::string>("PFBasicClusters");
0165 PFPreshowerClusterCollection_ = iConfig.getParameter<std::string>("PFPreshowerClusters");
0166 PFSuperClusterCollection_ = iConfig.getParameter<std::string>("PFSuperClusters");
0167 PFConversionCollection_ = iConfig.getParameter<std::string>("PFConversionCollection");
0168 PFPhotonCoreCollection_ = iConfig.getParameter<std::string>("PFPhotonCores");
0169 PFPhotonCollection_ = iConfig.getParameter<std::string>("PFPhotons");
0170
0171 EGPhotonCollection_ = iConfig.getParameter<std::string>("EGPhotons");
0172
0173 vertexProducer_ = iConfig.getParameter<std::string>("primaryVertexProducer");
0174
0175 barrelEcalHits_ = iConfig.getParameter<edm::InputTag>("barrelEcalHits");
0176 endcapEcalHits_ = iConfig.getParameter<edm::InputTag>("endcapEcalHits");
0177
0178 hcalTowers_ = iConfig.getParameter<edm::InputTag>("hcalTowers");
0179 hOverEConeSize_ = iConfig.getParameter<double>("hOverEConeSize");
0180
0181 if (iConfig.exists("emptyIsOk"))
0182 emptyIsOk_ = iConfig.getParameter<bool>("emptyIsOk");
0183 else
0184 emptyIsOk_ = false;
0185
0186 ecalMustacheSCParametersToken_ = esConsumes<EcalMustacheSCParameters, EcalMustacheSCParametersRcd>();
0187
0188 produces<reco::BasicClusterCollection>(PFBasicClusterCollection_);
0189 produces<reco::PreshowerClusterCollection>(PFPreshowerClusterCollection_);
0190 produces<reco::SuperClusterCollection>(PFSuperClusterCollection_);
0191 produces<reco::PhotonCoreCollection>(PFPhotonCoreCollection_);
0192 produces<reco::PhotonCollection>(PFPhotonCollection_);
0193 produces<reco::ConversionCollection>(PFConversionCollection_);
0194 }
0195
0196 PFPhotonTranslator::~PFPhotonTranslator() {}
0197
0198 void PFPhotonTranslator::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0199
0200 mustacheSCParams_ = &iSetup.getData(ecalMustacheSCParametersToken_);
0201
0202 auto basicClusters_p = std::make_unique<reco::BasicClusterCollection>();
0203
0204 auto psClusters_p = std::make_unique<reco::PreshowerClusterCollection>();
0205
0206
0207
0208
0209
0210 reco::SuperClusterCollection outputSuperClusterCollection;
0211 reco::ConversionCollection outputOneLegConversionCollection;
0212 reco::PhotonCoreCollection outputPhotonCoreCollection;
0213 reco::PhotonCollection outputPhotonCollection;
0214
0215 outputSuperClusterCollection.clear();
0216 outputOneLegConversionCollection.clear();
0217 outputPhotonCoreCollection.clear();
0218 outputPhotonCollection.clear();
0219
0220 edm::Handle<reco::PFCandidateCollection> pfCandidates;
0221 bool status = fetchCandidateCollection(pfCandidates, inputTagPFCandidates_, iEvent);
0222
0223 edm::Handle<reco::PhotonCollection> egPhotons;
0224 iEvent.getByLabel(EGPhotonCollection_, egPhotons);
0225
0226 Handle<reco::VertexCollection> vertexHandle;
0227
0228 IsolationValueMaps isolationValues(inputTagIsoVals_.size());
0229 for (size_t j = 0; j < inputTagIsoVals_.size(); ++j) {
0230 iEvent.getByLabel(inputTagIsoVals_[j], isolationValues[j]);
0231 }
0232
0233
0234 photPFCandidateIndex_.clear();
0235 basicClusters_.clear();
0236 pfClusters_.clear();
0237 preshowerClusters_.clear();
0238 superClusters_.clear();
0239 basicClusterPtr_.clear();
0240 preshowerClusterPtr_.clear();
0241 CandidatePtr_.clear();
0242 egSCRef_.clear();
0243 egPhotonRef_.clear();
0244 pfPhotonMva_.clear();
0245 energyRegression_.clear();
0246 energyRegressionError_.clear();
0247 pfConv_.clear();
0248 pfSingleLegConv_.clear();
0249 pfSingleLegConvMva_.clear();
0250 conv1legPFCandidateIndex_.clear();
0251 conv2legPFCandidateIndex_.clear();
0252
0253
0254
0255
0256
0257
0258 unsigned ncand = (status) ? pfCandidates->size() : 0;
0259
0260 unsigned iphot = 0;
0261 unsigned iconv1leg = 0;
0262 unsigned iconv2leg = 0;
0263
0264 for (unsigned i = 0; i < ncand; ++i) {
0265 const reco::PFCandidate &cand = (*pfCandidates)[i];
0266 if (cand.particleId() != reco::PFCandidate::gamma)
0267 continue;
0268
0269 if (cand.mva_nothing_gamma() > 0.001)
0270 {
0271
0272
0273
0274
0275 if (!cand.photonExtraRef()->conversionRef().empty()) {
0276 pfConv_.push_back(reco::ConversionRefVector());
0277
0278 const reco::ConversionRefVector &doubleLegConvColl = cand.photonExtraRef()->conversionRef();
0279 for (unsigned int iconv = 0; iconv < doubleLegConvColl.size(); iconv++) {
0280 pfConv_[iconv2leg].push_back(doubleLegConvColl[iconv]);
0281 }
0282
0283 conv2legPFCandidateIndex_.push_back(iconv2leg);
0284 iconv2leg++;
0285 } else
0286 conv2legPFCandidateIndex_.push_back(-1);
0287
0288 const std::vector<reco::TrackRef> &singleLegConvColl = cand.photonExtraRef()->singleLegConvTrackRef();
0289 const std::vector<float> &singleLegConvCollMva = cand.photonExtraRef()->singleLegConvMva();
0290
0291
0292
0293 if (!singleLegConvColl.empty()) {
0294 pfSingleLegConv_.push_back(std::vector<reco::TrackRef>());
0295 pfSingleLegConvMva_.push_back(std::vector<float>());
0296
0297
0298 for (unsigned int itk = 0; itk < singleLegConvColl.size(); itk++) {
0299
0300
0301 pfSingleLegConv_[iconv1leg].push_back(singleLegConvColl[itk]);
0302 pfSingleLegConvMva_[iconv1leg].push_back(singleLegConvCollMva[itk]);
0303 }
0304
0305 conv1legPFCandidateIndex_.push_back(iconv1leg);
0306
0307 iconv1leg++;
0308 } else
0309 conv1legPFCandidateIndex_.push_back(-1);
0310 }
0311
0312 photPFCandidateIndex_.push_back(i);
0313 pfPhotonMva_.push_back(cand.mva_nothing_gamma());
0314 energyRegression_.push_back(cand.photonExtraRef()->MVAGlobalCorrE());
0315 energyRegressionError_.push_back(cand.photonExtraRef()->MVAGlobalCorrEError());
0316 basicClusters_.push_back(reco::BasicClusterCollection());
0317 pfClusters_.push_back(std::vector<const reco::PFCluster *>());
0318 preshowerClusters_.push_back(reco::PreshowerClusterCollection());
0319 superClusters_.push_back(reco::SuperClusterCollection());
0320
0321 reco::PFCandidatePtr ptrToPFPhoton(pfCandidates, i);
0322 CandidatePtr_.push_back(ptrToPFPhoton);
0323 egSCRef_.push_back(cand.superClusterRef());
0324
0325
0326 int iegphot = 0;
0327 for (reco::PhotonCollection::const_iterator gamIter = egPhotons->begin(); gamIter != egPhotons->end(); ++gamIter) {
0328 if (cand.superClusterRef() == gamIter->superCluster()) {
0329 reco::PhotonRef PhotRef(reco::PhotonRef(egPhotons, iegphot));
0330 egPhotonRef_.push_back(PhotRef);
0331 }
0332 iegphot++;
0333 }
0334
0335
0336
0337 for (unsigned iele = 0; iele < cand.elementsInBlocks().size(); ++iele) {
0338
0339 reco::PFBlockRef blockRef = cand.elementsInBlocks()[iele].first;
0340
0341 unsigned elementIndex = cand.elementsInBlocks()[iele].second;
0342
0343 if (blockRef.isNull())
0344 continue;
0345
0346
0347 const edm::OwnVector<reco::PFBlockElement> &elements = (*blockRef).elements();
0348
0349 const reco::PFBlockElement &pfbe(elements[elementIndex]);
0350
0351 if (pfbe.type() == reco::PFBlockElement::ECAL) {
0352
0353
0354
0355
0356 createBasicCluster(pfbe, basicClusters_[iphot], pfClusters_[iphot], correspondingDaughterCandidate(cand, pfbe));
0357 }
0358 if (pfbe.type() == reco::PFBlockElement::PS1) {
0359
0360 createPreshowerCluster(pfbe, preshowerClusters_[iphot], 1);
0361 }
0362 if (pfbe.type() == reco::PFBlockElement::PS2) {
0363
0364 createPreshowerCluster(pfbe, preshowerClusters_[iphot], 2);
0365 }
0366
0367 }
0368
0369
0370 basicClusters_p->insert(basicClusters_p->end(), basicClusters_[iphot].begin(), basicClusters_[iphot].end());
0371
0372 psClusters_p->insert(psClusters_p->end(), preshowerClusters_[iphot].begin(), preshowerClusters_[iphot].end());
0373
0374 ++iphot;
0375
0376 }
0377
0378
0379
0380 const edm::OrphanHandle<reco::BasicClusterCollection> bcRefProd =
0381 iEvent.put(std::move(basicClusters_p), PFBasicClusterCollection_);
0382
0383
0384 const edm::OrphanHandle<reco::PreshowerClusterCollection> psRefProd =
0385 iEvent.put(std::move(psClusters_p), PFPreshowerClusterCollection_);
0386
0387
0388 createBasicClusterPtrs(bcRefProd);
0389
0390 createPreshowerClusterPtrs(psRefProd);
0391
0392
0393
0394 if (status)
0395 createSuperClusters(*pfCandidates, outputSuperClusterCollection);
0396
0397
0398
0399
0400 auto superClusters_p = std::make_unique<reco::SuperClusterCollection>(outputSuperClusterCollection);
0401 const edm::OrphanHandle<reco::SuperClusterCollection> scRefProd =
0402 iEvent.put(std::move(superClusters_p), PFSuperClusterCollection_);
0403
0404
0405
0406
0407
0408
0409
0410
0411
0412
0413
0414 if (status)
0415 createOneLegConversions(scRefProd, outputOneLegConversionCollection);
0416
0417 auto SingleLeg_p = std::make_unique<reco::ConversionCollection>(outputOneLegConversionCollection);
0418 const edm::OrphanHandle<reco::ConversionCollection> ConvRefProd =
0419 iEvent.put(std::move(SingleLeg_p), PFConversionCollection_);
0420
0421
0422
0423
0424
0425
0426
0427
0428
0429
0430
0431
0432
0433
0434
0435
0436 if (status)
0437 createPhotonCores(scRefProd, ConvRefProd, outputPhotonCoreCollection);
0438
0439
0440
0441
0442 auto photonCores_p = std::make_unique<reco::PhotonCoreCollection>(outputPhotonCoreCollection);
0443
0444 const edm::OrphanHandle<reco::PhotonCoreCollection> pcRefProd =
0445 iEvent.put(std::move(photonCores_p), PFPhotonCoreCollection_);
0446
0447
0448
0449
0450
0451
0452
0453
0454
0455
0456
0457
0458
0459
0460
0461
0462
0463
0464
0465
0466
0467
0468
0469
0470
0471
0472
0473
0474
0475
0476
0477 reco::VertexCollection vertexCollection;
0478 bool validVertex = true;
0479 iEvent.getByLabel(vertexProducer_, vertexHandle);
0480 if (!vertexHandle.isValid()) {
0481 edm::LogWarning("PhotonProducer") << "Error! Can't get the product primary Vertex Collection "
0482 << "\n";
0483 validVertex = false;
0484 }
0485 if (validVertex)
0486 vertexCollection = *(vertexHandle.product());
0487
0488
0489
0490
0491
0492
0493
0494
0495
0496
0497
0498
0499
0500
0501
0502
0503
0504
0505
0506
0507
0508
0509
0510
0511
0512
0513
0514
0515
0516
0517
0518
0519
0520
0521
0522
0523 if (status)
0524 createPhotons(vertexCollection, egPhotons, pcRefProd, isolationValues, outputPhotonCollection);
0525
0526
0527 auto photons_p = std::make_unique<reco::PhotonCollection>(outputPhotonCollection);
0528
0529 const edm::OrphanHandle<reco::PhotonCollection> photonRefProd = iEvent.put(std::move(photons_p), PFPhotonCollection_);
0530
0531
0532
0533
0534
0535
0536
0537
0538
0539
0540
0541
0542
0543
0544
0545
0546
0547
0548
0549
0550
0551
0552
0553
0554
0555
0556
0557
0558
0559
0560
0561
0562
0563
0564
0565
0566
0567
0568
0569
0570
0571
0572
0573
0574
0575
0576
0577
0578
0579 }
0580
0581 bool PFPhotonTranslator::fetchCandidateCollection(edm::Handle<reco::PFCandidateCollection> &c,
0582 const edm::InputTag &tag,
0583 const edm::Event &iEvent) const {
0584 bool found = iEvent.getByLabel(tag, c);
0585
0586 if (!found && !emptyIsOk_) {
0587 std::ostringstream err;
0588 err << " cannot get PFCandidates: " << tag << std::endl;
0589 edm::LogError("PFPhotonTranslator") << err.str();
0590 }
0591 return found;
0592 }
0593
0594
0595
0596
0597 void PFPhotonTranslator::createBasicCluster(const reco::PFBlockElement &PFBE,
0598 reco::BasicClusterCollection &basicClusters,
0599 std::vector<const reco::PFCluster *> &pfClusters,
0600 const reco::PFCandidate &coCandidate) const {
0601 const reco::PFClusterRef &myPFClusterRef = PFBE.clusterRef();
0602 if (myPFClusterRef.isNull())
0603 return;
0604
0605 const reco::PFCluster &myPFCluster(*myPFClusterRef);
0606 pfClusters.push_back(&myPFCluster);
0607
0608
0609
0610
0611 basicClusters.push_back(reco::CaloCluster(
0612 myPFCluster.energy(),
0613 myPFCluster.position(),
0614 myPFCluster.caloID(),
0615 myPFCluster.hitsAndFractions(),
0616 myPFCluster.algo(),
0617 myPFCluster.seed()));
0618 }
0619
0620 void PFPhotonTranslator::createPreshowerCluster(const reco::PFBlockElement &PFBE,
0621 reco::PreshowerClusterCollection &preshowerClusters,
0622 unsigned plane) const {
0623 const reco::PFClusterRef &myPFClusterRef = PFBE.clusterRef();
0624 preshowerClusters.push_back(reco::PreshowerCluster(
0625 myPFClusterRef->energy(), myPFClusterRef->position(), myPFClusterRef->hitsAndFractions(), plane));
0626 }
0627
0628 void PFPhotonTranslator::createBasicClusterPtrs(
0629 const edm::OrphanHandle<reco::BasicClusterCollection> &basicClustersHandle) {
0630 unsigned size = photPFCandidateIndex_.size();
0631 unsigned basicClusterCounter = 0;
0632 basicClusterPtr_.resize(size);
0633
0634 for (unsigned iphot = 0; iphot < size; ++iphot)
0635 {
0636 unsigned nbc = basicClusters_[iphot].size();
0637 for (unsigned ibc = 0; ibc < nbc; ++ibc)
0638 {
0639
0640 reco::CaloClusterPtr bcPtr(basicClustersHandle, basicClusterCounter);
0641 basicClusterPtr_[iphot].push_back(bcPtr);
0642 ++basicClusterCounter;
0643 }
0644 }
0645 }
0646
0647 void PFPhotonTranslator::createPreshowerClusterPtrs(
0648 const edm::OrphanHandle<reco::PreshowerClusterCollection> &preshowerClustersHandle) {
0649 unsigned size = photPFCandidateIndex_.size();
0650 unsigned psClusterCounter = 0;
0651 preshowerClusterPtr_.resize(size);
0652
0653 for (unsigned iphot = 0; iphot < size; ++iphot)
0654 {
0655 unsigned nbc = preshowerClusters_[iphot].size();
0656 for (unsigned ibc = 0; ibc < nbc; ++ibc)
0657 {
0658
0659 reco::CaloClusterPtr psPtr(preshowerClustersHandle, psClusterCounter);
0660 preshowerClusterPtr_[iphot].push_back(psPtr);
0661 ++psClusterCounter;
0662 }
0663 }
0664 }
0665
0666 void PFPhotonTranslator::createSuperClusters(const reco::PFCandidateCollection &pfCand,
0667 reco::SuperClusterCollection &superClusters) const {
0668 unsigned nphot = photPFCandidateIndex_.size();
0669 for (unsigned iphot = 0; iphot < nphot; ++iphot) {
0670
0671
0672
0673 double sclusterE = 0;
0674 double posX = 0.;
0675 double posY = 0.;
0676 double posZ = 0.;
0677
0678 unsigned nbasics = basicClusters_[iphot].size();
0679 for (unsigned ibc = 0; ibc < nbasics; ++ibc) {
0680
0681
0682 double e = basicClusters_[iphot][ibc].energy();
0683 sclusterE += e;
0684 posX += e * basicClusters_[iphot][ibc].position().X();
0685 posY += e * basicClusters_[iphot][ibc].position().Y();
0686 posZ += e * basicClusters_[iphot][ibc].position().Z();
0687 }
0688 posX /= sclusterE;
0689 posY /= sclusterE;
0690 posZ /= sclusterE;
0691
0692
0693
0694
0695
0696
0697
0698
0699
0700 PFClusterWidthAlgo pfwidth(pfClusters_[iphot]);
0701
0702 double correctedEnergy = pfCand[photPFCandidateIndex_[iphot]].ecalEnergy();
0703 reco::SuperCluster mySuperCluster(correctedEnergy, math::XYZPoint(posX, posY, posZ));
0704
0705 if (nbasics) {
0706
0707
0708
0709 mySuperCluster.setSeed(basicClusterPtr_[iphot][0]);
0710 } else {
0711
0712
0713
0714
0715
0716 mySuperCluster.setSeed(reco::CaloClusterPtr());
0717 }
0718
0719
0720 for (unsigned ibc = 0; ibc < nbasics; ++ibc) {
0721 mySuperCluster.addCluster(basicClusterPtr_[iphot][ibc]);
0722
0723 const std::vector<std::pair<DetId, float> > &v1 = basicClusters_[iphot][ibc].hitsAndFractions();
0724
0725 for (std::vector<std::pair<DetId, float> >::const_iterator diIt = v1.begin(); diIt != v1.end(); ++diIt) {
0726
0727 mySuperCluster.addHitAndFraction(diIt->first, diIt->second);
0728 }
0729 }
0730
0731 unsigned nps = preshowerClusterPtr_[iphot].size();
0732 for (unsigned ips = 0; ips < nps; ++ips) {
0733 mySuperCluster.addPreshowerCluster(preshowerClusterPtr_[iphot][ips]);
0734 }
0735
0736
0737 mySuperCluster.setPreshowerEnergy(pfCand[photPFCandidateIndex_[iphot]].pS1Energy() +
0738 pfCand[photPFCandidateIndex_[iphot]].pS2Energy());
0739
0740
0741 mySuperCluster.setEtaWidth(pfwidth.pflowEtaWidth());
0742 mySuperCluster.setPhiWidth(pfwidth.pflowPhiWidth());
0743
0744 mySuperCluster.rawEnergy();
0745
0746
0747
0748 superClusters.push_back(mySuperCluster);
0749
0750 }
0751 }
0752
0753 void PFPhotonTranslator::createOneLegConversions(
0754 const edm::OrphanHandle<reco::SuperClusterCollection> &superClustersHandle,
0755 reco::ConversionCollection &oneLegConversions) {
0756
0757
0758 unsigned nphot = photPFCandidateIndex_.size();
0759 for (unsigned iphot = 0; iphot < nphot; ++iphot) {
0760
0761
0762
0763 if (conv1legPFCandidateIndex_[iphot] > -1) {
0764 for (unsigned iConv = 0; iConv < pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]].size(); iConv++) {
0765 reco::CaloClusterPtrVector scPtrVec;
0766 std::vector<reco::CaloClusterPtr> matchingBC;
0767 math::Error<3>::type error;
0768 const reco::Vertex *convVtx =
0769 new reco::Vertex(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerPosition(), error);
0770
0771
0772
0773
0774 std::vector<reco::TrackRef> OneLegConvVector;
0775 OneLegConvVector.push_back(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]);
0776
0777 std::vector<reco::TrackRef> tr = OneLegConvVector;
0778 std::vector<math::XYZPointF> trackPositionAtEcalVec;
0779 std::vector<math::XYZPointF> innPointVec;
0780 std::vector<math::XYZVectorF> trackPinVec;
0781 std::vector<math::XYZVectorF> trackPoutVec;
0782 math::XYZPointF trackPositionAtEcal(
0783 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerPosition().X(),
0784 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerPosition().Y(),
0785 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerPosition().Z());
0786 math::XYZPointF innPoint(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerPosition().X(),
0787 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerPosition().Y(),
0788 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerPosition().Z());
0789 math::XYZVectorF trackPin(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerMomentum().X(),
0790 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerMomentum().Y(),
0791 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerMomentum().Z());
0792 math::XYZVectorF trackPout(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerMomentum().X(),
0793 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerMomentum().Y(),
0794 pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->outerMomentum().Z());
0795 float DCA = pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->d0();
0796 trackPositionAtEcalVec.push_back(trackPositionAtEcal);
0797 innPointVec.push_back(innPoint);
0798 trackPinVec.push_back(trackPin);
0799 trackPoutVec.push_back(trackPout);
0800 std::vector<float> OneLegMvaVector;
0801 reco::Conversion myOneLegConversion(scPtrVec,
0802 OneLegConvVector,
0803 trackPositionAtEcalVec,
0804 *convVtx,
0805 matchingBC,
0806 DCA,
0807 innPointVec,
0808 trackPinVec,
0809 trackPoutVec,
0810 pfSingleLegConvMva_[conv1legPFCandidateIndex_[iphot]][iConv],
0811 reco::Conversion::pflow);
0812 OneLegMvaVector.push_back(pfSingleLegConvMva_[conv1legPFCandidateIndex_[iphot]][iConv]);
0813 myOneLegConversion.setOneLegMVA(OneLegMvaVector);
0814
0815
0816
0817
0818
0819
0820
0821
0822
0823
0824
0825
0826
0827
0828
0829 oneLegConversions.push_back(myOneLegConversion);
0830
0831
0832 }
0833 }
0834 }
0835 }
0836
0837 void PFPhotonTranslator::createPhotonCores(const edm::OrphanHandle<reco::SuperClusterCollection> &superClustersHandle,
0838 const edm::OrphanHandle<reco::ConversionCollection> &oneLegConversionHandle,
0839 reco::PhotonCoreCollection &photonCores) {
0840
0841
0842 unsigned nphot = photPFCandidateIndex_.size();
0843
0844 unsigned i1legtot = 0;
0845
0846 for (unsigned iphot = 0; iphot < nphot; ++iphot) {
0847
0848
0849 reco::PhotonCore myPhotonCore;
0850
0851 reco::SuperClusterRef SCref(reco::SuperClusterRef(superClustersHandle, iphot));
0852
0853 myPhotonCore.setPFlowPhoton(true);
0854 myPhotonCore.setStandardPhoton(false);
0855 myPhotonCore.setParentSuperCluster(SCref);
0856 myPhotonCore.setSuperCluster(egSCRef_[iphot]);
0857
0858 reco::ElectronSeedRefVector pixelSeeds = egPhotonRef_[iphot]->electronPixelSeeds();
0859 for (unsigned iseed = 0; iseed < pixelSeeds.size(); iseed++) {
0860 myPhotonCore.addElectronPixelSeed(pixelSeeds[iseed]);
0861 }
0862
0863
0864
0865
0866
0867
0868 if (conv1legPFCandidateIndex_[iphot] > -1) {
0869 for (unsigned int iConv = 0; iConv < pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]].size(); iConv++) {
0870 const reco::ConversionRef &OneLegRef(reco::ConversionRef(oneLegConversionHandle, i1legtot));
0871 myPhotonCore.addOneLegConversion(OneLegRef);
0872
0873
0874
0875
0876
0877
0878
0879
0880
0881
0882
0883
0884
0885
0886
0887
0888
0889
0890 i1legtot++;
0891 }
0892 }
0893
0894 if (conv2legPFCandidateIndex_[iphot] > -1) {
0895 for (unsigned int iConv = 0; iConv < pfConv_[conv2legPFCandidateIndex_[iphot]].size(); iConv++) {
0896 const reco::ConversionRef &TwoLegRef(pfConv_[conv2legPFCandidateIndex_[iphot]][iConv]);
0897 myPhotonCore.addConversion(TwoLegRef);
0898 }
0899
0900
0901
0902
0903
0904
0905
0906
0907
0908
0909
0910
0911
0912
0913
0914
0915
0916 }
0917
0918 photonCores.push_back(myPhotonCore);
0919 }
0920
0921
0922 }
0923
0924 void PFPhotonTranslator::createPhotons(reco::VertexCollection &vertexCollection,
0925 edm::Handle<reco::PhotonCollection> &egPhotons,
0926 const edm::OrphanHandle<reco::PhotonCoreCollection> &photonCoresHandle,
0927 const IsolationValueMaps &isolationValues,
0928 reco::PhotonCollection &photons) {
0929
0930
0931 unsigned nphot = photPFCandidateIndex_.size();
0932
0933 for (unsigned iphot = 0; iphot < nphot; ++iphot) {
0934
0935
0936 reco::PhotonCoreRef PCref(reco::PhotonCoreRef(photonCoresHandle, iphot));
0937
0938 math::XYZPoint vtx(0., 0., 0.);
0939 if (!vertexCollection.empty())
0940 vtx = vertexCollection.begin()->position();
0941
0942
0943 math::XYZVector direction = PCref->parentSuperCluster()->position() - vtx;
0944
0945
0946 math::XYZVector P3 = direction.unit() * PCref->parentSuperCluster()->energy();
0947 LorentzVector P4(P3.x(), P3.y(), P3.z(), PCref->parentSuperCluster()->energy());
0948
0949 reco::Photon myPhoton(P4, PCref->parentSuperCluster()->position(), PCref, vtx);
0950
0951
0952 reco::Photon::ShowerShape showerShape;
0953 reco::Photon::SaturationInfo saturationInfo;
0954 reco::Photon::FiducialFlags fiducialFlags;
0955 reco::Photon::IsolationVariables isolationVariables03;
0956 reco::Photon::IsolationVariables isolationVariables04;
0957
0958 showerShape.e1x5 = egPhotonRef_[iphot]->e1x5();
0959 showerShape.e2x5 = egPhotonRef_[iphot]->e2x5();
0960 showerShape.e3x3 = egPhotonRef_[iphot]->e3x3();
0961 showerShape.e5x5 = egPhotonRef_[iphot]->e5x5();
0962 showerShape.maxEnergyXtal = egPhotonRef_[iphot]->maxEnergyXtal();
0963 showerShape.sigmaEtaEta = egPhotonRef_[iphot]->sigmaEtaEta();
0964 showerShape.sigmaIetaIeta = egPhotonRef_[iphot]->sigmaIetaIeta();
0965 for (uint id = 0; id < showerShape.hcalOverEcal.size(); ++id) {
0966 showerShape.hcalOverEcal[id] = egPhotonRef_[iphot]->hcalOverEcal(id + 1);
0967 }
0968 myPhoton.setShowerShapeVariables(showerShape);
0969
0970 saturationInfo.nSaturatedXtals = egPhotonRef_[iphot]->nSaturatedXtals();
0971 saturationInfo.isSeedSaturated = egPhotonRef_[iphot]->isSeedSaturated();
0972 myPhoton.setSaturationInfo(saturationInfo);
0973
0974 fiducialFlags.isEB = egPhotonRef_[iphot]->isEB();
0975 fiducialFlags.isEE = egPhotonRef_[iphot]->isEE();
0976 fiducialFlags.isEBEtaGap = egPhotonRef_[iphot]->isEBEtaGap();
0977 fiducialFlags.isEBPhiGap = egPhotonRef_[iphot]->isEBPhiGap();
0978 fiducialFlags.isEERingGap = egPhotonRef_[iphot]->isEERingGap();
0979 fiducialFlags.isEEDeeGap = egPhotonRef_[iphot]->isEEDeeGap();
0980 fiducialFlags.isEBEEGap = egPhotonRef_[iphot]->isEBEEGap();
0981 myPhoton.setFiducialVolumeFlags(fiducialFlags);
0982
0983 isolationVariables03.ecalRecHitSumEt = egPhotonRef_[iphot]->ecalRecHitSumEtConeDR03();
0984 for (uint id = 0; id < isolationVariables03.hcalRecHitSumEt.size(); ++id) {
0985 isolationVariables03.hcalRecHitSumEt[id] = egPhotonRef_[iphot]->hcalTowerSumEtConeDR03(id + 1);
0986 }
0987 isolationVariables03.trkSumPtSolidCone = egPhotonRef_[iphot]->trkSumPtSolidConeDR03();
0988 isolationVariables03.trkSumPtHollowCone = egPhotonRef_[iphot]->trkSumPtHollowConeDR03();
0989 isolationVariables03.nTrkSolidCone = egPhotonRef_[iphot]->nTrkSolidConeDR03();
0990 isolationVariables03.nTrkHollowCone = egPhotonRef_[iphot]->nTrkHollowConeDR03();
0991 isolationVariables04.ecalRecHitSumEt = egPhotonRef_[iphot]->ecalRecHitSumEtConeDR04();
0992 for (uint id = 0; id < isolationVariables04.hcalRecHitSumEt.size(); ++id) {
0993 isolationVariables04.hcalRecHitSumEt[id] = egPhotonRef_[iphot]->hcalTowerSumEtConeDR04(id + 1);
0994 }
0995 isolationVariables04.trkSumPtSolidCone = egPhotonRef_[iphot]->trkSumPtSolidConeDR04();
0996 isolationVariables04.trkSumPtHollowCone = egPhotonRef_[iphot]->trkSumPtHollowConeDR04();
0997 isolationVariables04.nTrkSolidCone = egPhotonRef_[iphot]->nTrkSolidConeDR04();
0998 isolationVariables04.nTrkHollowCone = egPhotonRef_[iphot]->nTrkHollowConeDR04();
0999 myPhoton.setIsolationVariables(isolationVariables04, isolationVariables03);
1000
1001 reco::Photon::PflowIsolationVariables myPFIso;
1002 myPFIso.chargedHadronIso = (*isolationValues[0])[CandidatePtr_[iphot]];
1003 myPFIso.photonIso = (*isolationValues[1])[CandidatePtr_[iphot]];
1004 myPFIso.neutralHadronIso = (*isolationValues[2])[CandidatePtr_[iphot]];
1005 myPhoton.setPflowIsolationVariables(myPFIso);
1006
1007 reco::Photon::PflowIDVariables myPFVariables;
1008
1009 reco::Mustache myMustache(mustacheSCParams_);
1010 myMustache.MustacheID(
1011 *(myPhoton.parentSuperCluster()), myPFVariables.nClusterOutsideMustache, myPFVariables.etOutsideMustache);
1012 myPFVariables.mva = pfPhotonMva_[iphot];
1013 myPhoton.setPflowIDVariables(myPFVariables);
1014
1015
1016
1017
1018 myPhoton.setCorrectedEnergy(
1019 reco::Photon::regression2, energyRegression_[iphot], energyRegressionError_[iphot], false);
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075 photons.push_back(myPhoton);
1076 }
1077
1078
1079 }
1080
1081 const reco::PFCandidate &PFPhotonTranslator::correspondingDaughterCandidate(const reco::PFCandidate &cand,
1082 const reco::PFBlockElement &pfbe) const {
1083 unsigned refindex = pfbe.index();
1084
1085 reco::PFCandidate::const_iterator myDaughterCandidate = cand.begin();
1086 reco::PFCandidate::const_iterator itend = cand.end();
1087
1088 for (; myDaughterCandidate != itend; ++myDaughterCandidate) {
1089 const reco::PFCandidate *myPFCandidate = (const reco::PFCandidate *)&*myDaughterCandidate;
1090 if (myPFCandidate->elementsInBlocks().size() != 1) {
1091
1092 return cand;
1093 }
1094 if (myPFCandidate->elementsInBlocks()[0].second == refindex) {
1095
1096 return *myPFCandidate;
1097 }
1098 }
1099 return cand;
1100 }