#include <DTDigitizer.h>

Inheritance diagram for DTDigitizer:

Classes
struct	hitLessT

Public Member Functions
	DTDigitizer (const edm::ParameterSet &)

virtual void	produce (edm::Event &, const edm::EventSetup &)

	~DTDigitizer ()

Public Member Functions inherited from edm::EDProducer
	EDProducer ()

virtual	~EDProducer ()

Public Member Functions inherited from edm::ProducerBase
	ProducerBase ()

void	registerProducts (ProducerBase , ProductRegistry , ModuleDescription const &)

boost::function< void(const BranchDescription &)>	registrationCallback () const
	used by the fwk to register list of products More...

virtual	~ProducerBase ()

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< hitAndT >	TDContainer

Private Member Functions
float	asymGausSmear (double mean, double sigmaLeft, double sigmaRight) const

std::pair< float, bool >	computeTime (const DTLayer layer, const DTWireId &wireId, const PSimHit hit, const LocalVector &BLoc)

std::pair< float, bool >	driftTimeFromParametrization (float x, float alpha, float By, float Bz) 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
std::string	collection_for_XF

float	deadTime

bool	debug

std::string	geometryType

bool	IdealModel

bool	interpolate

float	LinksTimeWindow

bool	MultipleLinks

bool	onlyMuHits

float	smearing

std::string	syncName

float	theConstVDrift

CLHEP::RandFlat *	theFlatDistribution

CLHEP::RandGaussQ *	theGaussianDistribution

DTDigiSyncBase *	theSync

double	vPropWire

Friends
class	DTDigitizerAnalysis

Additional Inherited Members
Public Types inherited from edm::EDProducer
typedef EDProducer	ModuleType

typedef WorkerT< EDProducer >	WorkerType

Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList	TypeLabelList

Static Public Member Functions inherited from edm::EDProducer
static const std::string &	baseType ()

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &descriptions)

Protected Member Functions inherited from edm::EDProducer
CurrentProcessingContext const *	currentContext () const

Protected Member Functions inherited from edm::ProducerBase
template<class TProducer , class TMethod >
void	callWhenNewProductsRegistered (TProducer *iProd, TMethod iMethod)

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.

Date:: 2008/07/01 07:41:35

Revision:: 1.12

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

Definition at line 39 of file DTDigitizer.h.

Member Typedef Documentation

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

private

Definition at line 51 of file DTDigitizer.h.

typedef DTWireIdMap::const_iterator DTDigitizer::DTWireIdMapConstIter

private

Definition at line 53 of file DTDigitizer.h.

typedef DTWireIdMap::iterator DTDigitizer::DTWireIdMapIter

private

Definition at line 52 of file DTDigitizer.h.

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

private

Definition at line 48 of file DTDigitizer.h.

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

private

Definition at line 49 of file DTDigitizer.h.

Constructor & Destructor Documentation

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

explicit

Definition at line 63 of file DTDigitizer.cc.

