CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFRecHitProducerPS.cc
Go to the documentation of this file.
2 
3 #include <memory>
4 
7 
11 
13 
17 
23 
26 
27 
28 
29 using namespace std;
30 using namespace edm;
31 
33  : PFRecHitProducer(iConfig) {
34 
35 
36 
37  // access to the collections of rechits
38 
40  iConfig.getParameter<InputTag>("ecalRecHitsES");
41 }
42 
43 
44 
46 
47 
48 
49 void PFRecHitProducerPS::createRecHits(vector<reco::PFRecHit>& rechits,
50  vector<reco::PFRecHit>& rechitsCleaned,
52  const edm::EventSetup& iSetup) {
53 
54  map<unsigned, unsigned > idSortedRecHits;
55  typedef map<unsigned, unsigned >::iterator IDH;
56 
57 
58  // get the ps geometry
60  iSetup.get<CaloGeometryRecord>().get(geoHandle);
61 
62  const CaloSubdetectorGeometry *psGeometry =
63  geoHandle->getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
64 
65 
66  // ShR 28 Jul 2008: check if geometry is NULL. If so we are using
67  // partial CMS gemoetry in Pilot1/2 scenarios which do not include the preshower
68  if(!psGeometry) {
69  LogDebug("PFRecHitProducerPS") << "No EcalPreshower geometry available. putting empty PS rechits collection in event";
70  return;
71  }
72 
73 
74  // get the ps topology
75  EcalPreshowerTopology psTopology(geoHandle);
76 
77  // process rechits
79 
80 
81 
83  pRecHits);
84 
85  if(!found) {
86  ostringstream err;
87  err<<"could not find rechits "<<inputTagEcalRecHitsES_;
88  LogError("PFRecHitProducerPS")<<err.str()<<endl;
89 
90  throw cms::Exception( "MissingProduct", err.str());
91  }
92  else {
93  assert( pRecHits.isValid() );
94 
95  const EcalRecHitCollection& psrechits = *( pRecHits.product() );
97 
98  for(IT i=psrechits.begin(); i!=psrechits.end(); i++) {
99  const EcalRecHit& hit = *i;
100 
101  double energy = hit.energy();
102  if( energy < thresh_Endcap_ ) continue;
103 
104  const ESDetId& detid = hit.detid();
105  const CaloCellGeometry *thisCell = psGeometry->getGeometry(detid);
106 
107  if(!thisCell) {
108  LogError("PFRecHitProducerPS")<<"warning detid "<<detid.rawId()
109  <<" not found in preshower geometry"
110  <<endl;
111  return;
112  }
113 
114  const GlobalPoint& position = thisCell->getPosition();
115 
117 
118  switch( detid.plane() ) {
119  case 1:
120  layer = PFLayer::PS1;
121  break;
122  case 2:
123  layer = PFLayer::PS2;
124  break;
125  default:
126  LogError("PFRecHitProducerPS")
127  <<"incorrect preshower plane !! plane number "
128  <<detid.plane()<<endl;
129  assert(0);
130  }
131 
132  reco::PFRecHit *pfrh
133  = new reco::PFRecHit( detid.rawId(), layer, energy,
134  position.x(), position.y(), position.z(),
135  0,0,0 );
136 
137  const CaloCellGeometry::CornersVec& corners = thisCell->getCorners();
138  assert( corners.size() == 8 );
139 
140  pfrh->setNECorner( corners[0].x(), corners[0].y(), corners[0].z() );
141  pfrh->setSECorner( corners[1].x(), corners[1].y(), corners[1].z() );
142  pfrh->setSWCorner( corners[2].x(), corners[2].y(), corners[2].z() );
143  pfrh->setNWCorner( corners[3].x(), corners[3].y(), corners[3].z() );
144 
145 
146  // if( !pfrh ) continue; // problem with this rechit. skip it
147 
148  rechits.push_back( *pfrh );
149  delete pfrh;
150  idSortedRecHits.insert( make_pair(detid.rawId(), rechits.size()-1 ) );
151  }
152  }
153 
154  // do navigation
155  for(unsigned i=0; i<rechits.size(); i++ ) {
156 
157  findRecHitNeighbours( rechits[i], idSortedRecHits,
158  psTopology,
159  *psGeometry,
160  psTopology,
161  *psGeometry);
162  }
163 }
164 
165 
166 
167 void
170  const map<unsigned,unsigned >& sortedHits,
171  const CaloSubdetectorTopology& barrelTopology,
172  const CaloSubdetectorGeometry& barrelGeometry,
173  const CaloSubdetectorTopology& endcapTopology,
174  const CaloSubdetectorGeometry& endcapGeometry ) {
175 
176  DetId detid( rh.detId() );
177 
178  const CaloSubdetectorTopology* topology = 0;
180 
181  switch( rh.layer() ) {
182  case PFLayer::ECAL_ENDCAP:
183  topology = &endcapTopology;
184  geometry = &endcapGeometry;
185  break;
186  case PFLayer::ECAL_BARREL:
187  topology = &barrelTopology;
188  geometry = &barrelGeometry;
189  break;
191  topology = &endcapTopology;
192  geometry = &endcapGeometry;
193  break;
195  topology = &barrelTopology;
196  geometry = &barrelGeometry;
197  break;
198  case PFLayer::PS1:
199  case PFLayer::PS2:
200  topology = &barrelTopology;
201  geometry = &barrelGeometry;
202  break;
203  default:
204  assert(0);
205  }
206 
207  assert( topology && geometry );
208 
209  CaloNavigator<DetId> navigator(detid, topology);
210 
211  DetId north = navigator.north();
212 
213  DetId northeast(0);
214  if( north != DetId(0) ) {
215  northeast = navigator.east();
216  }
217  navigator.home();
218 
219 
220  DetId south = navigator.south();
221 
222 
223 
224  DetId southwest(0);
225  if( south != DetId(0) ) {
226  southwest = navigator.west();
227  }
228  navigator.home();
229 
230 
231  DetId east = navigator.east();
232  DetId southeast;
233  if( east != DetId(0) ) {
234  southeast = navigator.south();
235  }
236  navigator.home();
237  DetId west = navigator.west();
238  DetId northwest;
239  if( west != DetId(0) ) {
240  northwest = navigator.north();
241  }
242  navigator.home();
243 
244  IDH i = sortedHits.find( north.rawId() );
245  if(i != sortedHits.end() )
246  rh.add4Neighbour( i->second );
247 
248  i = sortedHits.find( northeast.rawId() );
249  if(i != sortedHits.end() )
250  rh.add8Neighbour( i->second );
251 
252  i = sortedHits.find( south.rawId() );
253  if(i != sortedHits.end() )
254  rh.add4Neighbour( i->second );
255 
256  i = sortedHits.find( southwest.rawId() );
257  if(i != sortedHits.end() )
258  rh.add8Neighbour( i->second );
259 
260  i = sortedHits.find( east.rawId() );
261  if(i != sortedHits.end() )
262  rh.add4Neighbour( i->second );
263 
264  i = sortedHits.find( southeast.rawId() );
265  if(i != sortedHits.end() )
266  rh.add8Neighbour( i->second );
267 
268  i = sortedHits.find( west.rawId() );
269  if(i != sortedHits.end() )
270  rh.add4Neighbour( i->second );
271 
272  i = sortedHits.find( northwest.rawId() );
273  if(i != sortedHits.end() )
274  rh.add8Neighbour( i->second );
275 
276 
277 }
278 
#define LogDebug(id)
void setSECorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:226
T getParameter(std::string const &) const
void add4Neighbour(unsigned index)
Definition: PFRecHit.cc:176
int i
Definition: DBlmapReader.cc:9
void add8Neighbour(unsigned index)
Definition: PFRecHit.cc:181
void setNECorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:231
const DetId & detid() const
Definition: CaloRecHit.h:20
unsigned detId() const
rechit detId
Definition: PFRecHit.h:105
std::vector< EcalRecHit >::const_iterator const_iterator
T y() const
Definition: PV3DBase.h:63
float float float z
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
PFLayer::Layer layer() const
rechit layer
Definition: PFRecHit.h:108
void setSWCorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:221
void findRecHitNeighbours(reco::PFRecHit &rh, const std::map< unsigned, unsigned > &sortedHits, const CaloSubdetectorTopology &barrelTopo, const CaloSubdetectorGeometry &barrelGeom, const CaloSubdetectorTopology &endcapTopo, const CaloSubdetectorGeometry &endcapGeom)
int iEvent
Definition: GenABIO.cc:243
Particle flow rechit (rechit + geometry and topology information). See clustering algorithm in PFClus...
Definition: PFRecHit.h:31
Base producer for particle flow rechits (PFRecHit)
float energy() const
Definition: CaloRecHit.h:17
T z() const
Definition: PV3DBase.h:64
void createRecHits(std::vector< reco::PFRecHit > &rechits, std::vector< reco::PFRecHit > &rechitsCleaned, edm::Event &, const edm::EventSetup &)
edm::InputTag inputTagEcalRecHitsES_
std::vector< LinkConnSpec >::const_iterator IT
bool isValid() const
Definition: HandleBase.h:76
void setNWCorner(double posx, double posy, double posz)
search for pointers to neighbours, using neighbours&#39; DetId.
Definition: PFRecHit.cc:216
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
PFRecHitProducerPS(const edm::ParameterSet &)
Layer
layer definition
Definition: PFLayer.h:31
Definition: DetId.h:18
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: Handle.h:81
ESHandle< TrackerGeometry > geometry
std::map< unsigned, unsigned >::const_iterator IDH
int plane() const
Definition: ESDetId.h:46
Definition: DDAxes.h:10
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
T x() const
Definition: PV3DBase.h:62
virtual const CornersVec & getCorners() const =0
Returns the corner points of this cell&#39;s volume.