CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Classes | Public Member Functions | Private Types | Private Member Functions | Private Attributes | Friends
DTDigitizer Class Reference

#include <DTDigitizer.h>

Inheritance diagram for DTDigitizer:
edm::stream::EDProducer<> edm::stream::EDProducerBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Classes

struct  hitLessT
 

Public Member Functions

 DTDigitizer (const edm::ParameterSet &)
 
virtual void produce (edm::Event &, const edm::EventSetup &) override
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
- Public Member Functions inherited from edm::stream::EDProducerBase
 EDProducerBase ()
 
ModuleDescription const & moduleDescription () const
 
virtual ~EDProducerBase ()
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 ProducerBase ()
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription
const &)> 
registrationCallback () const
 used by the fwk to register list of products More...
 
virtual ~ProducerBase ()
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
 EDConsumerBase ()
 
ProductHolderIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
std::vector
< ProductHolderIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesDependentUpon (std::string const &iProcessName, std::string const &iModuleLabel, bool iPrint, std::vector< char const * > &oModuleLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
bool registeredToConsume (ProductHolderIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Private Types

typedef std::map< DTWireId,
std::vector< const PSimHit * > > 
DTWireIdMap
 
typedef DTWireIdMap::const_iterator DTWireIdMapConstIter
 
typedef DTWireIdMap::iterator DTWireIdMapIter
 
typedef std::pair< const
PSimHit *, float > 
hitAndT
 
typedef std::vector< hitAndTTDContainer
 

Private Member Functions

float asymGausSmear (double mean, double sigmaLeft, double sigmaRight, CLHEP::HepRandomEngine *) const
 
std::pair< float, bool > computeTime (const DTLayer *layer, const DTWireId &wireId, const PSimHit *hit, const LocalVector &BLoc, CLHEP::HepRandomEngine *)
 
std::pair< float, bool > driftTimeFromParametrization (float x, float alpha, float By, float Bz, CLHEP::HepRandomEngine *) const
 
std::pair< float, bool > driftTimeFromTimeMap () const
 
void dumpHit (const PSimHit *hit, float xEntry, float xExit, const DTTopology &topo)
 
float externalDelays (const DTLayer *layer, const DTWireId &wireId, const PSimHit *hit) const
 
void storeDigis (DTWireId &wireId, TDContainer &hits, DTDigiCollection &output, DTDigiSimLinkCollection &outputLinks)
 

Private Attributes

edm::EDGetTokenT
< CrossingFrame< PSimHit > > 
cf_token
 
std::string collection_for_XF
 
float deadTime
 
bool debug
 
std::string geometryType
 
bool IdealModel
 
bool interpolate
 
float LinksTimeWindow
 
std::string mix_
 
bool MultipleLinks
 
bool onlyMuHits
 
float smearing
 
std::string syncName
 
float theConstVDrift
 
std::unique_ptr< DTDigiSyncBasetheSync
 
double vPropWire
 

Friends

class DTDigitizerAnalysis
 

Additional Inherited Members

- Public Types inherited from edm::stream::EDProducer<>
typedef CacheContexts< T...> CacheTypes
 
typedef CacheTypes::GlobalCache GlobalCache
 
typedef AbilityChecker< T...> HasAbility
 
typedef
CacheTypes::LuminosityBlockCache 
LuminosityBlockCache
 
typedef
LuminosityBlockContextT
< LuminosityBlockCache,
RunCache, GlobalCache
LuminosityBlockContext
 
typedef
CacheTypes::LuminosityBlockSummaryCache 
LuminosityBlockSummaryCache
 
typedef CacheTypes::RunCache RunCache
 
typedef RunContextT< RunCache,
GlobalCache
RunContext
 
typedef CacheTypes::RunSummaryCache RunSummaryCache
 
- Public Types inherited from edm::stream::EDProducerBase
typedef EDProducerAdaptorBase ModuleType
 
- Public Types inherited from edm::ProducerBase
typedef
ProductRegistryHelper::TypeLabelList 
TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::stream::EDProducerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

Digitize the muon drift tubes. The parametrisation function in DTDriftTimeParametrization from P.G.Abia, J.Puerta is used in all cases where it is applicable.

Authors
: G. Bevilacqua, N. Amapane, G. Cerminara, R. Bellan

Definition at line 45 of file DTDigitizer.h.

Member Typedef Documentation

typedef std::map<DTWireId, std::vector<const PSimHit*> > DTDigitizer::DTWireIdMap
private

Definition at line 57 of file DTDigitizer.h.

typedef DTWireIdMap::const_iterator DTDigitizer::DTWireIdMapConstIter
private

Definition at line 59 of file DTDigitizer.h.

typedef DTWireIdMap::iterator DTDigitizer::DTWireIdMapIter
private

Definition at line 58 of file DTDigitizer.h.

typedef std::pair<const PSimHit*,float> DTDigitizer::hitAndT
private

Definition at line 54 of file DTDigitizer.h.

typedef std::vector<hitAndT> DTDigitizer::TDContainer
private

Definition at line 55 of file DTDigitizer.h.

Constructor & Destructor Documentation

DTDigitizer::DTDigitizer ( const edm::ParameterSet conf_)
explicit

Definition at line 57 of file DTDigitizer.cc.

References SurfaceDeformationFactory::create(), debug, Exception, reco::get(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), HLT_FULL_cff::InputTag, edm::Service< T >::isAvailable(), and AlCaHLTBitMon_QueryRunRegistry::string.

57  {
58 
59  // Set verbose output
60  debug=conf_.getUntrackedParameter<bool>("debug");
61 
62  if (debug) LogPrint("DTDigitizer")<<"Creating a DTDigitizer"<<endl;
63 
64  //register the Producer with a label
65  //produces<DTDigiCollection>("MuonDTDigis"); // FIXME: Do I pass it by ParameterSet?
66  produces<DTDigiCollection>(); // FIXME: Do I pass it by ParameterSet?
67  // produces<DTDigiSimLinkCollection>("MuonDTDigiSimLinks");
68  produces<DTDigiSimLinkCollection>();
69 
70  //Parameters:
71 
72  // build digis only for mu hits (for debug purposes)
73  onlyMuHits=conf_.getParameter<bool>("onlyMuHits");
74 
75  // interpolate parametrization function
76  interpolate=conf_.getParameter<bool>("interpolate");
77 
78  // Velocity of signal propagation along the wire (cm/ns)
79  // For the default value
80  // cfr. CMS-IN 2000-021: (2.56+-0.17)x1e8 m/s
81  // CMS NOTE 2003-17: (0.244) m/ns
82  vPropWire=conf_.getParameter<double>("vPropWire"); //24.4
83 
84  // Dead time for signals on the same wire (number from M. Pegoraro)
85  deadTime=conf_.getParameter<double>("deadTime"); //150
86 
87  // further configurable smearing
88  smearing=conf_.getParameter<double>("Smearing"); // 3.
89 
90  // Sync Algo
91  syncName = conf_.getParameter<string>("SyncName");
92  theSync.reset( DTDigiSyncFactory::get()->create(syncName,conf_.getParameter<ParameterSet>("pset")) );
93 
94  // Debug flag to switch to the Ideal model
95  // it uses a constant drift velocity and doesn't set any external delay
96  IdealModel = conf_.getParameter<bool>("IdealModel");
97 
98  // Constant drift velocity needed by the above flag
99  if(IdealModel)
100  theConstVDrift = conf_.getParameter<double>("IdealModelConstantDriftVelocity"); // 55 um/ns
101  else theConstVDrift = 55.;
102 
103  // get random engine
105  if ( ! rng.isAvailable()) {
106  throw cms::Exception("Configuration")
107  << "RandomNumberGeneratorService for DTDigitizer missing in cfg file";
108  }
109 
110  // MultipleLinks=false ==> one-to-one correspondence between digis and SimHits
111  MultipleLinks = conf_.getParameter<bool>("MultipleLinks");
112  // MultipleLinks=true ==> association of SimHits within a time window LinksTimeWindow
113  // (of the order of the resolution)
114  LinksTimeWindow = conf_.getParameter<double>("LinksTimeWindow"); // (10 ns)
115 
116  //Name of Collection used for create the XF
117  mix_ = conf_.getParameter<std::string>("mixLabel");
118  collection_for_XF = conf_.getParameter<std::string>("InputCollection");
119  cf_token = consumes<CrossingFrame<PSimHit> >( edm::InputTag(mix_, collection_for_XF) );
120 
121  //String to choice between ideal (the deafult) and (mis)aligned geometry for the digitization step
122  geometryType = conf_.getParameter<std::string>("GeometryType");
123 
124 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
bool onlyMuHits
Definition: DTDigitizer.h:113
float theConstVDrift
Definition: DTDigitizer.h:122
float LinksTimeWindow
Definition: DTDigitizer.h:126
bool isAvailable() const
Definition: Service.h:46
edm::EDGetTokenT< CrossingFrame< PSimHit > > cf_token
Definition: DTDigitizer.h:132
float smearing
Definition: DTDigitizer.h:110
bool interpolate
Definition: DTDigitizer.h:112
bool IdealModel
Definition: DTDigitizer.h:121
std::unique_ptr< DTDigiSyncBase > theSync
Definition: DTDigitizer.h:116
std::string mix_
Definition: DTDigitizer.h:129
std::string geometryType
Definition: DTDigitizer.h:118
std::string syncName
Definition: DTDigitizer.h:115
double vPropWire
Definition: DTDigitizer.h:108
float deadTime
Definition: DTDigitizer.h:109
SurfaceDeformation * create(int type, const std::vector< double > &params)
T get(const Candidate &c)
Definition: component.h:55
bool MultipleLinks
Definition: DTDigitizer.h:125
std::string collection_for_XF
Definition: DTDigitizer.h:130

Member Function Documentation

float DTDigitizer::asymGausSmear ( double  mean,
double  sigmaLeft,
double  sigmaRight,
CLHEP::HepRandomEngine *  engine 
) const
private

Definition at line 442 of file DTDigitizer.cc.

References f, and lumiQTWidget::t.

443  {
444 
445  double f = sigmaLeft/(sigmaLeft+sigmaRight);
446  double t;
447 
448  if (CLHEP::RandFlat::shoot(engine) <= f) {
449  t = CLHEP::RandGaussQ::shoot(engine, mean, sigmaLeft);
450  t = mean - fabs(t - mean);
451  } else {
452  t = CLHEP::RandGaussQ::shoot(engine, mean, sigmaRight);
453  t = mean + fabs(t - mean);
454  }
455  return static_cast<float>(t);
456 }
double f[11][100]
pair< float, bool > DTDigitizer::computeTime ( const DTLayer layer,
const DTWireId wireId,
const PSimHit hit,
const LocalVector BLoc,
CLHEP::HepRandomEngine *  engine 
)
private

Definition at line 229 of file DTDigitizer.cc.

References Vector3DBase< T, FrameTag >::cross(), ztail::d, debug, delta, dir, Vector3DBase< T, FrameTag >::dot(), PSimHit::entryPoint(), PSimHit::exitPoint(), PSimHit::localDirection(), PSimHit::localPosition(), M_PI, PV3DBase< T, PVType, FrameType >::mag(), PSimHit::momentumAtEntry(), DTTopology::none, DTTopology::onWhichBorder(), PSimHit::pabs(), PSimHit::particleType(), EnergyCorrector::pt, DTLayer::specificTopology(), mathSSE::sqrt(), theta(), Vector3DBase< T, FrameTag >::unit(), DTWireId::wire(), DTTopology::wirePosition(), PV3DBase< T, PVType, FrameType >::x(), x(), DTTopology::xMax, DTTopology::xMin, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

231  {
232 
233  LocalPoint entryP = hit->entryPoint();
234  LocalPoint exitP = hit->exitPoint();
235  int partType = hit->particleType();
236 
237  const DTTopology &topo = layer->specificTopology();
238 
239  // Pay attention: in CMSSW the rf of the SimHit is in the layer's rf
240 
241  if(debug) LogPrint("DTDigitizer")<<"Hit local entry point: "<<entryP<<endl
242  <<"Hit local exit point: "<<exitP<<endl;
243 
244  float xwire = topo.wirePosition(wireId.wire());
245  float xEntry = entryP.x() - xwire;
246  float xExit = exitP.x() - xwire;
247 
248  if(debug) LogPrint("DTDigitizer")<<"wire position: "<<xwire
249  <<" x entry in cell rf: "<<xEntry
250  <<" x exit in cell rf: "<<xExit<<endl;
251 
252  DTTopology::Side entrySide = topo.onWhichBorder(xEntry,entryP.y(),entryP.z());
253  DTTopology::Side exitSide = topo.onWhichBorder(xExit,exitP.y(),exitP.z());
254 
255  if (debug) dumpHit(hit, xEntry, xExit,topo);
256 
257  // The bolean is used to flag the drift time computation
258  pair<float,bool> driftTime(0.,false);
259 
260  // if delta in gas->ignore, since it is included in the parametrisation.
261  // FIXME: should check that it is actually a delta ray produced by a nearby
262  // muon hit.
263 
264  if (partType == 11 && entrySide == DTTopology::none) {
265  if (debug) LogPrint("DTDigitizer") << " e- hit in gas; discarding " << endl;
266  return driftTime;
267  }
268 
269  float By = BLoc.y();
270  float Bz = BLoc.z();
271 
272  // Radius and sagitta according to direction of momentum
273  // (just for printing)
274  // NOTE: in cmsim, d is always taken // pHat!
275  LocalVector d = (exitP-entryP);
276  LocalVector pHat = hit->localDirection().unit();
277  LocalVector hHat = (d.cross(pHat.cross(d))).unit();
278  float cosAlpha = hHat.dot(pHat);
279  float sinAlpha = sqrt(1.-cosAlpha*cosAlpha);
280  float radius_P = (d.mag())/(2.*cosAlpha);
281  float sagitta_P = radius_P*(1.-sinAlpha);
282 
283  // Radius, sagitta according to field bending
284  // (just for printing)
285  float halfd = d.mag()/2.;
286  float BMag = BLoc.mag();
287  LocalVector pT = (pHat - (BLoc.unit()*pHat.dot(BLoc.unit())))*(hit->pabs());
288  float radius_B = (pT.mag()/(0.3*BMag))*100.;
289  float sagitta_B;
290  if (radius_B > halfd) {
291  sagitta_B = radius_B - sqrt(radius_B*radius_B - halfd*halfd);
292  } else {
293  sagitta_B = radius_B;
294  }
295 
296  // cos(delta), delta= angle between direction at entry and hit segment
297  // (just for printing)
298  float delta = pHat.dot(d.unit());
299  if (debug) LogPrint("DTDigitizer") << " delta = " << delta << endl
300  << " cosAlpha = " << cosAlpha << endl
301  << " sinAlpha = " << sinAlpha << endl
302  << " pMag = " << pT.mag() << endl
303  << " bMag = " << BMag << endl
304  << " pT = " << pT << endl
305  << " halfd = " << halfd << endl
306  << " radius_P (cm) = " << radius_P << endl
307  << " sagitta_P (um) = " << sagitta_P*10000. << endl
308  << " radius_B (cm) = " << radius_B << endl
309  << " sagitta_B (um) = " << sagitta_B*10000. << endl;
310 
311  // Select cases where parametrization can not be used.
312  bool noParametrisation =
313  ( ( entrySide == DTTopology::none || exitSide == DTTopology::none ) // case # 2,3,8,9 or 11
314  || (entrySide == exitSide) // case # 4 or 10
315  || ((entrySide == DTTopology::xMin && exitSide == DTTopology::xMax) ||
316  (entrySide == DTTopology::xMax && exitSide == DTTopology::xMin)) // Hit is case # 7
317  );
318 
319  // FIXME: now, debug warning only; consider treating those
320  // with TM algo.
321  if ( delta < 0.99996 // Track is not straight. FIXME: use sagitta?
322  && (noParametrisation == false)) {
323  if (debug) LogPrint("DTDigitizer") << "*** WARNING: hit is not straight, type = " << partType << endl;
324  }
325 
326  //************ 5A ***************
327 
328  if (!noParametrisation) {
329 
330  LocalVector dir = hit->momentumAtEntry(); // ex Measurement3DVector dir = hit->measurementDirection(); //FIXME?
331  float theta = atan(dir.x()/-dir.z())*180/M_PI;
332 
333  // FIXME: use dir if M.S. is included as GARFIELD option...
334  // otherwise use hit segment dirction.
335  // LocalVector dir0 = (exitP-entryP).unit();
336  // float theta = atan(dir0.x()/-dir0.z())*180/M_PI;
337  float x;
338 
339  Local3DPoint pt = hit->localPosition(); //ex Measurement3DPoint pt = hit->measurementPosition(); // FIXME?
340 
341  if(fabs(pt.z()) < 0.002) {
342  // hit center within 20 um from z=0, no need to extrapolate.
343  x = pt.x() - xwire;
344  } else {
345  x = xEntry - (entryP.z()*(xExit-xEntry))/(exitP.z()-entryP.z());
346  }
347 
348  if(IdealModel) return make_pair(fabs(x)/theConstVDrift,true);
349  else driftTime = driftTimeFromParametrization(x, theta, By, Bz, engine);
350 
351  }
352 
353 
354  if ((driftTime.second)==false) {
355  // Parametrisation not applicable, or failed. Use time map.
356  driftTime = driftTimeFromTimeMap();
357  }
358 
359  //************ 5B ***************
360 
361  // Signal propagation, TOF etc.
362  if (driftTime.second) {
363  driftTime.first += externalDelays(layer,wireId,hit);
364  }
365  return driftTime;
366 }
dbl * delta
Definition: mlp_gen.cc:36
void dumpHit(const PSimHit *hit, float xEntry, float xExit, const DTTopology &topo)
Definition: DTDigitizer.cc:573
float wirePosition(int wireNumber) const
Returns the x position in the layer of a given wire number.
Definition: DTTopology.cc:86
Side onWhichBorder(float x, float y, float z) const
Definition: DTTopology.cc:109
LocalVector momentumAtEntry() const
The momentum of the track that produced the hit, at entry point.
Definition: PSimHit.h:47
float theConstVDrift
Definition: DTDigitizer.h:122
Geom::Theta< T > theta() const
T y() const
Definition: PV3DBase.h:63
PreciseFloatType< T, U >::Type dot(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:107
tuple d
Definition: ztail.py:151
const DTTopology & specificTopology() const
Definition: DTLayer.cc:42
Local3DPoint exitPoint() const
Exit point in the local Det frame.
Definition: PSimHit.h:38
T mag() const
Definition: PV3DBase.h:67
Local3DPoint localPosition() const
Definition: PSimHit.h:44
T sqrt(T t)
Definition: SSEVec.h:18
Vector3DBase< typename PreciseFloatType< T, U >::Type, FrameTag > cross(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:119
T z() const
Definition: PV3DBase.h:64
float pabs() const
fast and more accurate access to momentumAtEntry().mag()
Definition: PSimHit.h:63
Side
Sides of the cell.
Definition: DTTopology.h:88
bool IdealModel
Definition: DTDigitizer.h:121
#define M_PI
LocalVector localDirection() const
Obsolete. Same as momentumAtEntry().unit(), for backward compatibility.
Definition: PSimHit.h:52
int wire() const
Return the wire number.
Definition: DTWireId.h:56
Vector3DBase unit() const
Definition: Vector3DBase.h:57
std::pair< float, bool > driftTimeFromParametrization(float x, float alpha, float By, float Bz, CLHEP::HepRandomEngine *) const
Definition: DTDigitizer.cc:370
int particleType() const
Definition: PSimHit.h:85
float externalDelays(const DTLayer *layer, const DTWireId &wireId, const PSimHit *hit) const
Definition: DTDigitizer.cc:466
dbl *** dir
Definition: mlp_gen.cc:35
T x() const
Definition: PV3DBase.h:62
std::pair< float, bool > driftTimeFromTimeMap() const
Definition: DTDigitizer.cc:458
Local3DPoint entryPoint() const
Entry point in the local Det frame.
Definition: PSimHit.h:35
pair< float, bool > DTDigitizer::driftTimeFromParametrization ( float  x,
float  alpha,
float  By,
float  Bz,
CLHEP::HepRandomEngine *  engine 
) const
private

Definition at line 370 of file DTDigitizer.cc.

References debug, GeomDetEnumerators::DT, f, edm::max(), DTDriftTimeParametrization::MB_DT_drift_time(), DTDriftTimeParametrization::drift_time::t_drift, DTDriftTimeParametrization::drift_time::t_width_m, DTDriftTimeParametrization::drift_time::t_width_p, and theta().

371  {
372 
373  // Convert from CMSSW frame/units r.f. to parametrization ones.
374  x *= 10.; //cm -> mm
375 
376  // FIXME: Current parametrisation can extrapolate above 21 mm,
377  // however a detailed study is needed before using this.
378  if (fabs(x) > 21.) {
379  if (debug) LogPrint("DTDigitizer") << "*** WARNING: parametrisation: x out of range = "
380  << x << ", skipping" << endl;
381  return pair<float,bool>(0.f,false);
382  }
383 
384  // Different r.f. of the parametrization:
385  // X_par = X_ORCA; Y_par=Z_ORCA; Z_par = -Y_ORCA
386 
387  float By_par = Bz; // Bnorm
388  float Bz_par = -By; // Bwire
389  float theta_par = theta;
390 
391  // Parametrisation uses interpolation up to |theta|=45 deg,
392  // |Bwire|=0.4, |Bnorm|=0.75; extrapolation above.
393  if (fabs(theta_par)>45.) {
394  if (debug) LogPrint("DTDigitizer") << "*** WARNING: extrapolating theta > 45: "
395  << theta << endl;
396  // theta_par = min(fabs(theta_par),45.f)*((theta_par<0.)?-1.:1.);
397  }
398  if (fabs(By_par)>0.75) {
399  if (debug) LogPrint("DTDigitizer") << "*** WARNING: extrapolating Bnorm > 0.75: "
400  << By_par << endl;
401  // By_par = min(fabs(By_par),0.75f)*((By_par<0.)?-1.:1.);
402  }
403  if (fabs(Bz_par)>0.4) {
404  if (debug) LogPrint("DTDigitizer") << "*** WARNING: extrapolating Bwire >0.4: "
405  << Bz_par << endl;
406  // Bz_par = min(fabs(Bz_par),0.4)*((Bz_par<0.)?-1.:1.);
407  }
408 
410  static const DTDriftTimeParametrization par;
411  unsigned short flag = par.MB_DT_drift_time (x, theta_par, By_par, Bz_par, 0, &DT, interpolate);
412 
413  if (debug) {
414  LogPrint("DTDigitizer") << " Parametrisation: x, theta, Bnorm, Bwire = "
415  << x << " " << theta_par << " " << By_par << " " << Bz_par << endl
416  << " time=" << DT.t_drift
417  << " sigma_m=" << DT.t_width_m
418  << " sigma_p=" << DT.t_width_p << endl;
419  if (flag!=1) {
420  LogPrint("DTDigitizer") << "*** WARNING: call to parametrisation failed" << endl;
421  return pair<float,bool>(0.f,false);
422  }
423  }
424 
425  // Double half-gaussian smearing
426  float time = asymGausSmear(DT.t_drift, DT.t_width_m, DT.t_width_p, engine);
427 
428  // Do not allow the smearing to lead to negative values
429  time = max(time,0.f);
430 
431  // Apply a Gaussian smearing to account for electronic effects (cf. 2004 TB analysis)
432  // The width of the Gaussian can be configured with the "Smearing" parameter
433 
434  double u = CLHEP::RandGaussQ::shoot(engine, 0., smearing);
435  time += u;
436 
437  if (debug) LogPrint("DTDigitizer") << " drift time = " << time << endl;
438 
439  return pair<float,bool>(time,true);
440 }
Geom::Theta< T > theta() const
double f[11][100]
float smearing
Definition: DTDigitizer.h:110
bool interpolate
Definition: DTDigitizer.h:112
unsigned short MB_DT_drift_time(double x, double alpha, double by, double bz, short ifl, drift_time *DT, short interpolate) const
Calculate drift time and spread.
Structure used to return output values.
float asymGausSmear(double mean, double sigmaLeft, double sigmaRight, CLHEP::HepRandomEngine *) const
Definition: DTDigitizer.cc:442
pair< float, bool > DTDigitizer::driftTimeFromTimeMap ( ) const
private

Definition at line 458 of file DTDigitizer.cc.

References debug.

458  {
459  // FIXME: not yet implemented.
460  if (debug) LogPrint("DTDigitizer") << " TimeMap " << endl;
461  return pair<float,bool>(0.,false);
462 }
void DTDigitizer::dumpHit ( const PSimHit hit,
float  xEntry,
float  xExit,
const DTTopology topo 
)
private

Definition at line 573 of file DTDigitizer.cc.

References DTTopology::cellHeight(), DTTopology::cellLenght(), DTTopology::cellWidth(), PSimHit::energyLoss(), PSimHit::entryPoint(), PSimHit::exitPoint(), mag(), PSimHit::momentumAtEntry(), DTTopology::onWhichBorder(), PSimHit::pabs(), PSimHit::particleType(), PSimHit::processType(), PSimHit::trackId(), Vector3DBase< T, FrameTag >::unit(), x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

575  {
576 
577  LocalPoint entryP = hit->entryPoint();
578  LocalPoint exitP = hit->exitPoint();
579 
580  DTTopology::Side entrySide = topo.onWhichBorder(xEntry,entryP.y(),entryP.z());
581  DTTopology::Side exitSide = topo.onWhichBorder(xExit,exitP.y(),exitP.z());
582  // ProcessTypeEnumerator pTypes;
583 
584  LogPrint("DTDigitizer") << endl
585  << "------- SimHit: " << endl
586  << " Particle type = " << hit->particleType() << endl
587  << " process type = " << hit->processType() << endl
588  << " process type = " << hit->processType() << endl
589  // << " packedTrackId = " << hit->packedTrackId() << endl
590  << " trackId = " << hit->trackId() << endl // new,is the same as the
591  // previous?? FIXME-Check
592  << " |p| = " << hit->pabs() << endl
593  << " Energy loss = " << hit->energyLoss() << endl
594  // << " timeOffset = " << hit->timeOffset() << endl
595  // << " measurementPosition = " << hit->measurementPosition() << endl
596  // << " measurementDirection = " << hit->measurementDirection() << endl //FIXME
597  << " localDirection = " << hit->momentumAtEntry().unit() << endl //FIXME is it a versor?
598  << " Entry point = " << entryP << " cell x = " << xEntry << endl
599  << " Exit point = " << exitP << " cell x = " << xExit << endl
600  << " DR = = " << (exitP-entryP).mag() << endl
601  << " Dx = = " << (exitP-entryP).x() << endl
602  << " Cell w,h,l = (" << topo.cellWidth()
603  << " , " << topo.cellHeight()
604  << " , " << topo.cellLenght() << ") cm" << endl
605  << " DY entry from edge = " << topo.cellLenght()/2.-fabs(entryP.y())
606  << " DY exit from edge = " << topo.cellLenght()/2.-fabs(exitP.y())
607  << " entrySide = " << (int)entrySide
608  << " ; exitSide = " << (int)exitSide << endl;
609 
610 }
Side onWhichBorder(float x, float y, float z) const
Definition: DTTopology.cc:109
LocalVector momentumAtEntry() const
The momentum of the track that produced the hit, at entry point.
Definition: PSimHit.h:47
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
T y() const
Definition: PV3DBase.h:63
float cellWidth() const
Returns the cell width.
Definition: DTTopology.h:68
Local3DPoint exitPoint() const
Exit point in the local Det frame.
Definition: PSimHit.h:38
T z() const
Definition: PV3DBase.h:64
float cellHeight() const
Returns the cell height.
Definition: DTTopology.h:70
float pabs() const
fast and more accurate access to momentumAtEntry().mag()
Definition: PSimHit.h:63
Side
Sides of the cell.
Definition: DTTopology.h:88
Vector3DBase unit() const
Definition: Vector3DBase.h:57
unsigned short processType() const
Definition: PSimHit.h:118
float energyLoss() const
The energy deposit in the PSimHit, in ???.
Definition: PSimHit.h:75
int particleType() const
Definition: PSimHit.h:85
unsigned int trackId() const
Definition: PSimHit.h:102
float cellLenght() const
Definition: DTTopology.h:73
Local3DPoint entryPoint() const
Entry point in the local Det frame.
Definition: PSimHit.h:35
float DTDigitizer::externalDelays ( const DTLayer layer,
const DTWireId wireId,
const PSimHit hit 
) const
private

Definition at line 466 of file DTDigitizer.cc.

References DTTopology::cellLenght(), debug, PSimHit::localPosition(), DTLayer::specificTopology(), PSimHit::tof(), and PV3DBase< T, PVType, FrameType >::y().

468  {
469 
470  // Time of signal propagation along wire.
471 
472  float wireCoord = hit->localPosition().y();
473  float halfL = (layer->specificTopology().cellLenght())/2.;
474  float propgL = halfL - wireCoord; // the FE is always located at the pos coord.
475 
476  float propDelay = propgL/vPropWire;
477 
478  // Real TOF.
479  float tof = hit->tof();
480 
481  // Delays and t0 according to theSync
482 
483  double sync= theSync->digitizerOffset(&wireId,layer);
484 
485 
486  if (debug) {
487  LogPrint("DTDigitizer") << " propDelay =" << propDelay
488  << "; TOF=" << tof
489  << "; sync= " << sync
490  << endl;
491  }
492 
493  return propDelay + tof + sync;
494 }
float tof() const
deprecated name for timeOfFlight()
Definition: PSimHit.h:72
T y() const
Definition: PV3DBase.h:63
const DTTopology & specificTopology() const
Definition: DTLayer.cc:42
Local3DPoint localPosition() const
Definition: PSimHit.h:44
std::unique_ptr< DTDigiSyncBase > theSync
Definition: DTDigitizer.h:116
double vPropWire
Definition: DTDigitizer.h:108
float cellLenght() const
Definition: DTTopology.h:73
void DTDigitizer::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
)
overridevirtual

Implements edm::stream::EDProducerBase.

Definition at line 127 of file DTDigitizer.cc.

References debug, edm::EventID::event(), edm::EventSetup::get(), edm::Event::getByToken(), edm::RandomNumberGenerator::getEngine(), edm::EventBase::id(), DTWireId::layerId(), convertSQLitetoXML_cfg::output, edm::Handle< T >::product(), edm::Event::put(), edm::EventID::run(), trackerHits::simHits, edm::Event::streamID(), GeomDet::surface(), Surface::toGlobal(), and GloballyPositioned< T >::toLocal().

127  {
128 
130  CLHEP::HepRandomEngine* engine = &rng->getEngine(iEvent.streamID());
131 
132  if(debug)
133  LogPrint("DTDigitizer") << "--- Run: " << iEvent.id().run()
134  << " Event: " << iEvent.id().event() << endl;
135 
136  //************ 1 ***************
137  // create the container for the SimHits
138  // Handle<PSimHitContainer> simHits;
139  // iEvent.getByLabel("g4SimHits","MuonDTHits",simHits);
140 
141  // use MixCollection instead of the previous
143  iEvent.getByToken(cf_token, xFrame);
144 
145  auto_ptr<MixCollection<PSimHit> >
146  simHits( new MixCollection<PSimHit>(xFrame.product()) );
147 
148  // create the pointer to the Digi container
149  auto_ptr<DTDigiCollection> output(new DTDigiCollection());
150  // pointer to the DigiSimLink container
151  auto_ptr<DTDigiSimLinkCollection> outputLinks(new DTDigiSimLinkCollection());
152 
153  // Muon Geometry
154  ESHandle<DTGeometry> muonGeom;
155  iSetup.get<MuonGeometryRecord>().get(geometryType,muonGeom);
156 
157  // Magnetic Field
158  ESHandle<MagneticField> magnField;
159  iSetup.get<IdealMagneticFieldRecord>().get(magnField);
160 
161  //************ 2 ***************
162 
163  // These are sorted by DetId, i.e. by layer and then by wire #
164  // map<DTDetId, vector<const PSimHit*> > wireMap;
165  DTWireIdMap wireMap;
166 
167  for(MixCollection<PSimHit>::MixItr simHit = simHits->begin();
168  simHit != simHits->end(); simHit++){
169 
170  // Create the id of the wire, the simHits in the DT known also the wireId
171 
172  DTWireId wireId( (*simHit).detUnitId() );
173  // Fill the map
174  wireMap[wireId].push_back(&(*simHit));
175  }
176 
177  pair<float,bool> time(0.,false);
178 
179  //************ 3 ***************
180  // Loop over the wires
181  for(DTWireIdMapConstIter wire = wireMap.begin(); wire!=wireMap.end(); wire++){
182  // SimHit Container associated to the wire
183  const vector<const PSimHit*> & vhit = (*wire).second;
184  if(vhit.size()!=0) {
185  TDContainer tdCont; // It is a vector<pair<const PSimHit*,float> >;
186 
187  //************ 4 ***************
188  DTWireId wireId = (*wire).first;
189 
190  //const DTLayer* layer = dynamic_cast< const DTLayer* > (muonGeom->idToDet(wireId.layerId()));
191  const DTLayer* layer = muonGeom->layer(wireId.layerId());
192 
193  // Loop on the hits of this wire
194  for (vector<const PSimHit*>::const_iterator hit=vhit.begin();
195  hit != vhit.end(); hit++){
196  //************ 5 ***************
197  LocalPoint locPos = (*hit)->localPosition();
198 
199  const LocalVector BLoc=layer->surface().toLocal(magnField->inTesla(layer->surface().toGlobal(locPos)));
200 
201  time = computeTime(layer, wireId, *hit, BLoc, engine);
202 
203  //************ 6 ***************
204  if (time.second) {
205  tdCont.push_back(make_pair((*hit),time.first));
206  } else {
207  if (debug) LogPrint("DTDigitizer") << "hit discarded" << endl;
208  }
209  }
210 
211  //************ 7 ***************
212 
213  // the loading must be done by layer but
214  // the digitization must be done by wire (in order to take into account the dead time)
215 
216  storeDigis(wireId,tdCont,*output,*outputLinks);
217  }
218 
219  }
220 
221  //************ 8 ***************
222  // Load the Digi Container in the Event
223  //iEvent.put(output,"MuonDTDigis");
224  iEvent.put(output);
225  iEvent.put(outputLinks);
226 
227 }
RunNumber_t run() const
Definition: EventID.h:39
std::pair< float, bool > computeTime(const DTLayer *layer, const DTWireId &wireId, const PSimHit *hit, const LocalVector &BLoc, CLHEP::HepRandomEngine *)
Definition: DTDigitizer.cc:229
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:106
EventNumber_t event() const
Definition: EventID.h:41
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:462
virtual CLHEP::HepRandomEngine & getEngine(StreamID const &)=0
Use this engine in event methods.
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:40
void storeDigis(DTWireId &wireId, TDContainer &hits, DTDigiCollection &output, DTDigiSimLinkCollection &outputLinks)
Definition: DTDigitizer.cc:499
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:121
LocalPoint toLocal(const GlobalPoint &gp) const
edm::EDGetTokenT< CrossingFrame< PSimHit > > cf_token
Definition: DTDigitizer.h:132
DTWireIdMap::const_iterator DTWireIdMapConstIter
Definition: DTDigitizer.h:59
T const * product() const
Definition: Handle.h:81
std::string geometryType
Definition: DTDigitizer.h:118
const T & get() const
Definition: EventSetup.h:56
tuple simHits
Definition: trackerHits.py:16
std::vector< hitAndT > TDContainer
Definition: DTDigitizer.h:55
DTLayerId layerId() const
Return the corresponding LayerId.
Definition: DTWireId.h:62
edm::EventID id() const
Definition: EventBase.h:59
MuonDigiCollection< DTLayerId, DTDigi > DTDigiCollection
StreamID streamID() const
Definition: Event.h:80
MuonDigiCollection< DTLayerId, DTDigiSimLink > DTDigiSimLinkCollection
std::map< DTWireId, std::vector< const PSimHit * > > DTWireIdMap
Definition: DTDigitizer.h:57
void DTDigitizer::storeDigis ( DTWireId wireId,
TDContainer hits,
DTDigiCollection output,
DTDigiSimLinkCollection outputLinks 
)
private

Definition at line 499 of file DTDigitizer.cc.

References funct::abs(), DTDigi::countsTDC(), debug, MuonDigiCollection< IndexType, DigiType >::insertDigi(), DTWireId::layerId(), DTDigi::number(), funct::pow(), DTDigi::time(), and DTWireId::wire().

501  {
502 
503  //************ 7A ***************
504 
505  // sort signal times
506  sort(hits.begin(),hits.end(),hitLessT());
507 
508  //************ 7B ***************
509 
510  float wakeTime = -999999.0;
511  float resolTime = -999999.0;
512  int digiN = -1; // Digi identifier within the cell (for multiple digis)
513  DTDigi digi;
514 
515  // loop over signal times and drop signals inside dead time
516  for ( TDContainer::const_iterator hit = hits.begin() ; hit != hits.end() ;
517  hit++ ) {
518 
519  if (onlyMuHits && abs((*hit).first->particleType())!=13) continue;
520 
521  //************ 7C ***************
522 
523  float time = (*hit).second;
524  if ( time > wakeTime ) {
525  // Note that digi is constructed with a float value (in ns)
526  int wireN = wireId.wire();
527  digiN++;
528  digi = DTDigi(wireN, time, digiN);
529 
530  // Add association between THIS digi and the corresponding SimTrack
531  unsigned int SimTrackId = (*hit).first->trackId();
532  EncodedEventId evId = (*hit).first->eventId();
533  DTDigiSimLink digisimLink(wireN, digiN, time, SimTrackId, evId);
534 
535  if(debug) {
536  LogPrint("DTDigitizer")<<endl<<"---- DTDigitizer ----"<<endl;
537  LogPrint("DTDigitizer")<<"wireId: "<<wireId<<endl;
538  LogPrint("DTDigitizer")<<"sim. time = "<<time<<endl;
539  LogPrint("DTDigitizer")<<"digi number = "<< digi.number()<<", digi time = "<<digi.time()
540  <<", linked to SimTrack Id = "<<SimTrackId<<endl;
541  }
542 
543  //************ 7D ***************
544  if(digi.countsTDC() < pow(2.,16)){
545  DTLayerId layerID = wireId.layerId(); //taking the layer in which reside the wire
546  output.insertDigi(layerID, digi); // ordering Digis by layer
547  outputLinks.insertDigi(layerID, digisimLink);
548  wakeTime = time + deadTime;
549  resolTime = time + LinksTimeWindow;
550  }
551  else {
552  digiN--;
553  }
554  }
555  else if (MultipleLinks && time < resolTime){
556  int wireN = wireId.wire();
557  unsigned int SimTrackId = (*hit).first->trackId();
558  EncodedEventId evId = (*hit).first->eventId();
559  DTDigiSimLink digisimLink(wireN, digiN, time, SimTrackId, evId);
560  DTLayerId layerID = wireId.layerId();
561  outputLinks.insertDigi(layerID, digisimLink);
562 
563  if(debug) {
564  LogPrint("DTDigitizer")<<"\nAdded multiple link: \n"
565  <<"digi number = "<<digi.number()<<", digi time = "<<digi.time()<<" (sim. time = "<<time<<")"
566  <<", linked to SimTrack Id = "<<SimTrackId<<endl;
567  }
568  }
569  }
570 
571 }
bool onlyMuHits
Definition: DTDigitizer.h:113
float LinksTimeWindow
Definition: DTDigitizer.h:126
void insertDigi(const IndexType &index, const DigiType &digi)
insert a digi for a given DetUnit
double time() const
Get time in ns.
Definition: DTDigi.cc:63
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
Definition: DTDigi.h:17
uint32_t countsTDC() const
Get raw TDC count.
Definition: DTDigi.cc:65
int wire() const
Return the wire number.
Definition: DTWireId.h:56
DTLayerId layerId() const
Return the corresponding LayerId.
Definition: DTWireId.h:62
int number() const
Identifies different digis within the same.
Definition: DTDigi.cc:69
float deadTime
Definition: DTDigitizer.h:109
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
bool MultipleLinks
Definition: DTDigitizer.h:125

Friends And Related Function Documentation

friend class DTDigitizerAnalysis
friend

Definition at line 105 of file DTDigitizer.h.

Member Data Documentation

edm::EDGetTokenT<CrossingFrame<PSimHit> > DTDigitizer::cf_token
private

Definition at line 132 of file DTDigitizer.h.

std::string DTDigitizer::collection_for_XF
private

Definition at line 130 of file DTDigitizer.h.

float DTDigitizer::deadTime
private

Definition at line 109 of file DTDigitizer.h.

bool DTDigitizer::debug
private

Definition at line 111 of file DTDigitizer.h.

std::string DTDigitizer::geometryType
private

Definition at line 118 of file DTDigitizer.h.

bool DTDigitizer::IdealModel
private

Definition at line 121 of file DTDigitizer.h.

bool DTDigitizer::interpolate
private

Definition at line 112 of file DTDigitizer.h.

float DTDigitizer::LinksTimeWindow
private

Definition at line 126 of file DTDigitizer.h.

std::string DTDigitizer::mix_
private

Definition at line 129 of file DTDigitizer.h.

bool DTDigitizer::MultipleLinks
private

Definition at line 125 of file DTDigitizer.h.

bool DTDigitizer::onlyMuHits
private

Definition at line 113 of file DTDigitizer.h.

float DTDigitizer::smearing
private

Definition at line 110 of file DTDigitizer.h.

std::string DTDigitizer::syncName
private

Definition at line 115 of file DTDigitizer.h.

float DTDigitizer::theConstVDrift
private

Definition at line 122 of file DTDigitizer.h.

std::unique_ptr<DTDigiSyncBase> DTDigitizer::theSync
private

Definition at line 116 of file DTDigitizer.h.

double DTDigitizer::vPropWire
private

Definition at line 108 of file DTDigitizer.h.