References gather_cfg::cout, debug, edm::hlt::Exception, reco::get(), edm::RandomNumberGenerator::getEngine(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), and edm::Service< T >::isAvailable().

                                                   {
   
   // Set verbose output
   debug=conf_.getUntrackedParameter<bool>("debug"); 
     
   if (debug) cout<<"Creating a DTDigitizer"<<endl;
   
   //register the Producer with a label
   //produces<DTDigiCollection>("MuonDTDigis"); // FIXME: Do I pass it by ParameterSet?
   produces<DTDigiCollection>(); // FIXME: Do I pass it by ParameterSet?
   //  produces<DTDigiSimLinkCollection>("MuonDTDigiSimLinks");
   produces<DTDigiSimLinkCollection>();
   
   //Parameters:
 
   // build digis only for mu hits (for debug purposes) 
   onlyMuHits=conf_.getParameter<bool>("onlyMuHits");
   
   // interpolate parametrization function
   interpolate=conf_.getParameter<bool>("interpolate");
   
   // Velocity of signal propagation along the wire (cm/ns)
   // For the default value
   // cfr. CMS-IN 2000-021:   (2.56+-0.17)x1e8 m/s
   //      CMS NOTE 2003-17:  (0.244)  m/ns
   vPropWire=conf_.getParameter<double>("vPropWire"); //24.4
 
   // Dead time for signals on the same wire (number from M. Pegoraro)  
   deadTime=conf_.getParameter<double>("deadTime"); //150
   
   // further configurable smearing
   smearing=conf_.getParameter<double>("Smearing"); // 3.
 
   // Sync Algo
   syncName = conf_.getParameter<string>("SyncName");
   theSync = DTDigiSyncFactory::get()->create(syncName,conf_.getParameter<ParameterSet>("pset"));
 
   // Debug flag to switch to the Ideal model
   // it uses a constant drift velocity and doesn't set any external delay
   IdealModel = conf_.getParameter<bool>("IdealModel");
 
   // Constant drift velocity needed by the above flag
   if(IdealModel)
     theConstVDrift = conf_.getParameter<double>("IdealModelConstantDriftVelocity"); // 55 um/ns
   else theConstVDrift = 55.;
 
   // get random engine
   edm::Service<edm::RandomNumberGenerator> rng;
   if ( ! rng.isAvailable()) {
     throw cms::Exception("Configuration")
       << "RandomNumberGeneratorService for DTDigitizer missing in cfg file";
   }
   theGaussianDistribution = new CLHEP::RandGaussQ(rng->getEngine()); 
   theFlatDistribution = new CLHEP::RandFlat(rng->getEngine(), 0, 1); 
 
   // MultipleLinks=false ==> one-to-one correspondence between digis and SimHits 
   MultipleLinks = conf_.getParameter<bool>("MultipleLinks");
   // MultipleLinks=true ==> association of SimHits within a time window LinksTimeWindow 
   // (of the order of the resolution)
   LinksTimeWindow = conf_.getParameter<double>("LinksTimeWindow"); // (10 ns)
 
   //Name of Collection used for create the XF 
   collection_for_XF = conf_.getParameter<std::string>("InputCollection");
 
   //String to choice between ideal (the deafult) and (mis)aligned geometry for the digitization step 
   geometryType = conf_.getParameter<std::string>("GeometryType");
 }

DTDigitizer::~DTDigitizer ( )

Definition at line 132 of file DTDigitizer.cc.

                          {
   delete theGaussianDistribution;
   delete theFlatDistribution;
 }

Member Function Documentation

float DTDigitizer::asymGausSmear	(	double	mean,
		double	sigmaLeft,
		double	sigmaRight
	)		const

private

Definition at line 447 of file DTDigitizer.cc.

References f, and lumiQTWidget::t.

                                                                                        {
 
   double f = sigmaLeft/(sigmaLeft+sigmaRight);
   double t;
 
   if (theFlatDistribution->fire() <= f) {
     t = theGaussianDistribution->fire(mean,sigmaLeft);
     t = mean - fabs(t - mean);
   } else {
     t = theGaussianDistribution->fire(mean,sigmaRight);
     t = mean + fabs(t - mean);
   }
   return static_cast<float>(t);
 }

pair< float, bool > DTDigitizer::computeTime	(	const DTLayer *	layer,
		const DTWireId &	wireId,
		const PSimHit *	hit,
		const LocalVector &	BLoc
	)

private

Definition at line 236 of file DTDigitizer.cc.

