Skip to content

Commit 8062bf4

Browse files
davidlange6davidlange6
committed
Merge pull request cms-sw#2240 from ktf/rebase-1404
Calo towers speed-up
2 parents 469a069 + e76ebb0 commit 8062bf4

File tree

14 files changed

+487
-187
lines changed

14 files changed

+487
-187
lines changed

‎DataFormats/CaloTowers/interface/CaloTower.h

+17
Original file line numberDiff line numberDiff line change
@@ -44,10 +44,27 @@ class CaloTower : public reco::LeafCandidate {
4444
const LorentzVector& p4,
4545
const GlobalPoint& emPosition, const GlobalPoint& hadPosition);
4646

47+
CaloTower(CaloTowerDetId id,
48+
float emE, float hadE, float outerE,
49+
int ecal_tp, int hcal_tp,
50+
GlobalVector p3, float iEnergy, bool massless,
51+
GlobalPoint emPosition, GlobalPoint hadPosition);
52+
53+
CaloTower(CaloTowerDetId id,
54+
float emE, float hadE, float outerE,
55+
int ecal_tp, int hcal_tp,
56+
GlobalVector p3, float iEnergy, float imass,
57+
GlobalPoint emPosition, GlobalPoint hadPosition);
58+
59+
4760

4861
// setters
4962
void addConstituent( DetId id ) { constituents_.push_back( id ); }
5063
void addConstituents( const std::vector<DetId>& ids );
64+
#if defined(__GXX_EXPERIMENTAL_CXX0X__)
65+
void setConstituents( std::vector<DetId>&& ids ) { constituents_=std::move(ids);}
66+
#endif
67+
5168
void setEcalTime(int t) { ecalTime_ = t; };
5269
void setHcalTime(int t) { hcalTime_ = t; };
5370

‎DataFormats/CaloTowers/src/CaloTower.cc

