CMS 3D CMS Logo

PFPhotonTranslator.cc
Go to the documentation of this file.
6 
18 //#include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
22 
26 
27 //#include "Geometry/Records/interface/CaloGeometryRecord.h"
28 //#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
29 //#include "Geometry/CaloTopology/interface/CaloTopology.h"
30 //#include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
31 //#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
32 //#include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
33 //#include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
34 //#include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
35 //#include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
36 //#include "RecoEgamma/PhotonIdentification/interface/PhotonIsolationCalculator.h"
37 //#include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionFactory.h"
38 //#include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionBaseClass.h"
39 //#include "CondFormats/EcalObjects/interface/EcalFunctionParameters.h"
40 //#include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
42 
48 
49 #include <Math/VectorUtil.h>
50 #include <vector>
51 #include "TLorentzVector.h"
52 #include "TMath.h"
53 
54 using namespace edm;
55 using namespace std;
56 using namespace reco;
57 
58 using namespace ROOT::Math::VectorUtil;
62 
63 
65 
66  //std::cout << "PFPhotonTranslator" << std::endl;
67 
68  inputTagPFCandidates_
69  = iConfig.getParameter<edm::InputTag>("PFCandidate");
70 
71  edm::ParameterSet isoVals = iConfig.getParameter<edm::ParameterSet> ("isolationValues");
72  inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfChargedHadrons"));
73  inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfPhotons"));
74  inputTagIsoVals_.push_back(isoVals.getParameter<edm::InputTag>("pfNeutralHadrons"));
75 
76 
77  PFBasicClusterCollection_ = iConfig.getParameter<std::string>("PFBasicClusters");
78  PFPreshowerClusterCollection_ = iConfig.getParameter<std::string>("PFPreshowerClusters");
79  PFSuperClusterCollection_ = iConfig.getParameter<std::string>("PFSuperClusters");
80  PFConversionCollection_ = iConfig.getParameter<std::string>("PFConversionCollection");
81  PFPhotonCoreCollection_ = iConfig.getParameter<std::string>("PFPhotonCores");
82  PFPhotonCollection_ = iConfig.getParameter<std::string>("PFPhotons");
83 
84  EGPhotonCollection_ = iConfig.getParameter<std::string>("EGPhotons");
85 
86  vertexProducer_ = iConfig.getParameter<std::string>("primaryVertexProducer");
87 
88  barrelEcalHits_ = iConfig.getParameter<edm::InputTag>("barrelEcalHits");
89  endcapEcalHits_ = iConfig.getParameter<edm::InputTag>("endcapEcalHits");
90 
91  hcalTowers_ = iConfig.getParameter<edm::InputTag>("hcalTowers");
92  hOverEConeSize_ = iConfig.getParameter<double>("hOverEConeSize");
93 
94  if (iConfig.exists("emptyIsOk")) emptyIsOk_ = iConfig.getParameter<bool>("emptyIsOk");
95  else emptyIsOk_=false;
96 
97  produces<reco::BasicClusterCollection>(PFBasicClusterCollection_);
98  produces<reco::PreshowerClusterCollection>(PFPreshowerClusterCollection_);
99  produces<reco::SuperClusterCollection>(PFSuperClusterCollection_);
100  produces<reco::PhotonCoreCollection>(PFPhotonCoreCollection_);
101  produces<reco::PhotonCollection>(PFPhotonCollection_);
102  produces<reco::ConversionCollection>(PFConversionCollection_);
103 }
104 
106 
108  const edm::EventSetup& iSetup) {
109 
110  //cout << "NEW EVENT"<<endl;
111 
112  auto basicClusters_p = std::make_unique<reco::BasicClusterCollection>();
113 
114  auto psClusters_p = std::make_unique<reco::PreshowerClusterCollection>();
115 
116  /*
117  auto SingleLeg_p = std::make_unique<reco::ConversionCollection>();
118  */
119 
120  reco::SuperClusterCollection outputSuperClusterCollection;
121  reco::ConversionCollection outputOneLegConversionCollection;
122  reco::PhotonCoreCollection outputPhotonCoreCollection;
124 
125  outputSuperClusterCollection.clear();
126  outputOneLegConversionCollection.clear();
127  outputPhotonCoreCollection.clear();
128  outputPhotonCollection.clear();
129 
130 
132  bool status=fetchCandidateCollection(pfCandidates,
133  inputTagPFCandidates_,
134  iEvent );
135 
137  iEvent.getByLabel(EGPhotonCollection_, egPhotons);
138 
139 
140  Handle<reco::VertexCollection> vertexHandle;
141 
142 
143  IsolationValueMaps isolationValues(inputTagIsoVals_.size());
144  for (size_t j = 0; j<inputTagIsoVals_.size(); ++j) {
145  iEvent.getByLabel(inputTagIsoVals_[j], isolationValues[j]);
146  }
147 
148 
149  // clear the vectors
150  photPFCandidateIndex_.clear();
151  basicClusters_.clear();
152  pfClusters_.clear();
153  preshowerClusters_.clear();
154  superClusters_.clear();
155  basicClusterPtr_.clear();
156  preshowerClusterPtr_.clear();
157  CandidatePtr_.clear();
158  egSCRef_.clear();
159  egPhotonRef_.clear();
160  pfPhotonMva_.clear();
161  energyRegression_.clear();
162  energyRegressionError_.clear();
163  pfConv_.clear();
164  pfSingleLegConv_.clear();
165  pfSingleLegConvMva_.clear();
166  conv1legPFCandidateIndex_.clear();
167  conv2legPFCandidateIndex_.clear();
168 
169  // loop on the candidates
170  //CC@@
171  // we need first to create AND put the SuperCluster,
172  // basic clusters and presh clusters collection
173  // in order to get a working Handle
174  unsigned ncand=(status)?pfCandidates->size():0;
175 
176  unsigned iphot=0;
177  unsigned iconv1leg=0;
178  unsigned iconv2leg=0;
179 
180  for( unsigned i=0; i<ncand; ++i ) {
181 
182  const reco::PFCandidate& cand = (*pfCandidates)[i];
183  if(cand.particleId()!=reco::PFCandidate::gamma) continue;
184  //cout << "cand.mva_nothing_gamma()="<<cand. mva_nothing_gamma()<<endl;
185  if(cand. mva_nothing_gamma()>0.001)//Found PFPhoton with PFPhoton Extras saved
186  {
187 
188  //cout << "NEW PHOTON" << endl;
189 
190  //std::cout << "nDoubleLegConv="<<cand.photonExtraRef()->conversionRef().size()<<std::endl;
191 
192  if (!cand.photonExtraRef()->conversionRef().empty()){
193 
194  pfConv_.push_back(reco::ConversionRefVector());
195 
196  const reco::ConversionRefVector & doubleLegConvColl = cand.photonExtraRef()->conversionRef();
197  for (unsigned int iconv=0; iconv<doubleLegConvColl.size(); iconv++){
198  pfConv_[iconv2leg].push_back(doubleLegConvColl[iconv]);
199  }
200 
201  conv2legPFCandidateIndex_.push_back(iconv2leg);
202  iconv2leg++;
203  }
204  else conv2legPFCandidateIndex_.push_back(-1);
205 
206  const std::vector<reco::TrackRef> & singleLegConvColl = cand.photonExtraRef()->singleLegConvTrackRef();
207  const std::vector<float>& singleLegConvCollMva = cand.photonExtraRef()->singleLegConvMva();
208 
209  //std::cout << "nSingleLegConv=" <<singleLegConvColl.size() << std::endl;
210 
211  if (!singleLegConvColl.empty()){
212 
213  pfSingleLegConv_.push_back(std::vector<reco::TrackRef>());
214  pfSingleLegConvMva_.push_back(std::vector<float>());
215 
216  //cout << "nTracks="<< singleLegConvColl.size()<<endl;
217  for (unsigned int itk=0; itk<singleLegConvColl.size(); itk++){
218  //cout << "Track pt="<< singleLegConvColl[itk]->pt() <<endl;
219 
220  pfSingleLegConv_[iconv1leg].push_back(singleLegConvColl[itk]);
221  pfSingleLegConvMva_[iconv1leg].push_back(singleLegConvCollMva[itk]);
222  }
223 
224 
225  conv1legPFCandidateIndex_.push_back(iconv1leg);
226 
227  iconv1leg++;
228  }
229  else conv1legPFCandidateIndex_.push_back(-1);
230 
231  }
232 
233  photPFCandidateIndex_.push_back(i);
234  pfPhotonMva_.push_back(cand.mva_nothing_gamma());
235  energyRegression_.push_back(cand.photonExtraRef()->MVAGlobalCorrE());
236  energyRegressionError_.push_back(cand.photonExtraRef()->MVAGlobalCorrEError());
237  basicClusters_.push_back(reco::BasicClusterCollection());
238  pfClusters_.push_back(std::vector<const reco::PFCluster *>());
239  preshowerClusters_.push_back(reco::PreshowerClusterCollection());
240  superClusters_.push_back(reco::SuperClusterCollection());
241 
242  reco::PFCandidatePtr ptrToPFPhoton(pfCandidates,i);
243  CandidatePtr_.push_back(ptrToPFPhoton);
244  egSCRef_.push_back(cand.superClusterRef());
245  //std::cout << "PFPhoton cand " << iphot << std::endl;
246 
247  int iegphot=0;
248  for (reco::PhotonCollection::const_iterator gamIter = egPhotons->begin(); gamIter != egPhotons->end(); ++gamIter){
249  if (cand.superClusterRef()==gamIter->superCluster()){
250  reco::PhotonRef PhotRef(reco::PhotonRef(egPhotons, iegphot));
251  egPhotonRef_.push_back(PhotRef);
252  }
253  iegphot++;
254  }
255 
256 
257  //std::cout << "Cand elements in blocks : " << cand.elementsInBlocks().size() << std::endl;
258 
259  for(unsigned iele=0; iele<cand.elementsInBlocks().size(); ++iele) {
260  // first get the block
261  reco::PFBlockRef blockRef = cand.elementsInBlocks()[iele].first;
262  //
263  unsigned elementIndex = cand.elementsInBlocks()[iele].second;
264  // check it actually exists
265  if(blockRef.isNull()) continue;
266 
267  // then get the elements of the block
268  const edm::OwnVector< reco::PFBlockElement >& elements = (*blockRef).elements();
269 
270  const reco::PFBlockElement & pfbe (elements[elementIndex]);
271  // The first ECAL element should be the cluster associated to the GSF; defined as the seed
272  if(pfbe.type()==reco::PFBlockElement::ECAL)
273  {
274 
275  //std::cout << "BlockElement ECAL" << std::endl;
276  // the Brem photons are saved as daughter PFCandidate; this
277  // is convenient to access the corrected energy
278  // std::cout << " Found candidate " << correspondingDaughterCandidate(coCandidate,pfbe) << " " << coCandidate << std::endl;
279  createBasicCluster(pfbe,basicClusters_[iphot],pfClusters_[iphot],correspondingDaughterCandidate(cand,pfbe));
280  }
281  if(pfbe.type()==reco::PFBlockElement::PS1)
282  {
283  //std::cout << "BlockElement PS1" << std::endl;
284  createPreshowerCluster(pfbe,preshowerClusters_[iphot],1);
285  }
286  if(pfbe.type()==reco::PFBlockElement::PS2)
287  {
288  //std::cout << "BlockElement PS2" << std::endl;
289  createPreshowerCluster(pfbe,preshowerClusters_[iphot],2);
290  }
291 
292 
293  } // loop on the elements
294 
295  // save the basic clusters
296  basicClusters_p->insert(basicClusters_p->end(),basicClusters_[iphot].begin(), basicClusters_[iphot].end());
297  // save the preshower clusters
298  psClusters_p->insert(psClusters_p->end(),preshowerClusters_[iphot].begin(),preshowerClusters_[iphot].end());
299 
300  ++iphot;
301 
302  } // loop on PFCandidates
303 
304 
305  //Save the basic clusters and get an handle as to be able to create valid Refs (thanks to Claude)
306  // std::cout << " Number of basic clusters " << basicClusters_p->size() << std::endl;
308  iEvent.put(std::move(basicClusters_p),PFBasicClusterCollection_);
309 
310  //preshower clusters
312  iEvent.put(std::move(psClusters_p),PFPreshowerClusterCollection_);
313 
314  // now that the Basic clusters are in the event, the Ref can be created
315  createBasicClusterPtrs(bcRefProd);
316  // now that the preshower clusters are in the event, the Ref can be created
317  createPreshowerClusterPtrs(psRefProd);
318 
319  // and now the Super cluster can be created with valid references
320  //if(status) createSuperClusters(*pfCandidates,*superClusters_p);
321  if(status) createSuperClusters(*pfCandidates,outputSuperClusterCollection);
322 
323  //std::cout << "nb superclusters in collection : "<<outputSuperClusterCollection.size()<<std::endl;
324 
325  // Let's put the super clusters in the event
326  auto superClusters_p = std::make_unique<reco::SuperClusterCollection>(outputSuperClusterCollection);
327  const edm::OrphanHandle<reco::SuperClusterCollection> scRefProd = iEvent.put(std::move(superClusters_p),PFSuperClusterCollection_);
328 
329 
330  /*
331  int ipho=0;
332  for (reco::SuperClusterCollection::const_iterator gamIter = scRefProd->begin(); gamIter != scRefProd->end(); ++gamIter){
333  std::cout << "SC i="<<ipho<<" energy="<<gamIter->energy()<<std::endl;
334  ipho++;
335  }
336  */
337 
338 
339  //1-leg conversions
340 
341 
342  if (status) createOneLegConversions(scRefProd, outputOneLegConversionCollection);
343 
344 
345  auto SingleLeg_p = std::make_unique<reco::ConversionCollection>(outputOneLegConversionCollection);
346  const edm::OrphanHandle<reco::ConversionCollection> ConvRefProd = iEvent.put(std::move(SingleLeg_p),PFConversionCollection_);
347  /*
348  int iconv = 0;
349  for (reco::ConversionCollection::const_iterator convIter = ConvRefProd->begin(); convIter != ConvRefProd->end(); ++convIter){
350 
351  std::cout << "OneLegConv i="<<iconv<<" nTracks="<<convIter->nTracks()<<" EoverP="<<convIter->EoverP() <<std::endl;
352  std::vector<edm::RefToBase<reco::Track> > convtracks = convIter->tracks();
353  for (unsigned int itk=0; itk<convtracks.size(); itk++){
354  std::cout << "Track pt="<< convtracks[itk]->pt() << std::endl;
355  }
356 
357  iconv++;
358  }
359  */
360 
361  //create photon cores
362  //if(status) createPhotonCores(pfCandidates, scRefProd, *photonCores_p);
363  if(status) createPhotonCores(scRefProd, ConvRefProd, outputPhotonCoreCollection);
364 
365  //std::cout << "nb photoncores in collection : "<<outputPhotonCoreCollection.size()<<std::endl;
366 
367  // Put the photon cores in the event
368  auto photonCores_p = std::make_unique<reco::PhotonCoreCollection>(outputPhotonCoreCollection);
369  //std::cout << "photon core collection put in unique_ptr"<<std::endl;
370  const edm::OrphanHandle<reco::PhotonCoreCollection> pcRefProd = iEvent.put(std::move(photonCores_p),PFPhotonCoreCollection_);
371 
372  //std::cout << "photon core have been put in the event"<<std::endl;
373  /*
374  int ipho=0;
375  for (reco::PhotonCoreCollection::const_iterator gamIter = pcRefProd->begin(); gamIter != pcRefProd->end(); ++gamIter){
376  std::cout << "PhotonCore i="<<ipho<<" energy="<<gamIter->parentSuperCluster()->energy()<<std::endl;
377  //for (unsigned int i=0; i<)
378 
379  std::cout << "PhotonCore i="<<ipho<<" nconv2leg="<<gamIter->conversions().size()<<" nconv1leg="<<gamIter->conversionsOneLeg().size()<<std::endl;
380 
381  const reco::ConversionRefVector & conv = gamIter->conversions();
382  for (unsigned int iconv=0; iconv<conv.size(); iconv++){
383  cout << "2-leg iconv="<<iconv<<endl;
384  cout << "2-leg nTracks="<<conv[iconv]->nTracks()<<endl;
385  cout << "2-leg EoverP="<<conv[iconv]->EoverP()<<endl;
386  cout << "2-leg ConvAlgorithm="<<conv[iconv]->algo()<<endl;
387  }
388 
389  const reco::ConversionRefVector & convbis = gamIter->conversionsOneLeg();
390  for (unsigned int iconv=0; iconv<convbis.size(); iconv++){
391  cout << "1-leg iconv="<<iconv<<endl;
392  cout << "1-leg nTracks="<<convbis[iconv]->nTracks()<<endl;
393  cout << "1-leg EoverP="<<convbis[iconv]->EoverP()<<endl;
394  cout << "1-leg ConvAlgorithm="<<convbis[iconv]->algo()<<endl;
395  }
396 
397  ipho++;
398  }
399  */
400 
401  //load vertices
403  bool validVertex=true;
404  iEvent.getByLabel(vertexProducer_, vertexHandle);
405  if (!vertexHandle.isValid()) {
406  edm::LogWarning("PhotonProducer") << "Error! Can't get the product primary Vertex Collection "<< "\n";
407  validVertex=false;
408  }
409  if (validVertex) vertexCollection = *(vertexHandle.product());
410 
411  /*
412  //load Ecal rechits
413  bool validEcalRecHits=true;
414  Handle<EcalRecHitCollection> barrelHitHandle;
415  EcalRecHitCollection barrelRecHits;
416  iEvent.getByLabel(barrelEcalHits_, barrelHitHandle);
417  if (!barrelHitHandle.isValid()) {
418  edm::LogError("PhotonProducer") << "Error! Can't get the product "<<barrelEcalHits_.label();
419  validEcalRecHits=false;
420  }
421  if ( validEcalRecHits) barrelRecHits = *(barrelHitHandle.product());
422 
423  Handle<EcalRecHitCollection> endcapHitHandle;
424  iEvent.getByLabel(endcapEcalHits_, endcapHitHandle);
425  EcalRecHitCollection endcapRecHits;
426  if (!endcapHitHandle.isValid()) {
427  edm::LogError("PhotonProducer") << "Error! Can't get the product "<<endcapEcalHits_.label();
428  validEcalRecHits=false;
429  }
430  if( validEcalRecHits) endcapRecHits = *(endcapHitHandle.product());
431 
432  //load detector topology & geometry
433  iSetup.get<CaloGeometryRecord>().get(theCaloGeom_);
434 
435  edm::ESHandle<CaloTopology> pTopology;
436  iSetup.get<CaloTopologyRecord>().get(theCaloTopo_);
437  const CaloTopology *topology = theCaloTopo_.product();
438 
439  // get Hcal towers collection
440  Handle<CaloTowerCollection> hcalTowersHandle;
441  iEvent.getByLabel(hcalTowers_, hcalTowersHandle);
442  */
443 
444  //create photon collection
445  //if(status) createPhotons(vertexCollection, pcRefProd, topology, &barrelRecHits, &endcapRecHits, hcalTowersHandle, isolationValues, outputPhotonCollection);
446  if(status) createPhotons(vertexCollection, egPhotons, pcRefProd, isolationValues, outputPhotonCollection);
447 
448  // Put the photons in the event
449  auto photons_p = std::make_unique<reco::PhotonCollection>(outputPhotonCollection);
450  //std::cout << "photon collection put in unique_ptr"<<std::endl;
451  const edm::OrphanHandle<reco::PhotonCollection> photonRefProd = iEvent.put(std::move(photons_p),PFPhotonCollection_);
452  //std::cout << "photons have been put in the event"<<std::endl;
453 
454  /*
455  ipho=0;
456  for (reco::PhotonCollection::const_iterator gamIter = photonRefProd->begin(); gamIter != photonRefProd->end(); ++gamIter){
457  std::cout << "Photon i="<<ipho<<" pfEnergy="<<gamIter->parentSuperCluster()->energy()<<std::endl;
458 
459  const reco::ConversionRefVector & conv = gamIter->conversions();
460  cout << "conversions obtained : conv.size()="<< conv.size()<<endl;
461  for (unsigned int iconv=0; iconv<conv.size(); iconv++){
462  cout << "iconv="<<iconv<<endl;
463  cout << "nTracks="<<conv[iconv]->nTracks()<<endl;
464  cout << "EoverP="<<conv[iconv]->EoverP()<<endl;
465 
466  cout << "Vtx x="<<conv[iconv]->conversionVertex().x() << " y="<< conv[iconv]->conversionVertex().y()<<" z="<<conv[iconv]->conversionVertex().z()<< endl;
467  cout << "VtxError x=" << conv[iconv]->conversionVertex().xError() << endl;
468 
469  std::vector<edm::RefToBase<reco::Track> > convtracks = conv[iconv]->tracks();
470  //cout << "nTracks="<< convtracks.size()<<endl;
471  for (unsigned int itk=0; itk<convtracks.size(); itk++){
472  double convtrackpt = convtracks[itk]->pt();
473  const edm::RefToBase<reco::Track> & mytrack = convtracks[itk];
474  cout << "Track pt="<< convtrackpt <<endl;
475  cout << "Track origin="<<gamIter->conversionTrackOrigin(mytrack)<<endl;
476  }
477  }
478 
479  //1-leg
480  const reco::ConversionRefVector & convbis = gamIter->conversionsOneLeg();
481  cout << "conversions obtained : conv.size()="<< convbis.size()<<endl;
482  for (unsigned int iconv=0; iconv<convbis.size(); iconv++){
483  cout << "iconv="<<iconv<<endl;
484  cout << "nTracks="<<convbis[iconv]->nTracks()<<endl;
485  cout << "EoverP="<<convbis[iconv]->EoverP()<<endl;
486 
487  cout << "Vtx x="<<convbis[iconv]->conversionVertex().x() << " y="<< convbis[iconv]->conversionVertex().y()<<" z="<<convbis[iconv]->conversionVertex().z()<< endl;
488  cout << "VtxError x=" << convbis[iconv]->conversionVertex().xError() << endl;
489 
490  std::vector<edm::RefToBase<reco::Track> > convtracks = convbis[iconv]->tracks();
491  //cout << "nTracks="<< convtracks.size()<<endl;
492  for (unsigned int itk=0; itk<convtracks.size(); itk++){
493  double convtrackpt = convtracks[itk]->pt();
494  cout << "Track pt="<< convtrackpt <<endl;
495  cout << "Track origin="<<gamIter->conversionTrackOrigin((convtracks[itk]))<<endl;
496  }
497  }
498  ipho++;
499  }
500  */
501 
502 }
503 
505  const edm::InputTag& tag,
506  const edm::Event& iEvent) const {
507  bool found = iEvent.getByLabel(tag, c);
508 
509  if(!found && !emptyIsOk_)
510  {
511  std::ostringstream err;
512  err<<" cannot get PFCandidates: "
513  <<tag<<std::endl;
514  edm::LogError("PFPhotonTranslator")<<err.str();
515  }
516  return found;
517 
518 }
519 
520 // The basic cluster is a copy of the PFCluster -> the energy is not corrected
521 // It should be possible to get the corrected energy (including the associated PS energy)
522 // from the PFCandidate daugthers ; Needs some work
524  reco::BasicClusterCollection & basicClusters,
525  std::vector<const reco::PFCluster *> & pfClusters,
526  const reco::PFCandidate & coCandidate) const
527 {
528  const reco::PFClusterRef& myPFClusterRef = PFBE.clusterRef();
529  if(myPFClusterRef.isNull()) return;
530 
531  const reco::PFCluster & myPFCluster (*myPFClusterRef);
532  pfClusters.push_back(&myPFCluster);
533  //std::cout << " Creating BC " << myPFCluster.energy() << " " << coCandidate.ecalEnergy() <<" "<< coCandidate.rawEcalEnergy() <<std::endl;
534  //std::cout << " # hits " << myPFCluster.hitsAndFractions().size() << std::endl;
535 
536 // basicClusters.push_back(reco::CaloCluster(myPFCluster.energy(),
537  basicClusters.push_back(reco::CaloCluster(//coCandidate.rawEcalEnergy(),
538  myPFCluster.energy(),
539  myPFCluster.position(),
540  myPFCluster.caloID(),
541  myPFCluster.hitsAndFractions(),
542  myPFCluster.algo(),
543  myPFCluster.seed()));
544 }
545 
547 {
548  const reco::PFClusterRef& myPFClusterRef= PFBE.clusterRef();
549  preshowerClusters.push_back(reco::PreshowerCluster(myPFClusterRef->energy(),myPFClusterRef->position(),
550  myPFClusterRef->hitsAndFractions(),plane));
551 }
552 
554 {
555  unsigned size=photPFCandidateIndex_.size();
556  unsigned basicClusterCounter=0;
557  basicClusterPtr_.resize(size);
558 
559  for(unsigned iphot=0;iphot<size;++iphot) // loop on tracks
560  {
561  unsigned nbc=basicClusters_[iphot].size();
562  for(unsigned ibc=0;ibc<nbc;++ibc) // loop on basic clusters
563  {
564  // std::cout << "Track "<< iGSF << " ref " << basicClusterCounter << std::endl;
565  reco::CaloClusterPtr bcPtr(basicClustersHandle,basicClusterCounter);
566  basicClusterPtr_[iphot].push_back(bcPtr);
567  ++basicClusterCounter;
568  }
569  }
570 }
571 
573 {
574  unsigned size=photPFCandidateIndex_.size();
575  unsigned psClusterCounter=0;
576  preshowerClusterPtr_.resize(size);
577 
578  for(unsigned iphot=0;iphot<size;++iphot) // loop on tracks
579  {
580  unsigned nbc=preshowerClusters_[iphot].size();
581  for(unsigned ibc=0;ibc<nbc;++ibc) // loop on basic clusters
582  {
583  // std::cout << "Track "<< iGSF << " ref " << basicClusterCounter << std::endl;
584  reco::CaloClusterPtr psPtr(preshowerClustersHandle,psClusterCounter);
585  preshowerClusterPtr_[iphot].push_back(psPtr);
586  ++psClusterCounter;
587  }
588  }
589 }
590 
592  reco::SuperClusterCollection &superClusters) const
593 {
594  unsigned nphot=photPFCandidateIndex_.size();
595  for(unsigned iphot=0;iphot<nphot;++iphot)
596  {
597 
598  //cout << "SC iphot=" << iphot << endl;
599 
600  // Computes energy position a la e/gamma
601  double sclusterE=0;
602  double posX=0.;
603  double posY=0.;
604  double posZ=0.;
605 
606  unsigned nbasics=basicClusters_[iphot].size();
607  for(unsigned ibc=0;ibc<nbasics;++ibc)
608  {
609  //cout << "BC in SC : iphot="<<iphot<<endl;
610 
611  double e = basicClusters_[iphot][ibc].energy();
612  sclusterE += e;
613  posX += e * basicClusters_[iphot][ibc].position().X();
614  posY += e * basicClusters_[iphot][ibc].position().Y();
615  posZ += e * basicClusters_[iphot][ibc].position().Z();
616  }
617  posX /=sclusterE;
618  posY /=sclusterE;
619  posZ /=sclusterE;
620 
621  /*
622  if(pfCand[gsfPFCandidateIndex_[iphot]].gsfTrackRef()!=GsfTrackRef_[iphot])
623  {
624  edm::LogError("PFElectronTranslator") << " Major problem in PFElectron Translator" << std::endl;
625  }
626  */
627 
628  // compute the width
629  PFClusterWidthAlgo pfwidth(pfClusters_[iphot]);
630 
631  double correctedEnergy=pfCand[photPFCandidateIndex_[iphot]].ecalEnergy();
632  reco::SuperCluster mySuperCluster(correctedEnergy,math::XYZPoint(posX,posY,posZ));
633  // protection against empty basic cluster collection ; the value is -2 in this case
634  if(nbasics)
635  {
636 // std::cout << "SuperCluster creation; energy " << pfCand[gsfPFCandidateIndex_[iphot]].ecalEnergy();
637 // std::cout << " " << pfCand[gsfPFCandidateIndex_[iphot]].rawEcalEnergy() << std::endl;
638 // std::cout << "Seed energy from basic " << basicClusters_[iphot][0].energy() << std::endl;
639  mySuperCluster.setSeed(basicClusterPtr_[iphot][0]);
640  }
641  else
642  {
643  // std::cout << "SuperCluster creation ; seed energy " << 0 << std::endl;
644  //std::cout << "SuperCluster creation ; energy " << pfCand[photPFCandidateIndex_[iphot]].ecalEnergy();
645  //std::cout << " " << pfCand[photPFCandidateIndex_[iphot]].rawEcalEnergy() << std::endl;
646 // std::cout << " No seed found " << 0 << std::endl;
647 // std::cout << " MVA " << pfCand[gsfPFCandidateIndex_[iphot]].mva_e_pi() << std::endl;
648  mySuperCluster.setSeed(reco::CaloClusterPtr());
649  }
650  // the seed should be the first basic cluster
651 
652  for(unsigned ibc=0;ibc<nbasics;++ibc)
653  {
654  mySuperCluster.addCluster(basicClusterPtr_[iphot][ibc]);
655  // std::cout <<"Adding Ref to SC " << basicClusterPtr_[iphot][ibc].index() << std::endl;
656  const std::vector< std::pair<DetId, float> > & v1 = basicClusters_[iphot][ibc].hitsAndFractions();
657  // std::cout << " Number of cells " << v1.size() << std::endl;
658  for( std::vector< std::pair<DetId, float> >::const_iterator diIt = v1.begin();
659  diIt != v1.end();
660  ++diIt ) {
661  // std::cout << " Adding DetId " << (diIt->first).rawId() << " " << diIt->second << std::endl;
662  mySuperCluster.addHitAndFraction(diIt->first,diIt->second);
663  } // loop over rechits
664  }
665 
666  unsigned nps=preshowerClusterPtr_[iphot].size();
667  for(unsigned ips=0;ips<nps;++ips)
668  {
669  mySuperCluster.addPreshowerCluster(preshowerClusterPtr_[iphot][ips]);
670  }
671 
672 
673  // Set the preshower energy
674  mySuperCluster.setPreshowerEnergy(pfCand[photPFCandidateIndex_[iphot]].pS1Energy()+
675  pfCand[photPFCandidateIndex_[iphot]].pS2Energy());
676 
677  // Set the cluster width
678  mySuperCluster.setEtaWidth(pfwidth.pflowEtaWidth());
679  mySuperCluster.setPhiWidth(pfwidth.pflowPhiWidth());
680  // Force the computation of rawEnergy_ of the reco::SuperCluster
681  mySuperCluster.rawEnergy();
682 
683  //cout << "SC energy="<< mySuperCluster.energy()<<endl;
684 
685  superClusters.push_back(mySuperCluster);
686  //std::cout << "nb super clusters in collection : "<<superClusters.size()<<std::endl;
687  }
688 }
689 
691 {
692 
693  //std::cout << "createOneLegConversions" << std::endl;
694 
695  unsigned nphot=photPFCandidateIndex_.size();
696  for(unsigned iphot=0;iphot<nphot;++iphot)
697  {
698 
699  //if (conv1legPFCandidateIndex_[iphot]==-1) cout << "No OneLegConversions to add"<<endl;
700  //else std::cout << "Phot "<<iphot<< " nOneLegConversions to add : "<<pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]].size()<<endl;
701 
702 
703  if (conv1legPFCandidateIndex_[iphot]>-1){
704 
705  for (unsigned iConv=0; iConv<pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]].size(); iConv++){
706 
708  std::vector<reco::CaloClusterPtr>matchingBC;
710  const reco::Vertex * convVtx = new reco::Vertex(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->innerPosition(), error);
711 
712  //cout << "Vtx x="<<convVtx->x() << " y="<< convVtx->y()<<" z="<<convVtx->z()<< endl;
713  //cout << "VtxError x=" << convVtx->xError() << endl;
714 
715  std::vector<reco::TrackRef> OneLegConvVector;
716  OneLegConvVector.push_back(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]);
717 
718  std::vector<reco::TrackRef> tr=OneLegConvVector;
719  std::vector<math::XYZPointF>trackPositionAtEcalVec;
720  std::vector<math::XYZPointF>innPointVec;
721  std::vector<math::XYZVectorF>trackPinVec;
722  std::vector<math::XYZVectorF>trackPoutVec;
723  math::XYZPointF trackPositionAtEcal(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
724  outerPosition().X(),
725  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
726  outerPosition().Y(),
727  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
728  outerPosition().Z());
729  math::XYZPointF innPoint(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
730  innerPosition().X(),
731  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
732  innerPosition().Y(),
733  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
734  innerPosition().Z());
735  math::XYZVectorF trackPin(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
736  innerMomentum().X(),
737  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
738  innerMomentum().Y(),
739  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
740  innerMomentum().Z());
741  math::XYZVectorF trackPout(pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
742  outerMomentum().X(),
743  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
744  outerMomentum().Y(),
745  pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->
746  outerMomentum().Z());
747  float DCA=pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]][iConv]->d0();
748  trackPositionAtEcalVec.push_back(trackPositionAtEcal);
749  innPointVec.push_back(innPoint);
750  trackPinVec.push_back(trackPin);
751  trackPoutVec.push_back(trackPout);
752  std::vector< float > OneLegMvaVector;
753  reco::Conversion myOneLegConversion(scPtrVec,
754  OneLegConvVector,
755  trackPositionAtEcalVec,
756  *convVtx,
757  matchingBC,
758  DCA,
759  innPointVec,
760  trackPinVec,
761  trackPoutVec,
762  pfSingleLegConvMva_[conv1legPFCandidateIndex_[iphot]][iConv],
764  OneLegMvaVector.push_back(pfSingleLegConvMva_[conv1legPFCandidateIndex_[iphot]][iConv]);
765  myOneLegConversion.setOneLegMVA(OneLegMvaVector);
766  //reco::Conversion myOneLegConversion(scPtrVec,
767  //OneLegConvVector, *convVtx, reco::Conversion::pflow);
768 
769  /*
770  std::cout << "One leg conversion created" << endl;
771  std::vector<edm::RefToBase<reco::Track> > convtracks = myOneLegConversion.tracks();
772  const std::vector<float> mvalist = myOneLegConversion.oneLegMVA();
773 
774  cout << "nTracks="<< convtracks.size()<<endl;
775  for (unsigned int itk=0; itk<convtracks.size(); itk++){
776  //double convtrackpt = convtracks[itk]->pt();
777  std::cout << "Track pt="<< convtracks[itk]->pt() << std::endl;
778  std::cout << "Track mva="<< mvalist[itk] << std::endl;
779  }
780  */
781  oneLegConversions.push_back(myOneLegConversion);
782 
783  //cout << "OneLegConv added"<<endl;
784 
785  }
786  }
787  }
788 }
789 
790 
792 {
793 
794  //std::cout << "createPhotonCores" << std::endl;
795 
796  unsigned nphot=photPFCandidateIndex_.size();
797 
798  unsigned i1legtot = 0;
799 
800  for(unsigned iphot=0;iphot<nphot;++iphot)
801  {
802  //std::cout << "iphot="<<iphot<<std::endl;
803 
804  reco::PhotonCore myPhotonCore;
805 
806  reco::SuperClusterRef SCref(reco::SuperClusterRef(superClustersHandle, iphot));
807 
808  myPhotonCore.setPFlowPhoton(true);
809  myPhotonCore.setStandardPhoton(false);
810  myPhotonCore.setParentSuperCluster(SCref);
811  myPhotonCore.setSuperCluster(egSCRef_[iphot]);
812 
813  reco::ElectronSeedRefVector pixelSeeds = egPhotonRef_[iphot]->electronPixelSeeds();
814  for (unsigned iseed=0; iseed<pixelSeeds.size(); iseed++){
815  myPhotonCore.addElectronPixelSeed(pixelSeeds[iseed]);
816  }
817 
818  //cout << "PhotonCores : SC OK" << endl;
819 
820  //cout << "conv1legPFCandidateIndex_[iphot]="<<conv1legPFCandidateIndex_[iphot]<<endl;
821  //cout << "conv2legPFCandidateIndex_[iphot]="<<conv2legPFCandidateIndex_[iphot]<<endl;
822 
823  if (conv1legPFCandidateIndex_[iphot]>-1){
824  for (unsigned int iConv=0; iConv<pfSingleLegConv_[conv1legPFCandidateIndex_[iphot]].size(); iConv++){
825 
826  const reco::ConversionRef & OneLegRef(reco::ConversionRef(oneLegConversionHandle, i1legtot));
827  myPhotonCore.addOneLegConversion(OneLegRef);
828 
829  //cout << "PhotonCores : 1-leg OK" << endl;
830  /*
831  cout << "Testing 1-leg :"<<endl;
832  const reco::ConversionRefVector & conv = myPhotonCore.conversionsOneLeg();
833  for (unsigned int iconv=0; iconv<conv.size(); iconv++){
834  cout << "Testing 1-leg : iconv="<<iconv<<endl;
835  cout << "Testing 1-leg : nTracks="<<conv[iconv]->nTracks()<<endl;
836  cout << "Testing 1-leg : EoverP="<<conv[iconv]->EoverP()<<endl;
837  std::vector<edm::RefToBase<reco::Track> > convtracks = conv[iconv]->tracks();
838  for (unsigned int itk=0; itk<convtracks.size(); itk++){
839  //double convtrackpt = convtracks[itk]->pt();
840  std::cout << "Testing 1-leg : Track pt="<< convtracks[itk]->pt() << std::endl;
841  // std::cout << "Track mva="<< mvalist[itk] << std::endl;
842  }
843  }
844  */
845 
846  i1legtot++;
847  }
848  }
849 
850  if (conv2legPFCandidateIndex_[iphot]>-1){
851  for(unsigned int iConv=0; iConv<pfConv_[conv2legPFCandidateIndex_[iphot]].size(); iConv++) {
852 
853  const reco::ConversionRef & TwoLegRef(pfConv_[conv2legPFCandidateIndex_[iphot]][iConv]);
854  myPhotonCore.addConversion(TwoLegRef);
855 
856  }
857  //cout << "PhotonCores : 2-leg OK" << endl;
858 
859  /*
860  cout << "Testing 2-leg :"<<endl;
861  const reco::ConversionRefVector & conv = myPhotonCore.conversions();
862  for (unsigned int iconv=0; iconv<conv.size(); iconv++){
863  cout << "Testing 2-leg : iconv="<<iconv<<endl;
864  cout << "Testing 2-leg : nTracks="<<conv[iconv]->nTracks()<<endl;
865  cout << "Testing 2-leg : EoverP="<<conv[iconv]->EoverP()<<endl;
866  std::vector<edm::RefToBase<reco::Track> > convtracks = conv[iconv]->tracks();
867  for (unsigned int itk=0; itk<convtracks.size(); itk++){
868  //double convtrackpt = convtracks[itk]->pt();
869  std::cout << "Testing 2-leg : Track pt="<< convtracks[itk]->pt() << std::endl;
870  // std::cout << "Track mva="<< mvalist[itk] << std::endl;
871  }
872  }
873  */
874  }
875 
876  photonCores.push_back(myPhotonCore);
877 
878  }
879 
880  //std::cout << "end of createPhotonCores"<<std::endl;
881 }
882 
883 
885 {
886 
887  //cout << "createPhotons" << endl;
888 
889  unsigned nphot=photPFCandidateIndex_.size();
890 
891  for(unsigned iphot=0;iphot<nphot;++iphot)
892  {
893  //std::cout << "iphot="<<iphot<<std::endl;
894 
895  reco::PhotonCoreRef PCref(reco::PhotonCoreRef(photonCoresHandle, iphot));
896 
897  math::XYZPoint vtx(0.,0.,0.);
898  if (!vertexCollection.empty()) vtx = vertexCollection.begin()->position();
899  //std::cout << "vtx made" << std::endl;
900 
901  math::XYZVector direction = PCref->parentSuperCluster()->position() - vtx;
902 
903  //It could be that pfSC energy gives not the best resolution : use smaller agregates for some cases ?
904  math::XYZVector P3 = direction.unit() * PCref->parentSuperCluster()->energy();
905  LorentzVector P4(P3.x(), P3.y(), P3.z(), PCref->parentSuperCluster()->energy());
906 
907  reco::Photon myPhoton(P4, PCref->parentSuperCluster()->position(), PCref, vtx);
908  //cout << "photon created"<<endl;
909 
910 
911 
912  reco::Photon::ShowerShape showerShape;
913  reco::Photon::SaturationInfo saturationInfo;
914  reco::Photon::FiducialFlags fiducialFlags;
915  reco::Photon::IsolationVariables isolationVariables03;
916  reco::Photon::IsolationVariables isolationVariables04;
917 
918  showerShape.e1x5= egPhotonRef_[iphot]->e1x5();
919  showerShape.e2x5= egPhotonRef_[iphot]->e2x5();
920  showerShape.e3x3= egPhotonRef_[iphot]->e3x3();
921  showerShape.e5x5= egPhotonRef_[iphot]->e5x5();
922  showerShape.maxEnergyXtal = egPhotonRef_[iphot]->maxEnergyXtal();
923  showerShape.sigmaEtaEta = egPhotonRef_[iphot]->sigmaEtaEta();
924  showerShape.sigmaIetaIeta = egPhotonRef_[iphot]->sigmaIetaIeta();
925  showerShape.hcalDepth1OverEcal = egPhotonRef_[iphot]->hadronicDepth1OverEm();
926  showerShape.hcalDepth2OverEcal = egPhotonRef_[iphot]->hadronicDepth2OverEm();
927  myPhoton.setShowerShapeVariables ( showerShape );
928 
929 
930  saturationInfo.nSaturatedXtals = egPhotonRef_[iphot]->nSaturatedXtals();
931  saturationInfo.isSeedSaturated = egPhotonRef_[iphot]->isSeedSaturated();
932  myPhoton.setSaturationInfo(saturationInfo);
933 
934  fiducialFlags.isEB = egPhotonRef_[iphot]->isEB();
935  fiducialFlags.isEE = egPhotonRef_[iphot]->isEE();
936  fiducialFlags.isEBEtaGap = egPhotonRef_[iphot]->isEBEtaGap();
937  fiducialFlags.isEBPhiGap = egPhotonRef_[iphot]->isEBPhiGap();
938  fiducialFlags.isEERingGap = egPhotonRef_[iphot]->isEERingGap();
939  fiducialFlags.isEEDeeGap = egPhotonRef_[iphot]->isEEDeeGap();
940  fiducialFlags.isEBEEGap = egPhotonRef_[iphot]->isEBEEGap();
941  myPhoton.setFiducialVolumeFlags ( fiducialFlags );
942 
943  isolationVariables03.ecalRecHitSumEt = egPhotonRef_[iphot]->ecalRecHitSumEtConeDR03();
944  isolationVariables03.hcalTowerSumEt = egPhotonRef_[iphot]->hcalTowerSumEtConeDR03();
945  isolationVariables03.hcalDepth1TowerSumEt = egPhotonRef_[iphot]->hcalDepth1TowerSumEtConeDR03();
946  isolationVariables03.hcalDepth2TowerSumEt = egPhotonRef_[iphot]->hcalDepth2TowerSumEtConeDR03();
947  isolationVariables03.trkSumPtSolidCone = egPhotonRef_[iphot]->trkSumPtSolidConeDR03();
948  isolationVariables03.trkSumPtHollowCone = egPhotonRef_[iphot]->trkSumPtHollowConeDR03();
949  isolationVariables03.nTrkSolidCone = egPhotonRef_[iphot]->nTrkSolidConeDR03();
950  isolationVariables03.nTrkHollowCone = egPhotonRef_[iphot]->nTrkHollowConeDR03();
951  isolationVariables04.ecalRecHitSumEt = egPhotonRef_[iphot]->ecalRecHitSumEtConeDR04();
952  isolationVariables04.hcalTowerSumEt = egPhotonRef_[iphot]->hcalTowerSumEtConeDR04();
953  isolationVariables04.hcalDepth1TowerSumEt = egPhotonRef_[iphot]->hcalDepth1TowerSumEtConeDR04();
954  isolationVariables04.hcalDepth2TowerSumEt = egPhotonRef_[iphot]->hcalDepth2TowerSumEtConeDR04();
955  isolationVariables04.trkSumPtSolidCone = egPhotonRef_[iphot]->trkSumPtSolidConeDR04();
956  isolationVariables04.trkSumPtHollowCone = egPhotonRef_[iphot]->trkSumPtHollowConeDR04();
957  isolationVariables04.nTrkSolidCone = egPhotonRef_[iphot]->nTrkSolidConeDR04();
958  isolationVariables04.nTrkHollowCone = egPhotonRef_[iphot]->nTrkHollowConeDR04();
959  myPhoton.setIsolationVariables(isolationVariables04, isolationVariables03);
960 
961 
962 
963 
965  myPFIso.chargedHadronIso=(*isolationValues[0])[CandidatePtr_[iphot]];
966  myPFIso.photonIso=(*isolationValues[1])[CandidatePtr_[iphot]];
967  myPFIso.neutralHadronIso=(*isolationValues[2])[CandidatePtr_[iphot]];
968  myPhoton.setPflowIsolationVariables(myPFIso);
969 
970  reco::Photon::PflowIDVariables myPFVariables;
971 
972  reco::Mustache myMustache;
973  myMustache.MustacheID(*(myPhoton.parentSuperCluster()), myPFVariables.nClusterOutsideMustache, myPFVariables.etOutsideMustache );
974  myPFVariables.mva = pfPhotonMva_[iphot];
975  myPhoton.setPflowIDVariables(myPFVariables);
976 
977  //cout << "chargedHadronIso="<<myPhoton.chargedHadronIso()<<" photonIso="<<myPhoton.photonIso()<<" neutralHadronIso="<<myPhoton.neutralHadronIso()<<endl;
978 
979  // set PF-regression energy
980  myPhoton.setCorrectedEnergy(reco::Photon::regression2,energyRegression_[iphot],energyRegressionError_[iphot],false);
981 
982 
983  /*
984  if (basicClusters_[iphot].size()>0){
985  // Cluster shape variables
986  //Algorithms from EcalClusterTools could be adapted to PF photons ? (not based on 5x5 BC)
987  //It happens that energy computed in eg e5x5 is greater than pfSC energy (EcalClusterTools gathering energies from adjacent crystals even if not belonging to the SC)
988  const EcalRecHitCollection* hits = 0 ;
989  int subdet = PCref->parentSuperCluster()->seed()->hitsAndFractions()[0].first.subdetId();
990  if (subdet==EcalBarrel) hits = barrelRecHits;
991  else if (subdet==EcalEndcap) hits = endcapRecHits;
992  const CaloGeometry* geometry = theCaloGeom_.product();
993 
994  float maxXtal = EcalClusterTools::eMax( *(PCref->parentSuperCluster()->seed()), &(*hits) );
995  //cout << "maxXtal="<<maxXtal<<endl;
996  float e1x5 = EcalClusterTools::e1x5( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology));
997  //cout << "e1x5="<<e1x5<<endl;
998  float e2x5 = EcalClusterTools::e2x5Max( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology));
999  //cout << "e2x5="<<e2x5<<endl;
1000  float e3x3 = EcalClusterTools::e3x3( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology));
1001  //cout << "e3x3="<<e3x3<<endl;
1002  float e5x5 = EcalClusterTools::e5x5( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology));
1003  //cout << "e5x5="<<e5x5<<endl;
1004  std::vector<float> cov = EcalClusterTools::covariances( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology), geometry);
1005  float sigmaEtaEta = sqrt(cov[0]);
1006  //cout << "sigmaEtaEta="<<sigmaEtaEta<<endl;
1007  std::vector<float> locCov = EcalClusterTools::localCovariances( *(PCref->parentSuperCluster()->seed()), &(*hits), &(*topology));
1008  float sigmaIetaIeta = sqrt(locCov[0]);
1009  //cout << "sigmaIetaIeta="<<sigmaIetaIeta<<endl;
1010  //float r9 =e3x3/(PCref->parentSuperCluster()->rawEnergy());
1011 
1012 
1013  // calculate HoE
1014  const CaloTowerCollection* hcalTowersColl = hcalTowersHandle.product();
1015  EgammaTowerIsolation towerIso1(hOverEConeSize_,0.,0.,1,hcalTowersColl) ;
1016  EgammaTowerIsolation towerIso2(hOverEConeSize_,0.,0.,2,hcalTowersColl) ;
1017  double HoE1=towerIso1.getTowerESum(&(*PCref->parentSuperCluster()))/PCref->pfSuperCluster()->energy();
1018  double HoE2=towerIso2.getTowerESum(&(*PCref->parentSuperCluster()))/PCref->pfSuperCluster()->energy();
1019  //cout << "HoE1="<<HoE1<<endl;
1020  //cout << "HoE2="<<HoE2<<endl;
1021 
1022  reco::Photon::ShowerShape showerShape;
1023  showerShape.e1x5= e1x5;
1024  showerShape.e2x5= e2x5;
1025  showerShape.e3x3= e3x3;
1026  showerShape.e5x5= e5x5;
1027  showerShape.maxEnergyXtal = maxXtal;
1028  showerShape.sigmaEtaEta = sigmaEtaEta;
1029  showerShape.sigmaIetaIeta = sigmaIetaIeta;
1030  showerShape.hcalDepth1OverEcal = HoE1;
1031  showerShape.hcalDepth2OverEcal = HoE2;
1032  myPhoton.setShowerShapeVariables ( showerShape );
1033  //cout << "shower shape variables filled"<<endl;
1034  }
1035  */
1036 
1037 
1038  photons.push_back(myPhoton);
1039 
1040  }
1041 
1042  //std::cout << "end of createPhotons"<<std::endl;
1043 }
1044 
1045 
1047 {
1048  unsigned refindex=pfbe.index();
1049  // std::cout << " N daughters " << cand.numberOfDaughters() << std::endl;
1050  reco::PFCandidate::const_iterator myDaughterCandidate=cand.begin();
1052 
1053  for(;myDaughterCandidate!=itend;++myDaughterCandidate)
1054  {
1055  const reco::PFCandidate * myPFCandidate = (const reco::PFCandidate*)&*myDaughterCandidate;
1056  if(myPFCandidate->elementsInBlocks().size()!=1)
1057  {
1058  // std::cout << " Daughter with " << myPFCandidate.elementsInBlocks().size()<< " element in block " << std::endl;
1059  return cand;
1060  }
1061  if(myPFCandidate->elementsInBlocks()[0].second==refindex)
1062  {
1063  // std::cout << " Found it " << cand << std::endl;
1064  return *myPFCandidate;
1065  }
1066  }
1067  return cand;
1068 }
1069 
size
Write out results.
T getParameter(std::string const &) const
Abstract base class for a PFBlock element (track, cluster...)
const math::XYZPoint & position() const
cluster centroid position
Definition: CaloCluster.h:129
void setSuperCluster(const reco::SuperClusterRef &r)
set reference to SuperCluster
Definition: PhotonCore.h:50
void addOneLegConversion(const reco::ConversionRef &r)
add single ConversionRef to the vector of Refs
Definition: PhotonCore.h:56
void MustacheID(const CaloClusterPtrVector &clusters, int &nclusters, float &EoutsideMustache)
Definition: Mustache.cc:162
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:127
void addConversion(const reco::ConversionRef &r)
add single ConversionRef to the vector of Refs
Definition: PhotonCore.h:54
void createBasicClusterPtrs(const edm::OrphanHandle< reco::BasicClusterCollection > &basicClustersHandle)
void addHitAndFraction(DetId id, float fraction)
Definition: CaloCluster.h:188
Particle flow cluster, see clustering algorithm in PFClusterAlgo.
Definition: PFCluster.h:47
float mva_nothing_gamma() const
mva for gamma detection
Definition: PFCandidate.h:332
virtual const PFClusterRef & clusterRef() const
Type type() const
double pflowPhiWidth() const
bool exists(std::string const &parameterName) const
checks if a parameter exists
#define X(str)
Definition: MuonsGrabber.cc:48
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
const std::vector< std::pair< DetId, float > > & hitsAndFractions() const
Definition: CaloCluster.h:195
ErrorD< N >::type type
Definition: Error.h:33
void setSeed(const CaloClusterPtr &r)
list of used xtals by DetId // now inherited by CaloCluster
Definition: SuperCluster.h:96
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< float > > XYZPointF
point in space with cartesian internal representation
Definition: Point3D.h:10
void createSuperClusters(const reco::PFCandidateCollection &, reco::SuperClusterCollection &superClusters) const
reco::PFCandidatePhotonExtraRef photonExtraRef() const
return a reference to the photon extra
Definition: PFCandidate.cc:597
void setPhiWidth(double pw)
Definition: SuperCluster.h:62
double pflowEtaWidth() const
void createBasicCluster(const reco::PFBlockElement &, reco::BasicClusterCollection &basicClusters, std::vector< const reco::PFCluster * > &, const reco::PFCandidate &coCandidate) const
std::vector< Conversion > ConversionCollection
collectin of Conversion objects
Definition: ConversionFwd.h:9
AlgoId algo() const
algorithm identifier
Definition: CaloCluster.h:175
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
void setEtaWidth(double ew)
Definition: SuperCluster.h:63
int iEvent
Definition: GenABIO.cc:230
math::XYZTLorentzVector LorentzVector
void setParentSuperCluster(const reco::SuperClusterRef &r)
set reference to PFlow SuperCluster
Definition: PhotonCore.h:52
const CaloID & caloID() const
Definition: CaloCluster.h:186
unsigned index() const
std::vector< SuperCluster > SuperClusterCollection
collection of SuperCluser objectr
void setStandardPhoton(const bool prov)
Definition: PhotonCore.h:61
PFPhotonTranslator(const edm::ParameterSet &)
const_iterator end() const
last daughter const_iterator
Definition: Candidate.h:146
std::vector< PreshowerCluster > PreshowerClusterCollection
collection of PreshowerCluster objects
math::XYZTLorentzVector LorentzVector
math::XYZPoint Point
void createPhotonCores(const edm::OrphanHandle< reco::SuperClusterCollection > &superClustersHandle, const edm::OrphanHandle< reco::ConversionCollection > &oneLegConversionHandle, reco::PhotonCoreCollection &photonCores)
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< float > > XYZVectorF
spatial vector with cartesian internal representation
Definition: Vector3D.h:17
void createPreshowerCluster(const reco::PFBlockElement &PFBE, reco::PreshowerClusterCollection &preshowerClusters, unsigned plane) const
bool isValid() const
Definition: HandleBase.h:74
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:464
bool isNull() const
Checks for null.
Definition: Ref.h:250
double energy() const
cluster energy
Definition: PFCluster.h:82
void produce(edm::Event &, const edm::EventSetup &) override
void createPreshowerClusterPtrs(const edm::OrphanHandle< reco::PreshowerClusterCollection > &preshowerClustersHandle)
bool fetchCandidateCollection(edm::Handle< reco::PFCandidateCollection > &c, const edm::InputTag &tag, const edm::Event &iEvent) const
void setPFlowPhoton(const bool prov)
set the provenance
Definition: PhotonCore.h:60
const reco::PFCandidate & correspondingDaughterCandidate(const reco::PFCandidate &cand, const reco::PFBlockElement &pfbe) const
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
double rawEnergy() const
raw uncorrected energy (sum of energies of component BasicClusters)
Definition: SuperCluster.h:47
DetId seed() const
return DetId of seed
Definition: CaloCluster.h:205
int iseed
Definition: AMPTWrapper.h:124
T const * product() const
Definition: Handle.h:81
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
std::vector< Photon > PhotonCollection
collectin of Photon objects
Definition: PhotonFwd.h:9
void addPreshowerCluster(const CaloClusterPtr &r)
add reference to constituent BasicCluster
Definition: SuperCluster.h:117
math::XYZVector Vector
void createPhotons(reco::VertexCollection &vertexCollection, edm::Handle< reco::PhotonCollection > &egPhotons, const edm::OrphanHandle< reco::PhotonCoreCollection > &photonCoresHandle, const IsolationValueMaps &isolationValues, reco::PhotonCollection &photons)
std::vector< edm::Handle< edm::ValueMap< double > > > IsolationValueMaps
std::vector< BasicCluster > BasicClusterCollection
collection of BasicCluster objects
void createOneLegConversions(const edm::OrphanHandle< reco::SuperClusterCollection > &superClustersHandle, reco::ConversionCollection &oneLegConversions)
std::vector< PhotonCore > PhotonCoreCollection
collectin of PhotonCore objects
Definition: PhotonCoreFwd.h:9
void addCluster(const CaloClusterPtr &r)
add reference to constituent BasicCluster
Definition: SuperCluster.h:111
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:39
fixed size matrix
HLT enums.
void setOneLegMVA(const std::vector< float > &mva)
set the MVS output from PF for one leg conversions
Definition: Conversion.h:171
const_iterator begin() const
first daughter const_iterator
Definition: Candidate.h:144
size_type size() const
Size of the RefVector.
Definition: RefVector.h:107
void addElectronPixelSeed(const reco::ElectronSeedRef &r)
set electron pixel seed ref
Definition: PhotonCore.h:58
virtual ParticleType particleId() const
Definition: PFCandidate.h:373
const ElementsInBlocks & elementsInBlocks() const
Definition: PFCandidate.cc:687
def move(src, dest)
Definition: eostools.py:510
void setPreshowerEnergy(double preshowerEnergy)
Definition: SuperCluster.h:59
reco::SuperClusterRef superClusterRef() const
return a reference to the corresponding SuperCluster if any
Definition: PFCandidate.cc:605