CMS 3D CMS Logo

Classes | Functions
CommonUsefulStuff.h File Reference
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/EDAnalyzer.h"
#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
#include "DataFormats/DetId/interface/DetId.h"
#include "DataFormats/HcalDetId/interface/HcalDetId.h"
#include "DataFormats/EcalDetId/interface/EEDetId.h"
#include "DataFormats/EcalDetId/interface/EBDetId.h"
#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
#include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
#include "Geometry/Records/interface/CaloGeometryRecord.h"
#include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
#include "FWCore/Utilities/interface/Exception.h"
#include "DataFormats/GeometryVector/interface/GlobalPoint.h"
#include "DataFormats/GeometryVector/interface/GlobalVector.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "CommonTools/UtilAlgos/interface/DeltaR.h"
#include "CommonTools/UtilAlgos/interface/TFileService.h"
#include "TROOT.h"
#include "TFile.h"
#include "TTree.h"
#include "TH1F.h"
#include "TString.h"

Go to the source code of this file.

Classes

struct  MaxHit_struct
 

Functions

double ecalEnergyInCone (const GlobalPoint center, double radius, const EcalRecHitCollection ecalCol, const CaloGeometry *geo)
 
double ecalEnergyInCone (const GlobalVector trackMom, const GlobalPoint center, double radius, const EcalRecHitCollection ecalCol, const CaloGeometry *geo)
 
double getDistInPlane (const GlobalVector trackDirection, const GlobalPoint caloPoint, const GlobalPoint rechitPoint, double coneHeight)
 
double getDistInPlaneSimple (const GlobalPoint caloPoint, const GlobalPoint rechitPoint)
 
double getDistInPlaneTrackDir (const GlobalPoint caloPoint, const GlobalVector caloVector, const GlobalPoint rechitPoint)
 

Function Documentation

◆ ecalEnergyInCone() [1/2]

double ecalEnergyInCone ( const GlobalPoint  center,
double  radius,
const EcalRecHitCollection  ecalCol,
const CaloGeometry geo 
)
inline

Definition at line 176 of file CommonUsefulStuff.h.

179  {
180  double eECALcone = 0;
181  std::vector<int> usedHitsEcal;
182  usedHitsEcal.clear();
183 
184  for (std::vector<EcalRecHit>::const_iterator ehit = ecalCol.begin(); ehit != ecalCol.end(); ehit++) {
185  //This is a precaution for the case when hitCollection contains duplicats.
186  bool hitIsUsed = false;
187  int hitHashedIndex = -10000;
188  if (ehit->id().subdetId() == EcalBarrel) {
189  EBDetId did(ehit->id());
190  hitHashedIndex = did.hashedIndex();
191  }
192  if (ehit->id().subdetId() == EcalEndcap) {
193  EEDetId did(ehit->id());
194  hitHashedIndex = did.hashedIndex();
195  }
196  for (uint32_t i = 0; i < usedHitsEcal.size(); i++) {
197  if (usedHitsEcal[i] == hitHashedIndex)
198  hitIsUsed = true;
199  }
200  if (hitIsUsed)
201  continue;
202  usedHitsEcal.push_back(hitHashedIndex);
203  // -----------------------------------------------
204 
205  const GlobalPoint& pos = geo->getPosition((*ehit).detid());
206 
207  if (getDistInPlaneSimple(center, pos) < radius) {
208  eECALcone += ehit->energy();
209  }
210  }
211  return eECALcone;
212 }

References edm::SortedCollection< T, SORT >::begin(), EcalBarrel, EcalEndcap, edm::SortedCollection< T, SORT >::end(), getDistInPlaneSimple(), CaloGeometry::getPosition(), EBDetId::hashedIndex(), EEDetId::hashedIndex(), mps_fire::i, and CosmicsPD_Skims::radius.

◆ ecalEnergyInCone() [2/2]

double ecalEnergyInCone ( const GlobalVector  trackMom,
const GlobalPoint  center,
double  radius,
const EcalRecHitCollection  ecalCol,
const CaloGeometry geo 
)
inline

Definition at line 215 of file CommonUsefulStuff.h.

219  {
220  double eECALcone = 0;
221  std::vector<int> usedHitsEcal;
222  usedHitsEcal.clear();
223  for (std::vector<EcalRecHit>::const_iterator ehit = ecalCol.begin(); ehit != ecalCol.end(); ehit++) {
224  //This is a precaution for the case when hitCollection contains duplicats.
225  bool hitIsUsed = false;
226  int hitHashedIndex = -10000;
227  if (ehit->id().subdetId() == EcalBarrel) {
228  EBDetId did(ehit->id());
229  hitHashedIndex = did.hashedIndex();
230  }
231  if (ehit->id().subdetId() == EcalEndcap) {
232  EEDetId did(ehit->id());
233  hitHashedIndex = did.hashedIndex();
234  }
235  for (uint32_t i = 0; i < usedHitsEcal.size(); i++) {
236  if (usedHitsEcal[i] == hitHashedIndex)
237  hitIsUsed = true;
238  }
239  if (hitIsUsed)
240  continue;
241  usedHitsEcal.push_back(hitHashedIndex);
242  // -----------------------------------------------
243 
244  const GlobalPoint& pos = geo->getPosition((*ehit).detid());
245 
246  if (getDistInPlaneTrackDir(center, trackMom, pos) < radius) {
247  eECALcone += ehit->energy();
248  }
249  }
250  return eECALcone;
251 }