+24
Original file line numberDiff line numberDiff line change
@@ -21,6 +21,7 @@ CaloTower::CaloTower(const CaloTowerDetId& id,
2121
emLvl1_(ecal_tp), hadLvl1_(hcal_tp) {}
2222

2323

24+
2425
CaloTower::CaloTower(const CaloTowerDetId& id,
2526
double emE, double hadE, double outerE,
2627
int ecal_tp, int hcal_tp,
@@ -33,6 +34,29 @@ CaloTower::CaloTower(const CaloTowerDetId& id,
3334
emLvl1_(ecal_tp), hadLvl1_(hcal_tp) {}
3435

3536

37+
CaloTower::CaloTower(CaloTowerDetId id,
38+
float emE, float hadE, float outerE,
39+
int ecal_tp, int hcal_tp,
40+
GlobalVector p3, float iEnergy, bool massless,
41+
GlobalPoint emPos, GlobalPoint hadPos) :
42+
LeafCandidate(0, p3, iEnergy, massless, Point(0,0,0)),
43+
id_(id),
44+
emPosition_(emPos), hadPosition_(hadPos),
45+
emE_(emE), hadE_(hadE), outerE_(outerE),
46+
emLvl1_(ecal_tp), hadLvl1_(hcal_tp) {}
47+
48+
CaloTower::CaloTower(CaloTowerDetId id,
49+
float emE, float hadE, float outerE,
50+
int ecal_tp, int hcal_tp,
51+
GlobalVector p3, float iEnergy, float imass,
52+
GlobalPoint emPos, GlobalPoint hadPos) :
53+
LeafCandidate(0, p3, iEnergy, imass, Point(0,0,0)),
54+
id_(id),
55+
emPosition_(emPos), hadPosition_(hadPos),
56+
emE_(emE), hadE_(hadE), outerE_(outerE),
57+
emLvl1_(ecal_tp), hadLvl1_(hcal_tp) {}
58+
59+
3660
// recalculated momentum-related quantities wrt user provided vertex Z position
3761

3862

‎DataFormats/Candidate/interface/LeafCandidate.h

+30-5
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,8 @@
1313
#include "DataFormats/Candidate/interface/iterator_imp_specific.h"
1414

1515
#include "DataFormats/Math/interface/PtEtaPhiMass.h"
16+
#include "DataFormats/GeometryVector/interface/GlobalVector.h"
17+
1618

1719
namespace reco {
1820

@@ -33,6 +35,9 @@ namespace reco {
3335

3436
typedef unsigned int index;
3537

38+
static double magd(GlobalVector v) { return std::sqrt(double(v.x())*double(v.x()) + double(v.y())*double(v.y()) + double(v.z())*double(v.z()) );}
39+
static double dmass(GlobalVector v, double e) { double m2 = e*e-magd(v); return m2>0 ? std::sqrt(m2) : 0;}
40+
3641
/// default constructor
3742
LeafCandidate() :
3843
qx3_(0), pt_(0), eta_(0), phi_(0), mass_(0),
@@ -65,19 +70,39 @@ namespace reco {
6570
LeafCandidate( Charge q, const LorentzVector & p4, const Point & vtx = Point( 0, 0, 0 ),
6671
int pdgId = 0, int status = 0, bool integerCharge = true ) :
6772
qx3_( q ), pt_( p4.pt() ), eta_( p4.eta() ), phi_( p4.phi() ), mass_( p4.mass() ),
68-
vertex_( vtx ), pdgId_( pdgId ), status_( status ),
69-
cachePolarFixed_( false ), cacheCartesianFixed_( false ) {
73+
vertex_( vtx ), pdgId_( pdgId ), status_( status ), p4Cartesian_(p4),
74+
cachePolarFixed_( false ), cacheCartesianFixed_( true ) {
7075
if ( integerCharge ) qx3_ *= 3;
7176
}
7277
/// constructor from values
7378
LeafCandidate( Charge q, const PolarLorentzVector & p4, const Point & vtx = Point( 0, 0, 0 ),
7479
int pdgId = 0, int status = 0, bool integerCharge = true ) :
7580
qx3_( q ), pt_( p4.pt() ), eta_( p4.eta() ), phi_( p4.phi() ), mass_( p4.mass() ),
76-
vertex_( vtx ), pdgId_( pdgId ), status_( status ),
77-
cachePolarFixed_( false ), cacheCartesianFixed_( false ){
81+
vertex_( vtx ), pdgId_( pdgId ), status_( status ), p4Polar_(p4),
82+
cachePolarFixed_( true ), cacheCartesianFixed_( false ){
7883
if ( integerCharge ) qx3_ *= 3;
7984
}
8085

86+
/// constructor from values
87+
LeafCandidate( Charge q, const GlobalVector & p3, float iEnergy, bool massless, const Point & vtx = Point( 0, 0, 0 ),
88+
int pdgId = 0, int status = 0, bool integerCharge = true ) :
89+
qx3_( q ), pt_( p3.perp() ), eta_( p3.eta() ), phi_( p3.phi() ), mass_(massless ? 0. : dmass(p3,iEnergy) ),
90+
vertex_( vtx ), pdgId_( pdgId ), status_( status ), p4Polar_(pt_,eta_,phi_,mass_), p4Cartesian_(p3.x(),p3.y(),p3.z(), massless ? magd(p3) : iEnergy),
91+
cachePolarFixed_( true ), cacheCartesianFixed_( true ) {
92+
if ( integerCharge ) qx3_ *= 3;
93+
}
94+
95+
96+
/// constructor from values
97+
LeafCandidate( Charge q, const GlobalVector & p3, float iEnergy, float imass, const Point & vtx = Point( 0, 0, 0 ),
98+
int pdgId = 0, int status = 0, bool integerCharge = true ) :
99+
qx3_( q ), pt_( p3.perp() ), eta_( p3.eta() ), phi_( p3.phi() ), mass_(imass),
100+
vertex_( vtx ), pdgId_( pdgId ), status_( status ), p4Polar_(pt_,eta_,phi_,mass_), p4Cartesian_(p3.x(),p3.y(),p3.z(), iEnergy),
101+
cachePolarFixed_( true ), cacheCartesianFixed_( true ) {
102+
if ( integerCharge ) qx3_ *= 3;
103+
}
104+
105+
81106
/// destructor
82107
virtual ~LeafCandidate();
83108
/// first daughter const_iterator
@@ -133,7 +158,7 @@ namespace reco {
133158
/// energy
134159
virtual double energy() const GCC11_FINAL { cacheCartesian(); return p4Cartesian_.E(); }
135160
/// transverse energy
136-
virtual double et() const GCC11_FINAL { cachePolar(); return p4Polar_.Et(); }
161+
virtual double et() const GCC11_FINAL { cachePolar(); return (pt_<=0) ? 0 : p4Polar_.Et(); }
137162
/// mass
138163
virtual float mass() const GCC11_FINAL { return mass_; }
139164
/// mass squared

‎DataFormats/Candidate/test/testCandidate.cc

+42-1
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@
33
#include "DataFormats/Candidate/interface/Candidate.h"
44
#include "DataFormats/Candidate/interface/CandidateWithRef.h"
55
#include <memory>
6+
#include <iostream>
67

78
class testCandidate : public CppUnit::TestFixture {
89
CPPUNIT_TEST_SUITE(testCandidate);
@@ -59,7 +60,10 @@ namespace reco {
5960
}
6061

6162
void testCandidate::checkAll() {
62-
reco::Particle::LorentzVector p( 1.0, 2.0, 3.0, 4.0 );
63+
reco::Particle::LorentzVector p( 1.0, 2.0, 3.0, 5.0 );
64+
GlobalVector v(1.0, 2.0, 3.0);
65+
reco::LeafCandidate::PolarLorentzVector pl(p);
66+
6367
reco::Particle::Charge q( 1 );
6468
int x = 123, y0 = 111, y1 = 222;
6569
std::auto_ptr<reco::Candidate> c( new test::DummyCandidate1( p, q, x, y0, y1 ) );
@@ -78,4 +82,41 @@ void testCandidate::checkAll() {
7882
reco::Candidate::const_iterator b = c->begin(), e = c->end();
7983
CPPUNIT_ASSERT( b == e );
8084
CPPUNIT_ASSERT( e - b == 0 );
85+
86+
// test constructors
87+
88+
reco::LeafCandidate c1(q,p);
89+
reco::LeafCandidate c2(q,pl);
90+
reco::LeafCandidate c3(q,v,5.f,c1.mass());
91+
92+
auto ftoi = [](float x)->int { int i; memcpy(&i,&x,4); return i;};
93+
auto print = [](float a, float b)->bool { std::cout << "\nwhat? " << a <<' ' << b << std::endl; return false;};
94+
auto ok = [&](float a, float b)->bool { return std::abs(ftoi(a)-ftoi(b))<10 ? true : print(a,b); };
95+
96+
CPPUNIT_ASSERT(ok(c1.pt(),c2.pt()));
97+
CPPUNIT_ASSERT(ok(c1.eta(),c2.eta()));
98+
CPPUNIT_ASSERT(ok(c1.phi(),c2.phi()));
99+
CPPUNIT_ASSERT(ok(c1.mass(),c2.mass()));
100+
101+
CPPUNIT_ASSERT(ok(c1.y(),c2.y()));
102+
103+
CPPUNIT_ASSERT(ok(c1.pt(),c3.pt()));
104+
CPPUNIT_ASSERT(ok(c1.eta(),c3.eta()));
105+
CPPUNIT_ASSERT(ok(c1.phi(),c3.phi()));
106+
CPPUNIT_ASSERT(ok(c1.mass(),c3.mass()));
107+
108+
CPPUNIT_ASSERT(ok(c1.y(),c3.y()));
109+
110+
CPPUNIT_ASSERT(ok(c1.px(),c2.px()));
111+
CPPUNIT_ASSERT(ok(c1.py(),c2.py()));
112+
CPPUNIT_ASSERT(ok(c1.pz(),c2.pz()));
113+
CPPUNIT_ASSERT(ok(c1.energy(),c2.energy()));
114+
115+
CPPUNIT_ASSERT(ok(c1.px(),c3.px()));
116+
CPPUNIT_ASSERT(ok(c1.py(),c3.py()));
117+
CPPUNIT_ASSERT(ok(c1.pz(),c3.pz()));
118+
CPPUNIT_ASSERT(ok(c1.energy(),c3.energy()));
119+
120+
121+
81122
}

‎DataFormats/Common/interface/SortedCollection.h

+9
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,8 @@ unreliable if such duplicate entries are made.
3333
#include "DataFormats/Provenance/interface/ProductID.h"
3434
#include "FWCore/Utilities/interface/EDMException.h"
3535

36+
#include "FWCore/Utilities/interface/GCC11Compatibility.h"
37+
3638
#include <algorithm>
3739
#include <typeinfo>
3840
#include <vector>
@@ -99,6 +101,13 @@ namespace edm {
99101
//SortedCollection(InputIterator b, InputIterator e);
100102

101103
void push_back(T const& t);
104+
#if defined(__GXX_EXPERIMENTAL_CXX0X__)
105+
void push_back(T && t) { obj.push_back(t);}
106+
107+
template<typename... Args >
108+
void emplace_back( Args&&... args ) { obj.emplace_back(args...);}
109+
#endif
110+
void pop_back() { obj.pop_back(); }
102111

103112
void swap(SortedCollection& other);
104113

‎Geometry/CaloGeometry/interface/CaloCellGeometry.h

+14-1
Original file line numberDiff line numberDiff line change
@@ -73,7 +73,9 @@ class CaloCellGeometry
7373
virtual const CornersVec& getCorners() const = 0 ;
7474

7575
/// Returns the position of reference for this cell
76-
const GlobalPoint& getPosition() const {return m_refPoint ; }
76+
const GlobalPoint& getPosition() const {return m_refPoint;}
77+
const GlobalPoint& getBackPoint() const {return m_backPoint;}
78+
7779
float etaPos() const { return m_eta;}
7880
float phiPos() const { return m_phi;}
7981

@@ -123,10 +125,21 @@ class CaloCellGeometry
123125
getCorners()[2].eta());
124126
m_dPhi = std::abs(getCorners()[0].phi() -
125127
getCorners()[2].phi());
128+
initBack();
129+
130+
}
131+
132+
void initBack() const {
133+
// from CaloTower code
134+
CornersVec const & cv = getCorners();
135+
m_backPoint = GlobalPoint(0.25 * (cv[4].x() + cv[5].x() + cv[6].x() + cv[7].x()),
136+
0.25 * (cv[4].y() + cv[5].y() + cv[6].y() + cv[7].y()),
137+
0.25 * (cv[4].z() + cv[5].z() + cv[6].z() + cv[7].z()));
126138
}
127139

128140
private:
129141
GlobalPoint m_refPoint ;
142+
mutable GlobalPoint m_backPoint ;
130143
mutable CornersVec m_corners ;
131144
const CCGFloat* m_parms ;
132145
float m_eta, m_phi;

‎Geometry/CaloTopology/src/CaloTowerConstituentsMap.cc

+4
Original file line numberDiff line numberDiff line change
@@ -60,6 +60,10 @@ void CaloTowerConstituentsMap::assign(const DetId& cell, const CaloTowerDetId& t
6060

6161
void CaloTowerConstituentsMap::sort() {
6262
m_items.sort();
63+
64+
// for (auto const & it : m_items)
65+
// std::cout << std::hex << it.cell.rawId() << " " << it.tower.rawId() << std::dec << std::endl;
66+
6367
}
6468

6569
std::vector<DetId> CaloTowerConstituentsMap::constituentsOf(const CaloTowerDetId& id) const {

‎RecoJets/JetProducers/plugins/VirtualJetProducer.cc

+2
Original file line numberDiff line numberDiff line change
@@ -415,6 +415,7 @@ void VirtualJetProducer::inputTowers( )
415415
inEnd = inputs_.end(), i = inBegin;
416416
for (; i != inEnd; ++i ) {
417417
reco::CandidatePtr input = *i;
418+
// std::cout << "CaloTowerVI jets " << input->pt() << " " << input->et() << ' '<< input->energy() << ' ' << (isAnomalousTower(input) ? " bad" : " ok") << std::endl;
418419
if (edm::isNotFinite(input->pt())) continue;
419420
if (input->et() <inputEtMin_) continue;
420421
if (input->energy()<inputEMin_) continue;
@@ -697,6 +698,7 @@ void VirtualJetProducer::writeJets( edm::Event & iEvent, edm::EventSetup const&
697698

698699

699700
// std::cout << "area " << ijet << " " << jetArea << " " << Area<T>::get(jet) << std::endl;
701+
// std::cout << "JetVI " << ijet << jet.pt() << " " << jet.et() << ' '<< jet.energy() << ' '<< jet.mass() << std::endl;
700702

701703
// add to the list
702704
jets->push_back(jet);

‎RecoLocalCalo/CaloTowersCreator/interface/CaloTowersCreationAlgo.h

+23-10
Original file line numberDiff line numberDiff line change
@@ -23,6 +23,7 @@
2323

2424
// need if we want to store the handles
2525
#include "FWCore/Framework/interface/ESHandle.h"
26+
#include <tuple>
2627

2728

2829
#include <map>
@@ -46,6 +47,9 @@ class DetId;
4647

4748
class CaloTowersCreationAlgo {
4849
public:
50+
51+
int nalgo=-1;
52+
4953
CaloTowersCreationAlgo();
5054

5155
CaloTowersCreationAlgo(double EBthreshold, double EEthreshold,
@@ -111,6 +115,8 @@ class CaloTowersCreationAlgo {
111115
// The key is the calotower id.
112116
void makeHcalDropChMap();
113117

118+
void makeEcalBadChs();
119+
114120
void begin();
115121
void process(const HBHERecHitCollection& hbhe);
116122
void process(const HORecHitCollection& ho);
@@ -139,7 +145,8 @@ class CaloTowersCreationAlgo {
139145
// Called in assignHit to check if the energy should be added to
140146
// calotower, and how to flag the channel
141147
unsigned int hcalChanStatusForCaloTower(const CaloRecHit* hit);
142-
unsigned int ecalChanStatusForCaloTower(const CaloRecHit* hit);
148+
149+
std::tuple<unsigned int,bool> ecalChanStatusForCaloTower(const CaloRecHit* hit);
143150

144151
// Channel flagging is based on acceptable severity levels specified in the
145152
// configuration file. These methods are used to pass the values read in
@@ -192,31 +199,33 @@ class CaloTowersCreationAlgo {
192199
GlobalPoint hadSegmentShwrPos(DetId detId, float fracDepth);
193200
// "effective" point for the EM/HAD shower in CaloTower
194201
// position based on non-zero energy cells
195-
GlobalPoint hadShwrPos(const std::vector<std::pair<DetId,double> >& metaContains,
202+
GlobalPoint hadShwrPos(const std::vector<std::pair<DetId,float> >& metaContains,
196203
float fracDepth, double hadE);
197-
GlobalPoint emShwrPos(const std::vector<std::pair<DetId,double> >& metaContains,
204+
GlobalPoint emShwrPos(const std::vector<std::pair<DetId,float> >& metaContains,
198205
float fracDepth, double totEmE);
199206

200207
// overloaded function to get had position based on all had cells in the tower
201208
GlobalPoint hadShwrPos(CaloTowerDetId id, float fracDepth);
202209
GlobalPoint hadShwPosFromCells(DetId frontCell, DetId backCell, float fracDepth);
203210

204211
// for Chris
205-
GlobalPoint emShwrLogWeightPos(const std::vector<std::pair<DetId,double> >& metaContains,
212+
GlobalPoint emShwrLogWeightPos(const std::vector<std::pair<DetId,float> >& metaContains,
206213
float fracDepth, double totEmE);
207214

208215

209216
private:
210217

211218
struct MetaTower {
212-
MetaTower();
213-
double E, E_em, E_had, E_outer;
219+
MetaTower(){}
220+
bool empty() const { return metaConstituents.empty();}
214221
// contains also energy of RecHit
215-
std::vector< std::pair<DetId, double> > metaConstituents;
216-
double emSumTimeTimesE, hadSumTimeTimesE, emSumEForTime, hadSumEForTime; // Sum(Energy x Timing) : intermediate container
222+
std::vector< std::pair<DetId, float> > metaConstituents;
223+
CaloTowerDetId id;
224+
float E=0, E_em=0, E_had=0, E_outer=0;
225+
float emSumTimeTimesE=0, hadSumTimeTimesE=0, emSumEForTime=0, hadSumEForTime=0; // Sum(Energy x Timing) : intermediate container
217226

218227
// needed to set CaloTower status word
219-
int numBadEcalCells, numRecEcalCells, numProbEcalCells, numBadHcalCells, numRecHcalCells, numProbHcalCells;
228+
int numBadEcalCells=0, numRecEcalCells=0, numProbEcalCells=0, numBadHcalCells=0, numRecHcalCells=0, numProbHcalCells=0;
220229

221230
};
222231

@@ -315,14 +324,18 @@ class CaloTowersCreationAlgo {
315324

316325

317326
// internal map
318-
typedef std::map<CaloTowerDetId, MetaTower> MetaTowerMap;
327+
typedef std::vector<MetaTower> MetaTowerMap;
319328
MetaTowerMap theTowerMap;
329+
unsigned int theTowerMapSize=0;
320330

321331
// Number of channels in the tower that were not used in RecHit production (dead/off,...).
322332
// These channels are added to the other "bad" channels found in the recHit collection.
323333
typedef std::map<CaloTowerDetId, int> HcalDropChMap;
324334
HcalDropChMap hcalDropChMap;
325335

336+
// Number of bad Ecal channel in each tower
337+
unsigned short ecalBadChs[CaloTowerDetId::kSizeForDenseIndexing];
338+
326339
// clasification of channels in tower construction: the category definition is
327340
// affected by the setting in the configuration file
328341
//

0 commit comments

Comments
 (0)