References gather_cfg::cout, Vector3DBase< T, FrameTag >::cross(), 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(), DTLayer::specificTopology(), mathSSE::sqrt(), theta(), Vector3DBase< T, FrameTag >::unit(), DTWireId::wire(), DTTopology::wirePosition(), PV3DBase< T, PVType, FrameType >::x(), vdt::x, DTTopology::xMax, DTTopology::xMin, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

                                                                   { 
  
   LocalPoint entryP = hit->entryPoint();
   LocalPoint exitP = hit->exitPoint();
   int partType = hit->particleType();
 
   const DTTopology &topo = layer->specificTopology();
 
   // Pay attention: in CMSSW the rf of the SimHit is in the layer's rf
   
   if(debug)  cout<<"Hit local entry point: "<<entryP<<endl
          <<"Hit local exit point: "<<exitP<<endl;
 
   float xwire = topo.wirePosition(wireId.wire()); 
   float xEntry = entryP.x() - xwire;
   float xExit  = exitP.x() - xwire;
 
   if(debug) cout<<"wire position: "<<xwire
         <<" x entry in cell rf: "<<xEntry
         <<" x exit in cell rf: "<<xExit<<endl;
   
   DTTopology::Side entrySide = topo.onWhichBorder(xEntry,entryP.y(),entryP.z());
   DTTopology::Side exitSide  = topo.onWhichBorder(xExit,exitP.y(),exitP.z());
 
    if (debug) dumpHit(hit, xEntry, xExit,topo);
 
   // The bolean is used to flag the drift time computation
   pair<float,bool> driftTime(0.,false);  
 
   // if delta in gas->ignore, since it is included in the parametrisation.
   // FIXME: should check that it is actually a delta ray produced by a nearby
   // muon hit. 
 
   if (partType == 11 && entrySide == DTTopology::none) {
      if (debug) cout << "    e- hit in gas; discarding " << endl;
     return driftTime;
   }
 
   float By = BLoc.y();
   float Bz = BLoc.z();
 
   // Radius and sagitta according to direction of momentum
   // (just for printing)
   // NOTE: in cmsim, d is always taken // pHat!
   LocalVector d = (exitP-entryP);
   LocalVector pHat = hit->localDirection().unit();
   LocalVector hHat = (d.cross(pHat.cross(d))).unit();
   float cosAlpha = hHat.dot(pHat);
   float sinAlpha = sqrt(1.-cosAlpha*cosAlpha);
   float radius_P = (d.mag())/(2.*cosAlpha);
   float sagitta_P = radius_P*(1.-sinAlpha);
 
   // Radius, sagitta according to field bending
   // (just for printing)
   float halfd = d.mag()/2.;
   float BMag = BLoc.mag();
   LocalVector pT = (pHat - (BLoc.unit()*pHat.dot(BLoc.unit())))*(hit->pabs());
   float radius_B = (pT.mag()/(0.3*BMag))*100.;
   float sagitta_B;
   if (radius_B > halfd) {
     sagitta_B = radius_B - sqrt(radius_B*radius_B - halfd*halfd);
   } else {
     sagitta_B = radius_B;
   }
 
   // cos(delta), delta= angle between direction at entry and hit segment
   // (just for printing)
   float delta = pHat.dot(d.unit());
   if (debug) cout << "   delta                 = " << delta  << endl
           << "   cosAlpha              = " << cosAlpha << endl
           << "   sinAlpha              = " << sinAlpha << endl
           << "   pMag                  = " << pT.mag() << endl
           << "   bMag                  = " << BMag << endl
           << "   pT                    = " << pT << endl
           << "   halfd                 = " << halfd << endl
           << "   radius_P  (cm)        = " << radius_P << endl
           << "   sagitta_P (um)        = " << sagitta_P*10000. << endl
           << "   radius_B  (cm)        = " << radius_B << endl
           << "   sagitta_B (um)        = " << sagitta_B*10000. << endl;
 
   // Select cases where parametrization can not be used.
   bool noParametrisation = 
     ( ( entrySide == DTTopology::none || exitSide == DTTopology::none ) // case # 2,3,8,9 or 11
       || (entrySide == exitSide)                   // case # 4 or 10
       || ((entrySide == DTTopology::xMin && exitSide == DTTopology::xMax) || 
       (entrySide == DTTopology::xMax && exitSide == DTTopology::xMin)) // Hit is case # 7
       );
 
   // FIXME: now, debug warning only; consider treating those 
   // with TM algo. 
   if ( delta < 0.99996 // Track is not straight. FIXME: use sagitta?
        && (noParametrisation == false)) {
     if (debug) cout << "*** WARNING: hit is not straight, type = " << partType << endl;
   }
 
   //************ 5A ***************
 
   if (!noParametrisation) {
    
     LocalVector dir = hit->momentumAtEntry(); // ex Measurement3DVector dir = hit->measurementDirection(); //FIXME?
     float theta = atan(dir.x()/-dir.z())*180/M_PI;
 
     // FIXME: use dir if M.S. is included as GARFIELD option...
     //        otherwise use hit segment dirction.
     //    LocalVector dir0 = (exitP-entryP).unit();
     //    float theta = atan(dir0.x()/-dir0.z())*180/M_PI;
     float x;
 
     Local3DPoint pt = hit->localPosition();  //ex Measurement3DPoint pt = hit->measurementPosition(); // FIXME?
      
     if(fabs(pt.z()) < 0.002) { 
       // hit center within 20 um from z=0, no need to extrapolate.
       x = pt.x() - xwire;
     } else {
       x = xEntry - (entryP.z()*(xExit-xEntry))/(exitP.z()-entryP.z());
     }
 
     if(IdealModel) return make_pair(fabs(x)/theConstVDrift,true);
     else driftTime = driftTimeFromParametrization(x, theta, By, Bz);
 
   }
 
  
   if ((driftTime.second)==false) {
     // Parametrisation not applicable, or failed. Use time map.
     driftTime = driftTimeFromTimeMap();
   }
   
   //************ 5B ***************
 
   // Signal propagation, TOF etc.
   if (driftTime.second) {
     driftTime.first += externalDelays(layer,wireId,hit);
   }
   return driftTime;
 } 

