00001
00002
00003
00004 #include <sstream>
00005
00006 #include "FWCore/Utilities/interface/Exception.h"
00007 #include "DataFormats/RecoCandidate/interface/RecoCaloTowerCandidate.h"
00008
00009
00010 #include "DataFormats/JetReco/interface/CaloJet.h"
00011
00012 using namespace reco;
00013
00014 CaloJet::CaloJet (const LorentzVector& fP4, const Point& fVertex,
00015 const Specific& fSpecific)
00016 : Jet (fP4, fVertex),
00017 m_specific (fSpecific)
00018 {}
00019
00020 CaloJet::CaloJet (const LorentzVector& fP4, const Point& fVertex,
00021 const Specific& fSpecific,
00022 const Jet::Constituents& fConstituents)
00023 : Jet (fP4, fVertex, fConstituents),
00024 m_specific (fSpecific)
00025 {}
00026
00027 CaloJet::CaloJet (const LorentzVector& fP4,
00028 const Specific& fSpecific,
00029 const Jet::Constituents& fConstituents)
00030 : Jet (fP4, Point(0,0,0), fConstituents),
00031 m_specific (fSpecific)
00032 {}
00033
00034
00036 float CaloJet::physicsEtaQuick (float fZVertex) const {
00037 return Jet::physicsEta (fZVertex, eta());
00038 }
00039
00041 float CaloJet::physicsEtaDetailed (float fZVertex) const {
00042 Jet::LorentzVector correctedMomentum;
00043 std::vector<const Candidate*> towers = getJetConstituentsQuick ();
00044 for (unsigned i = 0; i < towers.size(); ++i) {
00045 const Candidate* c = towers[i];
00046 double etaRef = Jet::physicsEta (fZVertex, c->eta());
00047 math::PtEtaPhiMLorentzVectorD p4 (c->p()/cosh(etaRef), etaRef, c->phi(), c->mass());
00048 correctedMomentum += p4;
00049 }
00050 return correctedMomentum.eta();
00051 }
00052
00054 CaloJet::LorentzVector CaloJet::physicsP4 (float fZVertex) const {
00055 double physicsEta = Jet::physicsEta (fZVertex, eta());
00056 math::PtEtaPhiMLorentzVectorD p4 (p()/cosh(physicsEta), physicsEta, phi(), mass());
00057 return CaloJet::LorentzVector (p4);
00058 }
00059
00060
00061 CaloTowerPtr CaloJet::getCaloConstituent (unsigned fIndex) const {
00062 Constituent dau = daughterPtr (fIndex);
00063
00064 if ( dau.isNonnull() && dau.isAvailable() ) {
00065
00066 const CaloTower* towerCandidate = dynamic_cast <const CaloTower*> (dau.get());
00067
00068 if (towerCandidate) {
00069
00070
00071 return edm::Ptr<CaloTower> (dau.id(), towerCandidate, dau.key() );
00072 }
00073 else {
00074 throw cms::Exception("Invalid Constituent") << "CaloJet constituent is not of CaloTowere type";
00075 }
00076
00077 }
00078
00079 else {
00080 return CaloTowerPtr();
00081 }
00082 }
00083
00084
00085 std::vector <CaloTowerPtr > CaloJet::getCaloConstituents () const {
00086 std::vector <CaloTowerPtr> result;
00087 for (unsigned i = 0; i < numberOfDaughters (); i++) result.push_back (getCaloConstituent (i));
00088 return result;
00089 }
00090
00091
00092 CaloJet* CaloJet::clone () const {
00093 return new CaloJet (*this);
00094 }
00095
00096 bool CaloJet::overlap( const Candidate & ) const {
00097 return false;
00098 }
00099
00100 std::string CaloJet::print () const {
00101 std::ostringstream out;
00102 out << Jet::print ()
00103 << " CaloJet specific:" << std::endl
00104 << " energy fractions em/had: " << emEnergyFraction () << '/' << energyFractionHadronic () << std::endl
00105 << " em energy in EB/EE/HF: " << emEnergyInEB() << '/' << emEnergyInEE() << '/' << emEnergyInHF() << std::endl
00106 << " had energy in HB/HO/HE/HF: " << hadEnergyInHB() << '/' << hadEnergyInHO() << '/' << hadEnergyInHE() << '/' << hadEnergyInHF() << std::endl
00107 << " constituent towers area: " << towersArea() << std::endl;
00108 out << " Towers:" << std::endl;
00109 std::vector <CaloTowerPtr > towers = getCaloConstituents ();
00110 for (unsigned i = 0; i < towers.size (); i++) {
00111 if (towers[i].get ()) {
00112 out << " #" << i << " " << *(towers[i]) << std::endl;
00113 }
00114 else {
00115 out << " #" << i << " tower is not available in the event" << std::endl;
00116 }
00117 }
00118 return out.str ();
00119 }
00120
00121
00122
00123
00124
00125 std::vector<CaloTowerDetId> CaloJet::getTowerIndices() const {
00126 std::vector<CaloTowerDetId> result;
00127 std::vector <CaloTowerPtr> towers = getCaloConstituents ();
00128 for (unsigned i = 0; i < towers.size(); ++i) {
00129 result.push_back (towers[i]->id());
00130 }
00131 return result;
00132 }