References edm::SortedCollection< T, SORT >::begin(), EcalBarrel, EcalEndcap, edm::SortedCollection< T, SORT >::end(), getDistInPlaneTrackDir(), CaloGeometry::getPosition(), EBDetId::hashedIndex(), EEDetId::hashedIndex(), mps_fire::i, and CosmicsPD_Skims::radius.

◆ getDistInPlane()

double getDistInPlane ( const GlobalVector  trackDirection,
const GlobalPoint  caloPoint,
const GlobalPoint  rechitPoint,
double  coneHeight 
)
inline

Definition at line 106 of file CommonUsefulStuff.h.

109  {
110  // The iso track candidate hits the Calo (Ecal or Hcal) at "caloPoint"
111  // with direction "trackDirection".
112 
113  // "rechitPoint" is the position of the rechit. We only care about
114  // the direction of the rechit.
115 
116  // Consider the rechitPoint as characterized by angles theta and phi
117  // wrt the origin which points at the calo cell of the rechit. In
118  // some sense the distance along the line theta/phi is arbitrary. A
119  // simplified choice might be to put the rechit at the surface of
120  // the Hcal. Our approach is to see whether/where this line
121  // intersects a cone of height "coneHeight" with vertex at caloPoint
122  // aligned with trackDirection.
123  // To that end, this function returns the distance between the
124  // center of the base of the cone and the intersection of the rechit
125  // line and the plane that contains the base of the cone. This
126  // distance is compared with the radius of the cone outside this
127  // function.
128 
129  // Make vector of length cone height along track direction
130  const GlobalVector heightVector = trackDirection * coneHeight;
131 
132  // Make vector from origin to point where iso track intersects the
133  // calorimeter.
134  const GlobalVector caloIntersectVector(caloPoint.x(), caloPoint.y(), caloPoint.z());
135 
136  // Make vector from origin to point in center of base of cone
137  const GlobalVector coneBaseVector = heightVector + caloIntersectVector;
138 
139  // Make point in center of base of cone
140  const GlobalPoint coneBasePoint(coneBaseVector.x(), coneBaseVector.y(), coneBaseVector.z());
141 
142  // Make unit vector pointing at rechit.
143  const GlobalVector rechitVector(rechitPoint.x(), rechitPoint.y(), rechitPoint.z());
144  const GlobalVector rechitDirection = rechitVector.unit();
145 
146  // Find distance "r" along "rechit line" (with angles theta2 and
147  // phi2) where line intersects plane defined by base of cone.
148 
149  // Definition plane of that contains base of cone:
150  // trackDirection.x() (x - coneBaseVector.x()) +
151  // trackDirection.y() (y - coneBaseVector.y()) +
152  // trackDirection.z() (z - coneBaseVector.z()) = 0
153 
154  // Point P_{rh} where rechit line intersects plane:
155  // (rechitdist sin(theta2) cos(phi2),
156  // rechitdist sin(theta2) cos(phi2),
157  // rechitdist cos(theta2))
158 
159  // Substitute P_{rh} into equation for plane and solve for rechitdist.
160  // rechitDist turns out to be the ratio of dot products:
161 
162  double rechitdist = trackDirection.dot(coneBaseVector) / trackDirection.dot(rechitDirection);
163 
164  // Now find distance between point at center of cone base and point
165  // where the "rechit line" intersects the plane defined by the base
166  // of the cone; i.e. the effectiveRecHitPoint.
167  const GlobalVector effectiveRechitVector = rechitdist * rechitDirection;
168  const GlobalPoint effectiveRechitPoint(
169  effectiveRechitVector.x(), effectiveRechitVector.y(), effectiveRechitVector.z());
170 
171  GlobalVector distance_vector = effectiveRechitPoint - coneBasePoint;
172  return distance_vector.mag();
173 }

References Vector3DBase< T, FrameTag >::dot(), PV3DBase< T, PVType, FrameType >::mag(), Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

◆ getDistInPlaneSimple()

double getDistInPlaneSimple ( const GlobalPoint  caloPoint,
const GlobalPoint  rechitPoint 
)
inline

Definition at line 53 of file CommonUsefulStuff.h.