pair< float, bool > DTDigitizer::driftTimeFromParametrization	(	float	x,
		float	alpha,
		float	By,
		float	Bz
	)		const

private

Definition at line 376 of file DTDigitizer.cc.

References gather_cfg::cout, 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, theta(), and cond::rpcobgas::time.

                                                                                                          {
 
   // Convert from CMSSW frame/units r.f. to parametrization ones.
   x *= 10.;  //cm -> mm 
 
   // FIXME: Current parametrisation can extrapolate above 21 mm,
   // however a detailed study is needed before using this.
    if (fabs(x) > 21.) {
      if (debug) cout << "*** WARNING: parametrisation: x out of range = "
              << x << ", skipping" << endl;
      return pair<float,bool>(0.f,false);
   }
 
   // Different r.f. of the parametrization:
   // X_par = X_ORCA; Y_par=Z_ORCA; Z_par = -Y_ORCA  
 
   float By_par = Bz;  // Bnorm
   float Bz_par = -By; // Bwire
   float theta_par = theta;
 
   // Parametrisation uses interpolation up to |theta|=45 deg,
   // |Bwire|=0.4, |Bnorm|=0.75; extrapolation above.
   if (fabs(theta_par)>45.) {
     if (debug) cout << "*** WARNING: extrapolating theta > 45: "
             << theta << endl;
     // theta_par = min(fabs(theta_par),45.f)*((theta_par<0.)?-1.:1.);
   }
   if (fabs(By_par)>0.75) {
     if (debug) cout << "*** WARNING: extrapolating Bnorm > 0.75: "
             << By_par << endl;
     // By_par = min(fabs(By_par),0.75f)*((By_par<0.)?-1.:1.);
   }
   if (fabs(Bz_par)>0.4) {
     if (debug) cout << "*** WARNING: extrapolating Bwire >0.4: "
             << Bz_par << endl;
     // Bz_par = min(fabs(Bz_par),0.4)*((Bz_par<0.)?-1.:1.);
   }
 
   DTDriftTimeParametrization::drift_time DT;
   static DTDriftTimeParametrization par;
   unsigned short flag = par.MB_DT_drift_time (x, theta_par, By_par, Bz_par, 0, &DT, interpolate);
 
   if (debug) {
     cout << "    Parametrisation: x, theta, Bnorm, Bwire = "
      << x << " " <<  theta_par << " " << By_par << " " << Bz_par << endl
      << "  time=" << DT.t_drift
      << "  sigma_m=" <<  DT.t_width_m
      << "  sigma_p=" <<  DT.t_width_p << endl;
     if (flag!=1) {
       cout << "*** WARNING: call to parametrisation failed" << endl;
       return pair<float,bool>(0.f,false); 
     }
   }
 
   // Double half-gaussian smearing
   float time = asymGausSmear(DT.t_drift, DT.t_width_m, DT.t_width_p);
 
   // Do not allow the smearing to lead to negative values
   time = max(time,0.f);
 
   // Apply a Gaussian smearing to account for electronic effects (cf. 2004 TB analysis)
   // The width of the Gaussian can be configured with the "Smearing" parameter
 
   double u = theGaussianDistribution->fire(0.,smearing);
   time += u;
 
   if (debug) cout << "  drift time = " << time << endl;
   
   return pair<float,bool>(time,true); 
 }

pair< float, bool > DTDigitizer::driftTimeFromTimeMap ( ) const

private

Definition at line 462 of file DTDigitizer.cc.

References gather_cfg::cout, and debug.

                                                          {
   // FIXME: not yet implemented.
   if (debug) cout << "   TimeMap " << endl;
   return pair<float,bool>(0.,false);
 }

void DTDigitizer::dumpHit	(	const PSimHit *	hit,
		float	xEntry,
		float	xExit,
		const DTTopology &	topo
	)

private

Definition at line 577 of file DTDigitizer.cc.

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

                                                   {
   
   LocalPoint entryP = hit->entryPoint();
   LocalPoint exitP = hit->exitPoint();
   
   DTTopology::Side entrySide = topo.onWhichBorder(xEntry,entryP.y(),entryP.z());
   DTTopology::Side exitSide  = topo.onWhichBorder(xExit,exitP.y(),exitP.z());
   //  ProcessTypeEnumerator pTypes;
   
   cout << endl
        << "------- SimHit: " << endl
        << "   Particle type         = " << hit->particleType() << endl
        << "   process type          = " << hit->processType() << endl
        << "   process type          = " << hit->processType() << endl
     // << "   packedTrackId         = " << hit->packedTrackId() << endl
        << "   trackId               = " << hit->trackId() << endl // new,is the same as the
                                                                   // previous?? FIXME-Check
        << "   |p|                   = " << hit->pabs() << endl
        << "   Energy loss           = " << hit->energyLoss() << endl
     // << "   timeOffset            = " << hit->timeOffset() << endl
     // << "   measurementPosition   = " << hit->measurementPosition() << endl
     // << "   measurementDirection  = " << hit->measurementDirection() << endl      //FIXME
        << "   localDirection        = " << hit->momentumAtEntry().unit() << endl    //FIXME is it a versor?
        << "   Entry point           = " << entryP << " cell x = " << xEntry << endl
        << "   Exit point            = " << exitP << " cell x = " << xExit << endl
        << "   DR =                  = " << (exitP-entryP).mag() << endl
        << "   Dx =                  = " << (exitP-entryP).x() << endl
        << "   Cell w,h,l            = (" << topo.cellWidth()
        << " , " << topo.cellHeight() 
        << " , " << topo.cellLenght() << ") cm" << endl
        << "   DY entry from edge    = " << topo.cellLenght()/2.-fabs(entryP.y())
        << "   DY exit  from edge    = " << topo.cellLenght()/2.-fabs(exitP.y())
        << "   entrySide = "  << (int)entrySide
        << " ; exitSide = " << (int)exitSide << endl;
 
 }

float DTDigitizer::externalDelays	(	const DTLayer *	layer,
		const DTWireId &	wireId,
		const PSimHit *	hit
	)		const

private

Definition at line 470 of file DTDigitizer.cc.

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

                                             {
   
   // Time of signal propagation along wire.
   
   float wireCoord = hit->localPosition().y();
   float halfL     = (layer->specificTopology().cellLenght())/2.;
   float propgL = halfL - wireCoord; // the FE is always located at the pos coord.
 
   float propDelay = propgL/vPropWire;
 
   // Real TOF.
   float tof = hit->tof();  
 
   // Delays and t0 according to theSync
 
   double sync= theSync->digitizerOffset(&wireId,layer);
 
 
   if (debug) {
     cout << "    propDelay =" << propDelay
      << "; TOF=" << tof
      << "; sync= " << sync
      << endl;
   }
   
   return propDelay + tof + sync;
 }

void DTDigitizer::produce	(	edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)

virtual

Implements edm::EDProducer.

Definition at line 138 of file DTDigitizer.cc.