53  {
54  // Simplified version of getDistInPlane
55  // Assume track direction is origin -> point of hcal intersection
56 
57  const GlobalVector caloIntersectVector(caloPoint.x(), caloPoint.y(), caloPoint.z());
58 
59  const GlobalVector caloIntersectUnitVector = caloIntersectVector.unit();
60 
61  const GlobalVector rechitVector(rechitPoint.x(), rechitPoint.y(), rechitPoint.z());
62 
63  const GlobalVector rechitUnitVector = rechitVector.unit();
64  double dotprod = caloIntersectUnitVector.dot(rechitUnitVector);
65  double rechitdist = caloIntersectVector.mag() / dotprod;
66 
67  const GlobalVector effectiveRechitVector = rechitdist * rechitUnitVector;
68  const GlobalPoint effectiveRechitPoint(
69  effectiveRechitVector.x(), effectiveRechitVector.y(), effectiveRechitVector.z());
70 
71  GlobalVector distance_vector = effectiveRechitPoint - caloPoint;
72 
73  if (dotprod > 0.) {
74  return distance_vector.mag();
75  } else {
76  return 999999.;
77  }
78 }

References Vector3DBase< T, FrameTag >::dot(), PV3DBase< T, PVType, FrameType >::mag(), Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by ecalEnergyInCone().

◆ getDistInPlaneTrackDir()

double getDistInPlaneTrackDir ( const GlobalPoint  caloPoint,
const GlobalVector  caloVector,
const GlobalPoint  rechitPoint 
)
inline

Definition at line 80 of file CommonUsefulStuff.h.

82  {
83  // Simplified version of getDistInPlane : no cone "within" Hcal, but
84  // don't assume track direction is origin -> point of hcal
85  // intersection.
86  const GlobalVector caloIntersectVector(caloPoint.x(), caloPoint.y(),
87  caloPoint.z()); //p
88 
89  const GlobalVector caloUnitVector = caloVector.unit();
90  const GlobalVector rechitVector(rechitPoint.x(), rechitPoint.y(), rechitPoint.z());
91  const GlobalVector rechitUnitVector = rechitVector.unit();
92  double dotprod_denominator = caloUnitVector.dot(rechitUnitVector);
93  double dotprod_numerator = caloUnitVector.dot(caloIntersectVector);
94  double rechitdist = dotprod_numerator / dotprod_denominator;
95  const GlobalVector effectiveRechitVector = rechitdist * rechitUnitVector;
96  const GlobalPoint effectiveRechitPoint(
97  effectiveRechitVector.x(), effectiveRechitVector.y(), effectiveRechitVector.z());
98  GlobalVector distance_vector = effectiveRechitPoint - caloPoint;
99  if (dotprod_denominator > 0. && dotprod_numerator > 0.) {
100  return distance_vector.mag();
101  } else {
102  return 999999.;
103  }
104 }

References Vector3DBase< T, FrameTag >::dot(), PV3DBase< T, PVType, FrameType >::mag(), Vector3DBase< T, FrameTag >::unit(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by ecalEnergyInCone().

Vector3DBase
Definition: Vector3DBase.h:8
EBDetId::hashedIndex
int hashedIndex() const
get a compact index for arrays
Definition: EBDetId.h:82
mps_fire.i
i
Definition: mps_fire.py:355
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
CaloGeometry::getPosition
GlobalPoint getPosition(const DetId &id) const
Get the position of a given detector id.
Definition: CaloGeometry.cc:50
EBDetId
Definition: EBDetId.h:17
pos
Definition: PixelAliasList.h:18
EcalBarrel
Definition: EcalSubdetector.h:10
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
Vector3DBase::unit
Vector3DBase unit() const
Definition: Vector3DBase.h:54
edm::SortedCollection::begin
const_iterator begin() const
Definition: SortedCollection.h:262
Point3DBase< float, GlobalTag >
EEDetId
Definition: EEDetId.h:14
EcalEndcap
Definition: EcalSubdetector.h:10
getDistInPlaneSimple
double getDistInPlaneSimple(const GlobalPoint caloPoint, const GlobalPoint rechitPoint)
Definition: CommonUsefulStuff.h:53
edm::SortedCollection::end
const_iterator end() const
Definition: SortedCollection.h:267
EEDetId::hashedIndex
int hashedIndex() const
Definition: EEDetId.h:183
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
Vector3DBase::dot
PreciseFloatType< T, U >::Type dot(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:99
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:64
CosmicsPD_Skims.radius
radius
Definition: CosmicsPD_Skims.py:135
getDistInPlaneTrackDir
double getDistInPlaneTrackDir(const GlobalPoint caloPoint, const GlobalVector caloVector, const GlobalPoint rechitPoint)
Definition: CommonUsefulStuff.h:80