References MixCollection< T >::begin(), gather_cfg::cout, debug, edm::EventID::event(), edm::EventSetup::get(), edm::Event::getByLabel(), edm::EventBase::id(), DTWireId::layerId(), convertSQLitetoXML_cfg::output, edm::Handle< T >::product(), edm::Event::put(), edm::EventID::run(), trackerHits::simHits, GeomDet::surface(), cond::rpcobgas::time, Surface::toGlobal(), and GloballyPositioned< T >::toLocal().

                                                                 {
   if(debug)
     cout << "--- Run: " << iEvent.id().run()
      << " Event: " << iEvent.id().event() << endl;
   
   //************ 1 ***************
    // create the container for the SimHits
   //  Handle<PSimHitContainer> simHits; 
   //  iEvent.getByLabel("g4SimHits","MuonDTHits",simHits);
     
   // use MixCollection instead of the previous
   Handle<CrossingFrame<PSimHit> > xFrame;
   iEvent.getByLabel("mix",collection_for_XF,xFrame);
   
   auto_ptr<MixCollection<PSimHit> > 
     simHits( new MixCollection<PSimHit>(xFrame.product()) );
 
    // create the pointer to the Digi container
   auto_ptr<DTDigiCollection> output(new DTDigiCollection());
    // pointer to the DigiSimLink container
   auto_ptr<DTDigiSimLinkCollection> outputLinks(new DTDigiSimLinkCollection());
   
   // Muon Geometry
   ESHandle<DTGeometry> muonGeom;
   iSetup.get<MuonGeometryRecord>().get(geometryType,muonGeom);
 
   // Magnetic Field  
   ESHandle<MagneticField> magnField;
   iSetup.get<IdealMagneticFieldRecord>().get(magnField);
 
   //************ 2 ***************
 
   // These are sorted by DetId, i.e. by layer and then by wire #
   //  map<DTDetId, vector<const PSimHit*> > wireMap;     
   DTWireIdMap wireMap;     
   
   for(MixCollection<PSimHit>::MixItr simHit = simHits->begin();
        simHit != simHits->end(); simHit++){
     
     // Create the id of the wire, the simHits in the DT known also the wireId
  
     DTWireId wireId( (*simHit).detUnitId() );
     // Fill the map
     wireMap[wireId].push_back(&(*simHit));
   }
   
   pair<float,bool> time(0.,false);
 
   //************ 3 ***************
   // Loop over the wires
   for(DTWireIdMapConstIter wire = wireMap.begin(); wire!=wireMap.end(); wire++){
     // SimHit Container associated to the wire
     const vector<const PSimHit*> & vhit = (*wire).second; 
     if(vhit.size()!=0) {
       TDContainer tdCont; // It is a vector<pair<const PSimHit*,float> >;
       
       //************ 4 ***************
       DTWireId wireId = (*wire).first;
 
       //const DTLayer* layer = dynamic_cast< const DTLayer* > (muonGeom->idToDet(wireId.layerId()));
       const DTLayer* layer = muonGeom->layer(wireId.layerId()); 
 
       // Loop on the hits of this wire    
       for (vector<const PSimHit*>::const_iterator hit=vhit.begin();
        hit != vhit.end(); hit++){
     //************ 5 ***************
     LocalPoint locPos = (*hit)->localPosition();
 
     const LocalVector BLoc=layer->surface().toLocal(magnField->inTesla(layer->surface().toGlobal(locPos)));
 
     time = computeTime(layer, wireId, *hit, BLoc); 
 
     //************ 6 ***************
     if (time.second) {
       tdCont.push_back(make_pair((*hit),time.first));
     } else {
       if (debug) cout << "hit discarded" << endl;
     }
       }
 
       //************ 7 ***************
 
       // the loading must be done by layer but
       // the digitization must be done by wire (in order to take into account the dead time)
 
       storeDigis(wireId,tdCont,*output,*outputLinks);
     }
     
   }
 
   //************ 8 ***************  
   // Load the Digi Container in the Event
   //iEvent.put(output,"MuonDTDigis");
   iEvent.put(output);
   iEvent.put(outputLinks);
 
 }

void DTDigitizer::storeDigis	(	DTWireId &	wireId,
		TDContainer &	hits,
		DTDigiCollection &	output,
		DTDigiSimLinkCollection &	outputLinks
	)

private

Definition at line 503 of file DTDigitizer.cc.

References abs, DTDigi::countsTDC(), gather_cfg::cout, debug, MuonDigiCollection< IndexType, DigiType >::insertDigi(), DTWireId::layerId(), DTDigi::number(), funct::pow(), python.multivaluedict::sort(), DTDigi::time(), cond::rpcobgas::time, and DTWireId::wire().

                                                                                 {
 
   //************ 7A ***************
 
   // sort signal times
   sort(hits.begin(),hits.end(),hitLessT());
 
   //************ 7B ***************
 
   float wakeTime = -999999.0;
   float resolTime = -999999.0;
   int digiN = -1; // Digi identifier within the cell (for multiple digis)
   DTDigi digi;
 
   // loop over signal times and drop signals inside dead time
   for ( TDContainer::const_iterator hit = hits.begin() ; hit != hits.end() ; 
     hit++ ) {
 
     if (onlyMuHits && abs((*hit).first->particleType())!=13) continue;
 
     //************ 7C ***************
     
     float time = (*hit).second;
     if ( time > wakeTime ) {
       // Note that digi is constructed with a float value (in ns)
       int wireN = wireId.wire();
       digiN++;
       digi = DTDigi(wireN, time, digiN);
 
       // Add association between THIS digi and the corresponding SimTrack
       unsigned int SimTrackId = (*hit).first->trackId();
       EncodedEventId evId = (*hit).first->eventId();
       DTDigiSimLink digisimLink(wireN, digiN, time, SimTrackId, evId);
 
       if(debug) {
     cout<<endl<<"---- DTDigitizer ----"<<endl;
     cout<<"wireId: "<<wireId<<endl;
     cout<<"sim. time = "<<time<<endl;
     cout<<"digi number = "<< digi.number()<<", digi time = "<<digi.time()
         <<", linked to SimTrack Id = "<<SimTrackId<<endl;
       }
 
       //************ 7D ***************
       if(digi.countsTDC() < pow(2.,16)){
     DTLayerId layerID = wireId.layerId();  //taking the layer in which reside the wire
         output.insertDigi(layerID, digi); // ordering Digis by layer
         outputLinks.insertDigi(layerID, digisimLink);
     wakeTime = time + deadTime;
     resolTime = time + LinksTimeWindow;
       } 
       else {
     digiN--;
       }
     }
     else if (MultipleLinks && time < resolTime){
       int wireN = wireId.wire();
       unsigned int SimTrackId = (*hit).first->trackId();
       EncodedEventId evId = (*hit).first->eventId();
       DTDigiSimLink digisimLink(wireN, digiN, time, SimTrackId, evId);
       DTLayerId layerID = wireId.layerId(); 
       outputLinks.insertDigi(layerID, digisimLink);
 
       if(debug) { 
     cout<<"\nAdded multiple link: \n"
         <<"digi number = "<<digi.number()<<", digi time = "<<digi.time()<<" (sim. time = "<<time<<")"
         <<", linked to SimTrack Id = "<<SimTrackId<<endl;
       }
     }
   }
   
 }

Friends And Related Function Documentation

friend class DTDigitizerAnalysis

friend

Definition at line 98 of file DTDigitizer.h.

Member Data Documentation

std::string DTDigitizer::collection_for_XF

private

Definition at line 126 of file DTDigitizer.h.

float DTDigitizer::deadTime

private

Definition at line 102 of file DTDigitizer.h.

bool DTDigitizer::debug

private

Definition at line 104 of file DTDigitizer.h.

std::string DTDigitizer::geometryType

private

Definition at line 111 of file DTDigitizer.h.

bool DTDigitizer::IdealModel

private

Definition at line 114 of file DTDigitizer.h.

bool DTDigitizer::interpolate

private

Definition at line 105 of file DTDigitizer.h.

float DTDigitizer::LinksTimeWindow

private

Definition at line 123 of file DTDigitizer.h.

bool DTDigitizer::MultipleLinks

private

Definition at line 122 of file DTDigitizer.h.

bool DTDigitizer::onlyMuHits

private

Definition at line 106 of file DTDigitizer.h.

float DTDigitizer::smearing

private

Definition at line 103 of file DTDigitizer.h.

std::string DTDigitizer::syncName

private

Definition at line 108 of file DTDigitizer.h.

float DTDigitizer::theConstVDrift

private

Definition at line 115 of file DTDigitizer.h.

CLHEP::RandFlat* DTDigitizer::theFlatDistribution

private

Definition at line 119 of file DTDigitizer.h.

CLHEP::RandGaussQ* DTDigitizer::theGaussianDistribution

private

Definition at line 118 of file DTDigitizer.h.

DTDigiSyncBase* DTDigitizer::theSync

private

Definition at line 109 of file DTDigitizer.h.

double DTDigitizer::vPropWire

private

Definition at line 101 of file DTDigitizer.h.

Classes

Public Member Functions

Private Types

Private Member Functions

Private Attributes

Friends

Additional Inherited Members

Detailed Description

Member Typedef Documentation

Constructor & Destructor Documentation

Member Function Documentation

Friends And Related Function Documentation

Member Data Documentation