HitEff Class Reference

#include <HitEff.h>

Inheritance diagram for HitEff:

Public Member Functions
	HitEff (const edm::ParameterSet &conf)
	~HitEff () override=default
Public Member Functions inherited from edm::one::EDAnalyzer< edm::one::SharedResources >
	EDAnalyzer ()=default
	EDAnalyzer (const EDAnalyzer &)=delete
SerialTaskQueue *	globalLuminosityBlocksQueue () final
SerialTaskQueue *	globalRunsQueue () final
const EDAnalyzer &	operator= (const EDAnalyzer &)=delete
bool	wantsGlobalLuminosityBlocks () const final
bool	wantsGlobalRuns () const final
bool	wantsInputProcessBlocks () const final
bool	wantsProcessBlocks () const final
Public Member Functions inherited from edm::one::EDAnalyzerBase
void	callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
	EDAnalyzerBase ()
ModuleDescription const &	moduleDescription () const
bool	wantsStreamLuminosityBlocks () const
bool	wantsStreamRuns () const
	~EDAnalyzerBase () override
Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo >	consumesInfo () const
void	convertCurrentProcessAlias (std::string const &processName)
	Convert "@currentProcess" in InputTag process names to the actual current process name. More...
	EDConsumerBase ()
	EDConsumerBase (EDConsumerBase const &)=delete
	EDConsumerBase (EDConsumerBase &&)=default
ESResolverIndex const *	esGetTokenIndices (edm::Transition iTrans) const
std::vector< ESResolverIndex > const &	esGetTokenIndicesVector (edm::Transition iTrans) const
std::vector< ESRecordIndex > const &	esGetTokenRecordIndicesVector (edm::Transition iTrans) const
ProductResolverIndexAndSkipBit	indexFrom (EDGetToken, BranchType, TypeID const &) const
void	itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
void	itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
std::vector< ProductResolverIndexAndSkipBit > const &	itemsToGetFrom (BranchType iType) const
void	labelsForToken (EDGetToken iToken, Labels &oLabels) const
void	modulesWhoseProductsAreConsumed (std::array< std::vector< ModuleDescription const *> *, NumBranchTypes > &modulesAll, std::vector< ModuleProcessName > &modulesInPreviousProcesses, ProductRegistry const &preg, std::map< std::string, ModuleDescription const *> const &labelsToDesc, std::string const &processName) const
EDConsumerBase const &	operator= (EDConsumerBase const &)=delete
EDConsumerBase &	operator= (EDConsumerBase &&)=default
bool	registeredToConsume (ProductResolverIndex, bool, BranchType) const
void	selectInputProcessBlocks (ProductRegistry const &productRegistry, ProcessBlockHelperBase const &processBlockHelperBase)
ProductResolverIndexAndSkipBit	uncheckedIndexFrom (EDGetToken) const
void	updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
void	updateLookup (eventsetup::ESRecordsToProductResolverIndices const &)
virtual	~EDConsumerBase () noexcept(false)

Private Member Functions
void	analyze (const edm::Event &e, const edm::EventSetup &c) override
void	beginJob () override
void	endJob () override

Private Attributes
bool	addCommonMode_
bool	addLumi_
unsigned int	bunchx
const edm::ESGetToken< Chi2MeasurementEstimatorBase, TrackingComponentsRecord >	chi2MeasurementEstimatorToken_
float	ClusterLocErrX
float	ClusterLocErrY
float	ClusterLocX
float	ClusterLocY
const edm::EDGetTokenT< edmNew::DetSetVector< SiStripCluster > >	clusters_token_
float	ClusterStoN
const edm::EDGetTokenT< reco::TrackCollection >	combinatorialTracks_token_
float	commonMode
const edm::EDGetTokenT< edm::DetSetVector< SiStripRawDigi > >	commonModeToken_
int	compSettings
edm::ParameterSet	conf_
const edm::ESGetToken< StripClusterParameterEstimator, TkStripCPERecord >	cpeToken_
bool	cutOnTracks_
bool	DEBUG
const edm::EDGetTokenT< DetIdCollection >	digisCol_token_
const edm::EDGetTokenT< DetIdVector >	digisVec_token_
bool	doMissingHitsRecovery_
unsigned int	event
int	events
int	EventTrackCKF
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord >	geomToken_
bool	highPurity
std::vector< unsigned int >	hitRecoveryCounters
std::vector< unsigned int >	hitTotalCounters
unsigned int	Id
float	instLumi
unsigned int	layers
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord >	magFieldToken_
const edm::ESGetToken< MeasurementTracker, CkfComponentsRecord >	measurementTkToken_
const edm::EDGetTokenT< OnlineLuminosityRecord >	metaDataToken_
std::vector< int >	missHitPerLayer
unsigned int	ModIsBad
int	nHits
const edm::ESGetToken< Propagator, TrackingComponentsRecord >	propagatorToken_
float	pT
float	PU
float	ResX
float	ResXSig
unsigned int	run
const edm::EDGetTokenT< LumiScalersCollection >	scalerToken_
SiStripClusterInfo	siStripClusterInfo_
unsigned int	SiStripQualBad
const edm::ESGetToken< SiStripQuality, SiStripQualityRcd >	siStripQualityToken_
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd >	topoToken_
int	totalNbHits
int	tquality
const edm::EDGetTokenT< MeasurementTrackerEvent >	trackerEvent_token_
unsigned int	trackMultiplicityCut_
TTree *	traj
const edm::EDGetTokenT< std::vector< Trajectory > >	trajectories_token_
float	TrajGlbX
float	TrajGlbY
float	TrajGlbZ
unsigned int	trajHitValid
float	TrajLocAngleX
float	TrajLocAngleY
float	TrajLocErrX
float	TrajLocErrY
float	TrajLocX
float	TrajLocY
const edm::EDGetTokenT< TrajTrackAssociationCollection >	trajTrackAsso_token_
bool	useAllHitsFromTracksWithMissingHits_
bool	useFirstMeas_
bool	useLastMeas_
unsigned int	whatlayer
bool	withinAcceptance

Additional Inherited Members
Public Types inherited from edm::one::EDAnalyzerBase
typedef EDAnalyzerBase	ModuleType
Public Types inherited from edm::EDConsumerBase
typedef ProductLabels	Labels
Static Public Member Functions inherited from edm::one::EDAnalyzerBase
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)
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< B >	consumes (edm::InputTag tag) noexcept
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 ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto	esConsumes ()
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto	esConsumes (ESInputTag const &tag)
template<Transition Tr = Transition::Event>
constexpr auto	esConsumes ()
template<Transition Tr = Transition::Event>
auto	esConsumes (ESInputTag tag)
template<Transition Tr = Transition::Event>
ESGetTokenGeneric	esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
	Used with EventSetupRecord::doGet. More...
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)
void	resetItemsToGetFrom (BranchType iType)

Detailed Description

Definition at line 53 of file HitEff.h.

Constructor & Destructor Documentation

HitEff()

HitEff::HitEff ( const edm::ParameterSet &  conf )

explicit

Definition at line 74 of file HitEff.cc.

References addCommonMode_, addLumi_, compSettings, conf_, cutOnTracks_, DEBUG, doMissingHitsRecovery_, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), hitRecoveryCounters, hitTotalCounters, TFileService::kSharedResource, layers, trackMultiplicityCut_, useAllHitsFromTracksWithMissingHits_, useFirstMeas_, and useLastMeas_.

    : scalerToken_(consumes<LumiScalersCollection>(conf.getParameter<edm::InputTag>("lumiScalers"))),
      metaDataToken_(consumes<OnlineLuminosityRecord>(conf.getParameter<edm::InputTag>("metadata"))),
      commonModeToken_(mayConsume<edm::DetSetVector<SiStripRawDigi> >(conf.getParameter<edm::InputTag>("commonMode"))),
      siStripClusterInfo_(consumesCollector()),
      combinatorialTracks_token_(
          consumes<reco::TrackCollection>(conf.getParameter<edm::InputTag>("combinatorialTracks"))),
      trajectories_token_(consumes<std::vector<Trajectory> >(conf.getParameter<edm::InputTag>("trajectories"))),
      trajTrackAsso_token_(consumes<TrajTrackAssociationCollection>(conf.getParameter<edm::InputTag>("trajectories"))),
      clusters_token_(
          consumes<edmNew::DetSetVector<SiStripCluster> >(conf.getParameter<edm::InputTag>("siStripClusters"))),
      digisCol_token_(consumes(conf.getParameter<edm::InputTag>("siStripDigis"))),
      digisVec_token_(consumes(conf.getParameter<edm::InputTag>("siStripDigis"))),
      trackerEvent_token_(consumes<MeasurementTrackerEvent>(conf.getParameter<edm::InputTag>("trackerEvent"))),
      topoToken_(esConsumes()),
      geomToken_(esConsumes()),
      cpeToken_(esConsumes(edm::ESInputTag("", "StripCPEfromTrackAngle"))),
      siStripQualityToken_(esConsumes()),
      magFieldToken_(esConsumes()),
      measurementTkToken_(esConsumes()),
      chi2MeasurementEstimatorToken_(esConsumes(edm::ESInputTag("", "Chi2"))),
      propagatorToken_(esConsumes(edm::ESInputTag("", "PropagatorWithMaterial"))),
      conf_(conf) {
  usesResource(TFileService::kSharedResource);
  compSettings = conf_.getUntrackedParameter<int>("CompressionSettings", -1);
  layers = conf_.getParameter<int>("Layer");
  DEBUG = conf_.getParameter<bool>("Debug");
  addLumi_ = conf_.getUntrackedParameter<bool>("addLumi", false);
  addCommonMode_ = conf_.getUntrackedParameter<bool>("addCommonMode", false);
  cutOnTracks_ = conf_.getUntrackedParameter<bool>("cutOnTracks", false);
  trackMultiplicityCut_ = conf.getUntrackedParameter<unsigned int>("trackMultiplicity", 100);
  useFirstMeas_ = conf_.getUntrackedParameter<bool>("useFirstMeas", false);
  useLastMeas_ = conf_.getUntrackedParameter<bool>("useLastMeas", false);
  useAllHitsFromTracksWithMissingHits_ =
      conf_.getUntrackedParameter<bool>("useAllHitsFromTracksWithMissingHits", false);
  doMissingHitsRecovery_ = conf_.getUntrackedParameter<bool>("doMissingHitsRecovery", false);

  hitRecoveryCounters.resize(k_END_OF_LAYERS, 0);
  hitTotalCounters.resize(k_END_OF_LAYERS, 0);
}

~HitEff()

HitEff::~HitEff ( )

overridedefault

Member Function Documentation

analyze()

void HitEff::analyze ( const edm::Event &  e, const edm::EventSetup &  c  )

overrideprivatevirtual

Implements edm::one::EDAnalyzerBase.

Definition at line 177 of file HitEff.cc.

References funct::abs(), addCommonMode_, addLumi_, AnalyticalPropagator_cfi::AnalyticalPropagator, anyDirection, OnlineLuminosityRecord::avgPileUp(), edm::AssociationMap< Tag >::begin(), Surface::bounds(), bunchx, nano_mu_local_reco_cff::chi2, chi2MeasurementEstimatorToken_, ClusterLocErrX, ClusterLocErrY, ClusterLocX, ClusterLocY, clusters_token_, ClusterStoN, combinatorialTracks_token_, commonMode, commonModeToken_, GeometricSearchDet::compatibleDets(), edm::AssociationMap< edm::OneToOne< std::vector< Trajectory >, reco::TrackCollection, unsigned short > >::const_iterator, cpeToken_, cutOnTracks_, digisCol_token_, digisVec_token_, doMissingHitsRecovery_, MillePedeFileConverter_cfg::e, edm::AssociationMap< Tag >::end(), events, EventTrackCKF, Exception, geomToken_, SiStripQuality::getBadApvs(), edm::EventSetup::getData(), edm::EventSetup::getHandle(), highPurity, muons_cff::highPurity, hitRecoveryCounters, hitTotalCounters, mps_fire::i, Id, TrackerGeometry::idToDetUnit(), cuy::ii, SiStripClusterInfo::initEvent(), input, OnlineLuminosityRecord::instLumi(), instLumi, createfilelist::int, sistrip::SpyUtilities::isValid(), edm::HandleBase::isValid(), ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, nano_mu_digi_cff::layer, layers, StripClusterParameterEstimator::localParameters(), LogDebug, V0Monitor_cfi::lumiScalers, magFieldToken_, measurementTkToken_, HLT_2024v14_cff::measurementTrackerEvent, metaDataToken_, missHitPerLayer, combine::missing, ModIsBad, heavyionUCCDQM_cfi::nClusters, reco::MuonTime::nDof, nHits, StripTopology::nstrips(), AlCaHLTBitMon_ParallelJobs::p, TrapezoidalPlaneBounds::parameters(), edm::Handle< T >::product(), TrackCandidateProducer_cfi::propagator, propagatorToken_, pT, PU, nano_mu_digi_cff::rawId, rpcPointValidation_cfi::recHit, TrajectoryMeasurement::recHit(), ResX, ResXSig, run, scalerToken_, SiStripClusterInfo::setCluster(), SiStripClusterInfo::signalOverNoise(), siStripClusterInfo_, SiStripQualBad, siStripQualityToken_, edm::DetSetVector< T >::size(), StripGeomDetUnit::specificTopology(), mathSSE::sqrt(), subdetector, GeomDet::surface(), StripSubdetector::TEC, TrackerTopology::tecSide(), TrackerTopology::tidSide(), reco::MuonTime::timeAtIpInOut, reco::MuonTime::timeAtIpInOutErr, createJobs::tmp, topoToken_, totalNbHits, PbPb_ZMuSkimMuonDPG_cff::tracker, trackerEvent_token_, trackMultiplicityCut_, traj, trajectories_token_, TrajGlbX, TrajGlbY, TrajGlbZ, trajHitValid, TrajLocAngleX, TrajLocAngleY, TrajLocErrX, TrajLocErrY, TrajLocX, TrajLocY, trajTrackAsso_token_, useAllHitsFromTracksWithMissingHits_, useFirstMeas_, useLastMeas_, whatlayer, Bounds::width(), and withinAcceptance.

                                                               {
  //Retrieve tracker topology from geometry
  const TrackerTopology* tTopo = &es.getData(topoToken_);
  siStripClusterInfo_.initEvent(es);

  //  bool DEBUG = false;

  LogDebug("SiStripHitEfficiency:HitEff") << "beginning analyze from HitEff" << endl;

  using namespace edm;
  using namespace reco;
  // Step A: Get Inputs

  int run_nr = e.id().run();
  int ev_nr = e.id().event();
  int bunch_nr = e.bunchCrossing();

  // Luminosity informations
  edm::Handle<LumiScalersCollection> lumiScalers = e.getHandle(scalerToken_);
  edm::Handle<OnlineLuminosityRecord> metaData = e.getHandle(metaDataToken_);

  instLumi = 0;
  PU = 0;
  if (addLumi_) {
    if (lumiScalers.isValid() && !lumiScalers->empty()) {
      if (lumiScalers->begin() != lumiScalers->end()) {
        instLumi = lumiScalers->begin()->instantLumi();
        PU = lumiScalers->begin()->pileup();
      }
    } else if (metaData.isValid()) {
      instLumi = metaData->instLumi();
      PU = metaData->avgPileUp();
    } else {
      edm::LogWarning("SiStripHitEfficiencyWorker") << "could not find a source for the Luminosity and PU";
    }
  }

  // CM
  edm::Handle<edm::DetSetVector<SiStripRawDigi> > commonModeDigis;
  if (addCommonMode_)
    e.getByToken(commonModeToken_, commonModeDigis);

  //CombinatoriaTrack
  edm::Handle<reco::TrackCollection> trackCollectionCKF;
  //edm::InputTag TkTagCKF = conf_.getParameter<edm::InputTag>("combinatorialTracks");
  e.getByToken(combinatorialTracks_token_, trackCollectionCKF);

  edm::Handle<std::vector<Trajectory> > TrajectoryCollectionCKF;
  //edm::InputTag TkTrajCKF = conf_.getParameter<edm::InputTag>("trajectories");
  e.getByToken(trajectories_token_, TrajectoryCollectionCKF);

  edm::Handle<TrajTrackAssociationCollection> trajTrackAssociationHandle;
  e.getByToken(trajTrackAsso_token_, trajTrackAssociationHandle);

  // Clusters
  // get the SiStripClusters from the event
  edm::Handle<edmNew::DetSetVector<SiStripCluster> > theClusters;
  //e.getByLabel("siStripClusters", theClusters);
  e.getByToken(clusters_token_, theClusters);

  //get tracker geometry
  edm::ESHandle<TrackerGeometry> tracker = es.getHandle(geomToken_);
  const TrackerGeometry* tkgeom = &(*tracker);

  //get Cluster Parameter Estimator
  //std::string cpe = conf_.getParameter<std::string>("StripCPE");
  edm::ESHandle<StripClusterParameterEstimator> parameterestimator = es.getHandle(cpeToken_);
  const StripClusterParameterEstimator& stripcpe(*parameterestimator);

  // get the SiStripQuality records
  edm::ESHandle<SiStripQuality> SiStripQuality_ = es.getHandle(siStripQualityToken_);

  const MagneticField* magField_ = &es.getData(magFieldToken_);

  // get the list of module IDs with FED-detected errors
  //  - In Aug-2023, the data format was changed from DetIdCollection to DetIdVector.
  //  - To provide some level of backward-compatibility,
  //    the plugin checks for both types giving preference to the new format.
  //  - If only the old format is available, the collection is
  //    converted to the new format, then used downstream.
  auto const& fedErrorIdsCol_h = e.getHandle(digisCol_token_);
  auto const& fedErrorIdsVec_h = e.getHandle(digisVec_token_);
  if (not fedErrorIdsCol_h.isValid() and not fedErrorIdsVec_h.isValid()) {
    throw cms::Exception("InvalidProductSiStripDetIdsWithFEDErrors")
        << "no valid product for SiStrip DetIds with FED errors (see parameter \"siStripDigis\"), "
           "neither for new format (DetIdVector) nor old format (DetIdCollection)";
  }
  auto const& fedErrorIds = fedErrorIdsVec_h.isValid() ? *fedErrorIdsVec_h : fedErrorIdsCol_h->as_vector();

  edm::ESHandle<MeasurementTracker> measurementTrackerHandle = es.getHandle(measurementTkToken_);

  edm::Handle<MeasurementTrackerEvent> measurementTrackerEvent;
  //e.getByLabel("MeasurementTrackerEvent", measurementTrackerEvent);
  e.getByToken(trackerEvent_token_, measurementTrackerEvent);

  const MeasurementEstimator* estimator = &es.getData(chi2MeasurementEstimatorToken_);
  const Propagator* thePropagator = &es.getData(propagatorToken_);

  events++;

  // *************** SiStripCluster Collection
  const edmNew::DetSetVector<SiStripCluster>& input = *theClusters;

  //go through clusters to write out global position of good clusters for the layer understudy for comparison
  // Loop through clusters just to print out locations
  // Commented out to avoid discussion, should really be deleted.
  /*
  for (edmNew::DetSetVector<SiStripCluster>::const_iterator DSViter = input.begin(); DSViter != input.end(); DSViter++) {
    // DSViter is a vector of SiStripClusters located on a single module
    unsigned int ClusterId = DSViter->id();
    DetId ClusterDetId(ClusterId);
    const StripGeomDetUnit * stripdet=(const StripGeomDetUnit*)tkgeom->idToDetUnit(ClusterDetId);
    
    edmNew::DetSet<SiStripCluster>::const_iterator begin=DSViter->begin();
    edmNew::DetSet<SiStripCluster>::const_iterator end  =DSViter->end();
    for(edmNew::DetSet<SiStripCluster>::const_iterator iter=begin;iter!=end;++iter) {
      //iter is a single SiStripCluster
      StripClusterParameterEstimator::LocalValues parameters=stripcpe.localParameters(*iter,*stripdet);
      
      const Surface* surface;
      surface = &(tracker->idToDet(ClusterDetId)->surface());
      LocalPoint lp = parameters.first;
      GlobalPoint gp = surface->toGlobal(lp);
      unsigned int layer = ::checkLayer(ClusterId, tTopo);
            if(DEBUG) LOGPRINT << "Found hit in cluster collection layer = " << layer << " with id = " << ClusterId << "   local X position = " << lp.x() << " +- " << sqrt(parameters.second.xx()) << "   matched/stereo/rphi = " << ((ClusterId & 0x3)==0) << "/" << ((ClusterId & 0x3)==1) << "/" << ((ClusterId & 0x3)==2) << endl;
    }
  }
  */

  // Tracking
  const reco::TrackCollection* tracksCKF = trackCollectionCKF.product();
  LogDebug("SiStripHitEfficiency:HitEff") << "number ckf tracks found = " << tracksCKF->size() << endl;
  //if (tracksCKF->size() == 1 ){
  if (!tracksCKF->empty()) {
    if (cutOnTracks_ && (tracksCKF->size() >= trackMultiplicityCut_))
      return;
    if (cutOnTracks_)
      LogDebug("SiStripHitEfficiency:HitEff")
          << "starting checking good event with < " << trackMultiplicityCut_ << " tracks" << endl;

    EventTrackCKF++;

#ifdef ExtendedCALIBTree
    //get dEdx info if available
    edm::Handle<ValueMap<DeDxData> > dEdxUncalibHandle;
    if (e.getByLabel("dedxMedianCTF", dEdxUncalibHandle)) {
      const ValueMap<DeDxData> dEdxTrackUncalib = *dEdxUncalibHandle.product();

      reco::TrackRef itTrack = reco::TrackRef(trackCollectionCKF, 0);
      dedx = dEdxTrackUncalib[itTrack].dEdx();
      dedxNOM = dEdxTrackUncalib[itTrack].numberOfMeasurements();
    } else {
      dedx = -999.0;
      dedxNOM = -999;
    }

    //get muon and ecal timing info if available
    edm::Handle<MuonCollection> muH;
    if (e.getByLabel("muonsWitht0Correction", muH)) {
      const MuonCollection& muonsT0 = *muH.product();
      if (!muonsT0.empty()) {
        MuonTime mt0 = muonsT0[0].time();
        timeDT = mt0.timeAtIpInOut;
        timeDTErr = mt0.timeAtIpInOutErr;
        timeDTDOF = mt0.nDof;

        bool hasCaloEnergyInfo = muonsT0[0].isEnergyValid();
        if (hasCaloEnergyInfo)
          timeECAL = muonsT0[0].calEnergy().ecal_time;
      }
    } else {
      timeDT = -999.0;
      timeDTErr = -999.0;
      timeDTDOF = -999;
      timeECAL = -999.0;
    }

#endif
    // actually should do a loop over all the tracks in the event here

    // Looping over traj-track associations to be able to get traj & track informations
    for (TrajTrackAssociationCollection::const_iterator it = trajTrackAssociationHandle->begin();
         it != trajTrackAssociationHandle->end();
         it++) {
      edm::Ref<std::vector<Trajectory> > itraj = it->key;
      reco::TrackRef itrack = it->val;

      // for each track, fill some variables such as number of hits and momentum
      nHits = itraj->foundHits();
#ifdef ExtendedCALIBTree
      nLostHits = itraj->lostHits();
      chi2 = (itraj->chiSquared() / itraj->ndof());
      p = itraj->lastMeasurement().updatedState().globalMomentum().mag();
#endif
      pT = sqrt((itraj->lastMeasurement().updatedState().globalMomentum().x() *
                 itraj->lastMeasurement().updatedState().globalMomentum().x()) +
                (itraj->lastMeasurement().updatedState().globalMomentum().y() *
                 itraj->lastMeasurement().updatedState().globalMomentum().y()));

      // track quality
      highPurity = itrack->quality(reco::TrackBase::TrackQuality::highPurity);

      std::vector<TrajectoryMeasurement> TMeas = itraj->measurements();
      totalNbHits += int(TMeas.size());
      vector<TrajectoryMeasurement>::iterator itm;
      double xloc = 0.;
      double yloc = 0.;
      double xErr = 0.;
      double yErr = 0.;
      double angleX = -999.;
      double angleY = -999.;
      double xglob, yglob, zglob;

      // Check whether the trajectory has some missing hits
      bool hasMissingHits = false;
      int previous_layer = 999;
      vector<unsigned int> missedLayers;
      for (const auto& itm : TMeas) {
        auto theHit = itm.recHit();
        unsigned int iidd = theHit->geographicalId().rawId();
        int layer = ::checkLayer(iidd, tTopo);
        int missedLayer = (layer + 1);
        int previousMissedLayer = (layer + 2);
        int diffPreviousLayer = (layer - previous_layer);
        if (doMissingHitsRecovery_) {
          //Layers from TIB + TOB
          if (diffPreviousLayer == -2 && missedLayer > k_LayersStart && missedLayer < k_LayersAtTOBEnd) {
            missHitPerLayer[missedLayer] += 1;
            hasMissingHits = true;
          }
          //Layers from TID
          else if (diffPreviousLayer == -2 && (missedLayer > k_LayersAtTOBEnd + 1 && missedLayer <= k_LayersAtTIDEnd)) {
            missHitPerLayer[missedLayer] += 1;
            hasMissingHits = true;
          }
          //Layers from TEC
          else if (diffPreviousLayer == -2 && missedLayer > k_LayersAtTIDEnd && missedLayer <= k_LayersAtTECEnd) {
            missHitPerLayer[missedLayer] += 1;
            hasMissingHits = true;
          }

          //##### TID Layer 11 in transition TID -> TIB : layer is in TIB, previous layer  = 12
          if ((layer > k_LayersStart && layer <= k_LayersAtTIBEnd) && (previous_layer == 12)) {
            missHitPerLayer[11] += 1;
            hasMissingHits = true;
          }

          //##### TEC Layer 14 in transition TEC -> TOB : layer is in TOB, previous layer =  15
          if ((layer > k_LayersAtTIBEnd && layer <= k_LayersAtTOBEnd) && (previous_layer == 15)) {
            missHitPerLayer[14] += 1;
            hasMissingHits = true;
          }

          //####### Consecutive missing hits case #######

          //##### Layers from TIB + TOB
          if (diffPreviousLayer == -3 && missedLayer > k_LayersStart && missedLayer < k_LayersAtTOBEnd &&
              previousMissedLayer > k_LayersStart && previousMissedLayer < k_LayersAtTOBEnd) {
            missHitPerLayer[missedLayer] += 1;
            missHitPerLayer[previousMissedLayer] += 1;
            hasMissingHits = true;
          }

          //##### Layers from TEC
          else if (diffPreviousLayer == -3 && missedLayer > k_LayersAtTIDEnd && missedLayer <= k_LayersAtTECEnd &&
                   previousMissedLayer > k_LayersAtTIDEnd && previousMissedLayer <= k_LayersAtTECEnd) {
            missHitPerLayer[missedLayer] += 1;
            missHitPerLayer[previousMissedLayer] += 1;
            hasMissingHits = true;
          }
        }

        if (theHit->getType() == TrackingRecHit::Type::missing)
          hasMissingHits = true;

        if (hasMissingHits)
          missedLayers.push_back(layer);
        previous_layer = layer;
      }

      // Loop on each measurement and take it into consideration
      //--------------------------------------------------------
      unsigned int prev_TKlayers = 0;
      for (itm = TMeas.begin(); itm != TMeas.end(); itm++) {
        auto theInHit = (*itm).recHit();

        LogDebug("SiStripHitEfficiency:HitEff") << "theInHit is valid = " << theInHit->isValid() << endl;

        unsigned int iidd = theInHit->geographicalId().rawId();
        bool foundConsMissingHits = false;
        unsigned int TKlayers = ::checkLayer(iidd, tTopo);
        LogDebug("SiStripHitEfficiency:HitEff") << "TKlayer from trajectory: " << TKlayers << "  from module = " << iidd
                                                << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/"
                                                << ((iidd & 0x3) == 1) << "/" << ((iidd & 0x3) == 2) << endl;

        // Test first and last points of the trajectory
        // the list of measurements starts from outer layers  !!! This could change -> should add a check
        bool isFirstMeas = (itm == (TMeas.end() - 1));
        bool isLastMeas = (itm == (TMeas.begin()));

        if (!useFirstMeas_ && isFirstMeas)
          continue;
        if (!useLastMeas_ && isLastMeas)
          continue;

        // In case of missing hit in the track, check whether to use the other hits or not.
        if (hasMissingHits && theInHit->getType() != TrackingRecHit::Type::missing &&
            !useAllHitsFromTracksWithMissingHits_)
          continue;

        // If Trajectory measurement from TOB 6 or TEC 9, skip it because it's always valid they are filled later
        if (TKlayers == 10 || TKlayers == 22) {
          LogDebug("SiStripHitEfficiency:HitEff") << "skipping original TM for TOB 6 or TEC 9" << endl;
          continue;
        }

        // Make vector of TrajectoryAtInvalidHits to hold the trajectories
        std::vector<TrajectoryAtInvalidHit> TMs;

        // Make AnalyticalPropagator to use in TAVH constructor
        AnalyticalPropagator propagator(magField_, anyDirection);

        // for double sided layers check both sensors--if no hit was found on either sensor surface,
        // the trajectory measurements only have one invalid hit entry on the matched surface
        // so get the TrajectoryAtInvalidHit for both surfaces and include them in the study
        if (::isDoubleSided(iidd, tTopo) && ((iidd & 0x3) == 0)) {
          // do hit eff check twice--once for each sensor
          //add a TM for each surface
          TMs.push_back(TrajectoryAtInvalidHit(*itm, tTopo, tkgeom, propagator, 1));
          TMs.push_back(TrajectoryAtInvalidHit(*itm, tTopo, tkgeom, propagator, 2));
        } else if (::isDoubleSided(iidd, tTopo) && (!::check2DPartner(iidd, TMeas))) {
          // if only one hit was found the trajectory measurement is on that sensor surface, and the other surface from
          // the matched layer should be added to the study as well
          TMs.push_back(TrajectoryAtInvalidHit(*itm, tTopo, tkgeom, propagator, 1));
          TMs.push_back(TrajectoryAtInvalidHit(*itm, tTopo, tkgeom, propagator, 2));
          LogDebug("SiStripHitEfficiency:HitEff") << " found a hit with a missing partner";
        } else {
          //only add one TM for the single surface and the other will be added in the next iteration
          TMs.push_back(TrajectoryAtInvalidHit(*itm, tTopo, tkgeom, propagator));
        }
        bool missingHitAdded = false;

        vector<TrajectoryMeasurement> tmpTmeas, prev_tmpTmeas;
        unsigned int misLayer = TKlayers + 1;
        unsigned int previousMisLayer = TKlayers + 2;
        //Use bool doMissingHitsRecovery to add possible missing hits based on actual/previous hit
        if (doMissingHitsRecovery_) {
          if (int(TKlayers) - int(prev_TKlayers) == -2) {
            const DetLayer* detlayer = itm->layer();
            const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
            const TrajectoryStateOnSurface tsos = itm->updatedState();
            std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, *thePropagator, *estimator);

            if (misLayer > k_LayersAtTIDEnd && misLayer < k_LayersAtTECEnd) {  //TEC
              std::vector<ForwardDetLayer const*> negTECLayers =
                  measurementTrackerHandle->geometricSearchTracker()->negTecLayers();
              std::vector<ForwardDetLayer const*> posTECLayers =
                  measurementTrackerHandle->geometricSearchTracker()->posTecLayers();
              const DetLayer* tecLayerneg = negTECLayers[misLayer - k_LayersAtTIDEnd - 1];
              const DetLayer* tecLayerpos = posTECLayers[misLayer - k_LayersAtTIDEnd - 1];
              if (tTopo->tecSide(iidd) == 1) {
                tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, *thePropagator, *estimator);
              } else if (tTopo->tecSide(iidd) == 2) {
                tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, *thePropagator, *estimator);
              }
            }

            else if (misLayer == (k_LayersAtTIDEnd - 1) ||
                     misLayer == k_LayersAtTIDEnd) {  // This is for TID layers 12 and 13

              std::vector<ForwardDetLayer const*> negTIDLayers =
                  measurementTrackerHandle->geometricSearchTracker()->negTidLayers();
              std::vector<ForwardDetLayer const*> posTIDLayers =
                  measurementTrackerHandle->geometricSearchTracker()->posTidLayers();
              const DetLayer* tidLayerneg = negTIDLayers[misLayer - k_LayersAtTOBEnd - 1];
              const DetLayer* tidLayerpos = posTIDLayers[misLayer - k_LayersAtTOBEnd - 1];

              if (tTopo->tidSide(iidd) == 1) {
                tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, *thePropagator, *estimator);
              } else if (tTopo->tidSide(iidd) == 2) {
                tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, *thePropagator, *estimator);
              }
            }

            if (misLayer > k_LayersStart && misLayer < k_LayersAtTOBEnd) {  // Barrel

              std::vector<BarrelDetLayer const*> barrelTIBLayers =
                  measurementTrackerHandle->geometricSearchTracker()->tibLayers();
              std::vector<BarrelDetLayer const*> barrelTOBLayers =
                  measurementTrackerHandle->geometricSearchTracker()->tobLayers();

              if (misLayer > k_LayersStart && misLayer <= k_LayersAtTIBEnd) {
                const DetLayer* tibLayer = barrelTIBLayers[misLayer - k_LayersStart - 1];
                tmpTmeas = layerMeasurements.measurements(*tibLayer, tsos, *thePropagator, *estimator);
              } else if (misLayer > k_LayersAtTIBEnd && misLayer < k_LayersAtTOBEnd) {
                const DetLayer* tobLayer = barrelTOBLayers[misLayer - k_LayersAtTIBEnd - 1];
                tmpTmeas = layerMeasurements.measurements(*tobLayer, tsos, *thePropagator, *estimator);
              }
            }
          }
          if ((int(TKlayers) > k_LayersStart && int(TKlayers) <= k_LayersAtTIBEnd) && int(prev_TKlayers) == 12) {
            const DetLayer* detlayer = itm->layer();
            const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
            const TrajectoryStateOnSurface tsos = itm->updatedState();
            std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, *thePropagator, *estimator);
            std::vector<ForwardDetLayer const*> negTIDLayers =
                measurementTrackerHandle->geometricSearchTracker()->negTidLayers();
            std::vector<ForwardDetLayer const*> posTIDLayers =
                measurementTrackerHandle->geometricSearchTracker()->posTidLayers();

            const DetLayer* tidLayerneg = negTIDLayers[k_LayersStart];
            const DetLayer* tidLayerpos = posTIDLayers[k_LayersStart];
            if (tTopo->tidSide(iidd) == 1) {
              tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, *thePropagator, *estimator);
            } else if (tTopo->tidSide(iidd) == 2) {
              tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, *thePropagator, *estimator);
            }
          }

          if ((int(TKlayers) > k_LayersAtTIBEnd && int(TKlayers) <= k_LayersAtTOBEnd) && int(prev_TKlayers) == 15) {
            const DetLayer* detlayer = itm->layer();
            const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
            const TrajectoryStateOnSurface tsos = itm->updatedState();
            std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, *thePropagator, *estimator);

            std::vector<ForwardDetLayer const*> negTECLayers =
                measurementTrackerHandle->geometricSearchTracker()->negTecLayers();
            std::vector<ForwardDetLayer const*> posTECLayers =
                measurementTrackerHandle->geometricSearchTracker()->posTecLayers();

            const DetLayer* tecLayerneg = negTECLayers[k_LayersStart];
            const DetLayer* tecLayerpos = posTECLayers[k_LayersStart];
            if (tTopo->tecSide(iidd) == 1) {
              tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, *thePropagator, *estimator);
            } else if (tTopo->tecSide(iidd) == 2) {
              tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, *thePropagator, *estimator);
            }
          }

          //Test for two consecutive missing hits
          if (int(TKlayers) - int(prev_TKlayers) == -3) {
            foundConsMissingHits = true;
            const DetLayer* detlayer = itm->layer();
            const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
            const TrajectoryStateOnSurface tsos = itm->updatedState();
            std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, *thePropagator, *estimator);

            if (misLayer > k_LayersStart && misLayer <= k_LayersAtTOBEnd && previousMisLayer > k_LayersStart &&
                previousMisLayer <= k_LayersAtTOBEnd) {  //Barrel case
              std::vector<BarrelDetLayer const*> barrelTIBLayers =
                  measurementTrackerHandle->geometricSearchTracker()->tibLayers();
              std::vector<BarrelDetLayer const*> barrelTOBLayers =
                  measurementTrackerHandle->geometricSearchTracker()->tobLayers();
              if (misLayer > k_LayersStart && misLayer < k_LayersAtTIBEnd) {
                const DetLayer* tibLayer = barrelTIBLayers[misLayer - k_LayersStart - 1];
                const DetLayer* prevTibLayer = barrelTIBLayers[previousMisLayer - k_LayersStart - 1];

                tmpTmeas = layerMeasurements.measurements(*tibLayer, tsos, *thePropagator, *estimator);
                prev_tmpTmeas = layerMeasurements.measurements(*prevTibLayer, tsos, *thePropagator, *estimator);
              } else if (misLayer > k_LayersAtTIBEnd && misLayer < k_LayersAtTOBEnd) {
                const DetLayer* tobLayer = barrelTOBLayers[misLayer - k_LayersAtTIBEnd - 1];
                const DetLayer* prevTobLayer = barrelTOBLayers[previousMisLayer - k_LayersAtTIBEnd - 1];
                tmpTmeas = layerMeasurements.measurements(*tobLayer, tsos, *thePropagator, *estimator);
                prev_tmpTmeas = layerMeasurements.measurements(*prevTobLayer, tsos, *thePropagator, *estimator);
              }
            } else if (misLayer > k_LayersAtTIDEnd && misLayer < k_LayersAtTECEnd &&
                       previousMisLayer > k_LayersAtTIDEnd && previousMisLayer < k_LayersAtTECEnd) {  //TEC
              std::vector<ForwardDetLayer const*> negTECLayers =
                  measurementTrackerHandle->geometricSearchTracker()->negTecLayers();
              std::vector<ForwardDetLayer const*> posTECLayers =
                  measurementTrackerHandle->geometricSearchTracker()->posTecLayers();

              const DetLayer* tecLayerneg = negTECLayers[misLayer - k_LayersAtTIDEnd - 1];
              const DetLayer* prevTecLayerneg = negTECLayers[previousMisLayer - k_LayersAtTIDEnd - 1];

              const DetLayer* tecLayerpos = posTECLayers[misLayer - k_LayersAtTIDEnd - 1];
              const DetLayer* prevTecLayerpos = posTECLayers[previousMisLayer - k_LayersAtTIDEnd - 1];

              if (tTopo->tecSide(iidd) == 1) {
                tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, *thePropagator, *estimator);
                prev_tmpTmeas = layerMeasurements.measurements(*prevTecLayerneg, tsos, *thePropagator, *estimator);
              } else if (tTopo->tecSide(iidd) == 2) {
                tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, *thePropagator, *estimator);
                prev_tmpTmeas = layerMeasurements.measurements(*prevTecLayerpos, tsos, *thePropagator, *estimator);
              }
            }
          }
          if (!tmpTmeas.empty() && !foundConsMissingHits) {
            TrajectoryMeasurement TM_tmp(tmpTmeas.back());
            unsigned int iidd_tmp = TM_tmp.recHit()->geographicalId().rawId();
            if (iidd_tmp != 0) {
              LogDebug("SiStripHitEfficiency:HitEff") << " hit actually being added to TM vector";
              if ((!useAllHitsFromTracksWithMissingHits_ || (!useFirstMeas_ && isFirstMeas)))
                TMs.clear();
              if (::isDoubleSided(iidd_tmp, tTopo)) {
                TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator, 1));
                TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator, 2));
              } else
                TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator));
              missingHitAdded = true;
              hitRecoveryCounters[misLayer] += 1;
            }
          }

          if (!tmpTmeas.empty() && !prev_tmpTmeas.empty() &&
              foundConsMissingHits) {  //Found two consecutive missing hits
            TrajectoryMeasurement TM_tmp1(tmpTmeas.back());
            TrajectoryMeasurement TM_tmp2(prev_tmpTmeas.back());
            //Inner and outer hits module IDs
            unsigned int modIdInner = TM_tmp1.recHit()->geographicalId().rawId();
            unsigned int modIdOuter = TM_tmp2.recHit()->geographicalId().rawId();
            bool innerModInactive = false, outerModInactive = false;
            for (const auto& tm : tmpTmeas) {  //Check if inner module is inactive
              unsigned int tmModId = tm.recHit()->geographicalId().rawId();
              if (tmModId == modIdInner && tm.recHit()->getType() == 2) {
                innerModInactive = true;
                break;
              }
            }
            for (const auto& tm : prev_tmpTmeas) {  //Check if outer module is inactive
              unsigned int tmModId = tm.recHit()->geographicalId().rawId();
              if (tmModId == modIdOuter && tm.recHit()->getType() == 2) {
                outerModInactive = true;
                break;  //Found the inactive module
              }
            }

            if (outerModInactive) {  //If outer missing hit is in inactive module, recover the inner one
              if (modIdInner != 0) {
                LogDebug("SiStripHitEfficiency:HitEff") << " hit actually being added to TM vector";
                if ((!useAllHitsFromTracksWithMissingHits_ || (!useFirstMeas_ && isFirstMeas)))
                  TMs.clear();
                if (::isDoubleSided(modIdInner, tTopo)) {
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp1, tTopo, tkgeom, propagator, 1));
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp1, tTopo, tkgeom, propagator, 2));
                } else
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp1, tTopo, tkgeom, propagator));
                missingHitAdded = true;
                hitRecoveryCounters[misLayer] += 1;
              }
            }
            if (innerModInactive) {  //If inner missing hit is in inactive module, recover the outer one
              if (modIdOuter != 0) {
                LogDebug("SiStripHitEfficiency:HitEff") << " hit actually being added to TM vector";
                if ((!useAllHitsFromTracksWithMissingHits_ || (!useFirstMeas_ && isFirstMeas)))
                  TMs.clear();
                if (::isDoubleSided(modIdOuter, tTopo)) {
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp2, tTopo, tkgeom, propagator, 1));
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp2, tTopo, tkgeom, propagator, 2));
                } else
                  TMs.push_back(TrajectoryAtInvalidHit(TM_tmp2, tTopo, tkgeom, propagator));
                missingHitAdded = true;
                hitRecoveryCounters[previousMisLayer] += 1;
              }
            }
          }
        }

        prev_TKlayers = TKlayers;
        if (!useFirstMeas_ && isFirstMeas && !missingHitAdded)
          continue;
        if (!useLastMeas_ && isLastMeas)
          continue;
        bool hitsWithBias = false;
        for (auto ilayer : missedLayers) {
          if (ilayer < TKlayers)
            hitsWithBias = true;
        }
        if (hasMissingHits && theInHit->getType() != TrackingRecHit::Type::missing && !missingHitAdded &&
            hitsWithBias && !useAllHitsFromTracksWithMissingHits_) {
          continue;
        }
        //Now check for tracks at TOB6 and TEC9

        // to make sure we only propagate on the last TOB5 hit check the next entry isn't also in TOB5
        // to avoid bias, make sure the TOB5 hit is valid (an invalid hit on TOB5 could only exist with a valid hit on TOB6)

        bool isValid = theInHit->isValid();
        bool isLast = (itm == (TMeas.end() - 1));
        bool isLastTOB5 = true;
        if (!isLast) {
          if (::checkLayer((++itm)->recHit()->geographicalId().rawId(), tTopo) == 9)
            isLastTOB5 = false;
          else
            isLastTOB5 = true;
          --itm;
        }

        if (TKlayers == 9 && isValid && isLastTOB5) {
          //      if ( TKlayers==9 && itm==TMeas.rbegin()) {
          //      if ( TKlayers==9 && (itm==TMeas.back()) ) {     // to check for only the last entry in the trajectory for propagation
          std::vector<BarrelDetLayer const*> barrelTOBLayers =
              measurementTrackerHandle->geometricSearchTracker()->tobLayers();
          const DetLayer* tob6 = barrelTOBLayers[barrelTOBLayers.size() - 1];
          const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
          const TrajectoryStateOnSurface tsosTOB5 = itm->updatedState();
          auto tmp = layerMeasurements.measurements(*tob6, tsosTOB5, *thePropagator, *estimator);

          if (!tmp.empty()) {
            LogDebug("SiStripHitEfficiency:HitEff") << "size of TM from propagation = " << tmp.size() << endl;

            // take the last of the TMs, which is always an invalid hit
            // if no detId is available, ie detId==0, then no compatible layer was crossed
            // otherwise, use that TM for the efficiency measurement
            TrajectoryMeasurement tob6TM(tmp.back());
            const auto& tob6Hit = tob6TM.recHit();

            if (tob6Hit->geographicalId().rawId() != 0) {
              LogDebug("SiStripHitEfficiency:HitEff") << "tob6 hit actually being added to TM vector" << endl;
              TMs.push_back(TrajectoryAtInvalidHit(tob6TM, tTopo, tkgeom, propagator));
            }
          }
        }

        bool isLastTEC8 = true;
        if (!isLast) {
          if (::checkLayer((++itm)->recHit()->geographicalId().rawId(), tTopo) == 21)
            isLastTEC8 = false;
          else
            isLastTEC8 = true;
          --itm;
        }

        if (TKlayers == 21 && isValid && isLastTEC8) {
          std::vector<const ForwardDetLayer*> posTecLayers =
              measurementTrackerHandle->geometricSearchTracker()->posTecLayers();
          const DetLayer* tec9pos = posTecLayers[posTecLayers.size() - 1];
          std::vector<const ForwardDetLayer*> negTecLayers =
              measurementTrackerHandle->geometricSearchTracker()->negTecLayers();
          const DetLayer* tec9neg = negTecLayers[negTecLayers.size() - 1];
          const LayerMeasurements layerMeasurements{*measurementTrackerHandle, *measurementTrackerEvent};
          const TrajectoryStateOnSurface tsosTEC9 = itm->updatedState();

          // check if track on positive or negative z
          if (!(iidd == StripSubdetector::TEC))
            LogDebug("SiStripHitEfficiency:HitEff") << "there is a problem with TEC 9 extrapolation" << endl;

          //LOGPRINT << " tec9 id = " << iidd << " and side = " << tTopo->tecSide(iidd) << endl;
          vector<TrajectoryMeasurement> tmp;
          if (tTopo->tecSide(iidd) == 1) {
            tmp = layerMeasurements.measurements(*tec9neg, tsosTEC9, *thePropagator, *estimator);
            //LOGPRINT << "on negative side" << endl;
          }
          if (tTopo->tecSide(iidd) == 2) {
            tmp = layerMeasurements.measurements(*tec9pos, tsosTEC9, *thePropagator, *estimator);
            //LOGPRINT << "on positive side" << endl;
          }

          if (!tmp.empty()) {
            // take the last of the TMs, which is always an invalid hit
            // if no detId is available, ie detId==0, then no compatible layer was crossed
            // otherwise, use that TM for the efficiency measurement
            TrajectoryMeasurement tec9TM(tmp.back());
            const auto& tec9Hit = tec9TM.recHit();

            unsigned int tec9id = tec9Hit->geographicalId().rawId();
            LogDebug("SiStripHitEfficiency:HitEff")
                << "tec9id = " << tec9id << " is Double sided = " << ::isDoubleSided(tec9id, tTopo)
                << "  and 0x3 = " << (tec9id & 0x3) << endl;

            if (tec9Hit->geographicalId().rawId() != 0) {
              LogDebug("SiStripHitEfficiency:HitEff") << "tec9 hit actually being added to TM vector" << endl;
              // in tec the hit can be single or doubled sided. whenever the invalid hit at the end of vector of TMs is
              // double sided it is always on the matched surface, so we need to split it into the true sensor surfaces
              if (::isDoubleSided(tec9id, tTopo)) {
                TMs.push_back(TrajectoryAtInvalidHit(tec9TM, tTopo, tkgeom, propagator, 1));
                TMs.push_back(TrajectoryAtInvalidHit(tec9TM, tTopo, tkgeom, propagator, 2));
              } else
                TMs.push_back(TrajectoryAtInvalidHit(tec9TM, tTopo, tkgeom, propagator));
            }
          }  //else LOGPRINT << "tec9 tmp empty" << endl;
        }
        hitTotalCounters[TKlayers] += 1;


        // Modules Constraints

        for (std::vector<TrajectoryAtInvalidHit>::const_iterator TM = TMs.begin(); TM != TMs.end(); ++TM) {
          // --> Get trajectory from combinatedState
          iidd = TM->monodet_id();
          LogDebug("SiStripHitEfficiency:HitEff") << "setting iidd = " << iidd << " before checking efficiency and ";

          xloc = TM->localX();
          yloc = TM->localY();

          angleX = atan(TM->localDxDz());
          angleY = atan(TM->localDyDz());

          TrajLocErrX = 0.0;
          TrajLocErrY = 0.0;

          xglob = TM->globalX();
          yglob = TM->globalY();
          zglob = TM->globalZ();
          xErr = TM->localErrorX();
          yErr = TM->localErrorY();

          TrajGlbX = 0.0;
          TrajGlbY = 0.0;
          TrajGlbZ = 0.0;
          withinAcceptance = TM->withinAcceptance();

          trajHitValid = TM->validHit();
          int TrajStrip = -1;

          // reget layer from iidd here, to account for TOB 6 and TEC 9 TKlayers being off
          TKlayers = ::checkLayer(iidd, tTopo);

          if ((layers == TKlayers) || (layers == 0)) {  // Look at the layer not used to reconstruct the track
            whatlayer = TKlayers;
            LogDebug("SiStripHitEfficiency:HitEff") << "Looking at layer under study" << endl;
            ModIsBad = 2;
            Id = 0;
            SiStripQualBad = 0;
            run = 0;
            event = 0;
            TrajLocX = 0.0;
            TrajLocY = 0.0;
            TrajLocAngleX = -999.0;
            TrajLocAngleY = -999.0;
            ResX = 0.0;
            ResXSig = 0.0;
            ClusterLocX = 0.0;
            ClusterLocY = 0.0;
            ClusterLocErrX = 0.0;
            ClusterLocErrY = 0.0;
            ClusterStoN = 0.0;
            bunchx = 0;
            commonMode = -100;

            // RPhi RecHit Efficiency

            if (!input.empty()) {
              LogDebug("SiStripHitEfficiency:HitEff") << "Checking clusters with size = " << input.size() << endl;
              int nClusters = 0;
              std::vector<std::vector<float> >
                  VCluster_info;  //fill with X residual, X residual pull, local X, sig(X), local Y, sig(Y), StoN
              for (edmNew::DetSetVector<SiStripCluster>::const_iterator DSViter = input.begin(); DSViter != input.end();
                   DSViter++) {
                // DSViter is a vector of SiStripClusters located on a single module
                //if (DEBUG)      LOGPRINT << "the ID from the DSViter = " << DSViter->id() << endl;
                unsigned int ClusterId = DSViter->id();
                if (ClusterId == iidd) {
                  LogDebug("SiStripHitEfficiency:HitEff")
                      << "found  (ClusterId == iidd) with ClusterId = " << ClusterId << " and iidd = " << iidd << endl;
                  DetId ClusterDetId(ClusterId);
                  const StripGeomDetUnit* stripdet = (const StripGeomDetUnit*)tkgeom->idToDetUnit(ClusterDetId);
                  const StripTopology& Topo = stripdet->specificTopology();

                  float hbedge = 0.0;
                  float htedge = 0.0;
                  float hapoth = 0.0;
                  float uylfac = 0.0;
                  float uxlden = 0.0;
                  if (TKlayers >= 11) {
                    const BoundPlane& plane = stripdet->surface();
                    const TrapezoidalPlaneBounds* trapezoidalBounds(
                        dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
                    std::array<const float, 4> const& parameterTrap =
                        (*trapezoidalBounds).parameters();  // el bueno aqui
                    hbedge = parameterTrap[0];
                    htedge = parameterTrap[1];
                    hapoth = parameterTrap[3];
                    uylfac = (htedge - hbedge) / (htedge + hbedge) / hapoth;
                    uxlden = 1 + yloc * uylfac;
                  }

                  // Need to know position of trajectory in strip number for selecting the right APV later
                  if (TrajStrip == -1) {
                    int nstrips = Topo.nstrips();
                    float pitch = stripdet->surface().bounds().width() / nstrips;
                    TrajStrip = xloc / pitch + nstrips / 2.0;
                    // Need additionnal corrections for endcap
                    if (TKlayers >= 11) {
                      float TrajLocXMid = xloc / (1 + (htedge - hbedge) * yloc / (htedge + hbedge) /
                                                          hapoth);  // radialy extrapolated x loc position at middle
                      TrajStrip = TrajLocXMid / pitch + nstrips / 2.0;
                    }
                    //LOGPRINT<<" Layer "<<TKlayers<<" TrajStrip: "<<nstrips<<" "<<pitch<<" "<<TrajStrip<<endl;
                  }

                  for (edmNew::DetSet<SiStripCluster>::const_iterator iter = DSViter->begin(); iter != DSViter->end();
                       ++iter) {
                    //iter is a single SiStripCluster
                    StripClusterParameterEstimator::LocalValues parameters = stripcpe.localParameters(*iter, *stripdet);
                    float res = (parameters.first.x() - xloc);
                    float sigma = ::checkConsistency(parameters, xloc, xErr);
                    // The consistency is probably more accurately measured with the Chi2MeasurementEstimator. To use it
                    // you need a TransientTrackingRecHit instead of the cluster
                    //theEstimator=       new Chi2MeasurementEstimator(30);
                    //const Chi2MeasurementEstimator *theEstimator(100);
                    //theEstimator->estimate(TM->tsos(), TransientTrackingRecHit);

                    if (TKlayers >= 11) {
                      res = parameters.first.x() - xloc / uxlden;  // radialy extrapolated x loc position at middle
                      sigma = abs(res) /
                              sqrt(parameters.second.xx() + xErr * xErr / uxlden / uxlden +
                                   yErr * yErr * xloc * xloc * uylfac * uylfac / uxlden / uxlden / uxlden / uxlden);
                    }

                    siStripClusterInfo_.setCluster(*iter, ClusterId);
                    // signal to noise from SiStripClusterInfo not working in 225. I'll fix this after the interface
                    // redesign in 300 -ku
                    //float cluster_info[7] = {res, sigma, parameters.first.x(), sqrt(parameters.second.xx()), parameters.first.y(), sqrt(parameters.second.yy()), signal_to_noise};
                    std::vector<float> cluster_info;
                    cluster_info.push_back(res);
                    cluster_info.push_back(sigma);
                    cluster_info.push_back(parameters.first.x());
                    cluster_info.push_back(sqrt(parameters.second.xx()));
                    cluster_info.push_back(parameters.first.y());
                    cluster_info.push_back(sqrt(parameters.second.yy()));
                    cluster_info.push_back(siStripClusterInfo_.signalOverNoise());
                    VCluster_info.push_back(cluster_info);
                    nClusters++;
                    LogDebug("SiStripHitEfficiency:HitEff") << "Have ID match. residual = " << VCluster_info.back()[0]
                                                            << "  res sigma = " << VCluster_info.back()[1] << endl;
                    LogDebug("SiStripHitEfficiency:HitEff")
                        << "trajectory measurement compatability estimate = " << (*itm).estimate() << endl;
                    LogDebug("SiStripHitEfficiency:HitEff")
                        << "hit position = " << parameters.first.x() << "  hit error = " << sqrt(parameters.second.xx())
                        << "  trajectory position = " << xloc << "  traj error = " << xErr << endl;
                  }
                }
              }
              float FinalResSig = 1000.0;
              float FinalCluster[7] = {1000.0, 1000.0, 0.0, 0.0, 0.0, 0.0, 0.0};
              if (nClusters > 0) {
                LogDebug("SiStripHitEfficiency:HitEff") << "found clusters > 0" << endl;
                if (nClusters > 1) {
                  //get the smallest one
                  vector<vector<float> >::iterator ires;
                  for (ires = VCluster_info.begin(); ires != VCluster_info.end(); ires++) {
                    if (abs((*ires)[1]) < abs(FinalResSig)) {
                      FinalResSig = (*ires)[1];
                      for (unsigned int i = 0; i < ires->size(); i++) {
                        LogDebug("SiStripHitEfficiency:HitEff")
                            << "filling final cluster. i = " << i << " before fill FinalCluster[i]=" << FinalCluster[i]
                            << " and (*ires)[i] =" << (*ires)[i] << endl;
                        FinalCluster[i] = (*ires)[i];
                        LogDebug("SiStripHitEfficiency:HitEff")
                            << "filling final cluster. i = " << i << " after fill FinalCluster[i]=" << FinalCluster[i]
                            << " and (*ires)[i] =" << (*ires)[i] << endl;
                      }
                    }
                    LogDebug("SiStripHitEfficiency:HitEff")
                        << "iresidual = " << (*ires)[0] << "  isigma = " << (*ires)[1]
                        << "  and FinalRes = " << FinalCluster[0] << endl;
                  }
                } else {
                  FinalResSig = VCluster_info.at(0)[1];
                  for (unsigned int i = 0; i < VCluster_info.at(0).size(); i++) {
                    FinalCluster[i] = VCluster_info.at(0)[i];
                  }
                }
                VCluster_info.clear();
              }

              LogDebug("SiStripHitEfficiency:HitEff")
                  << "Final residual in X = " << FinalCluster[0] << "+-" << (FinalCluster[0] / FinalResSig) << endl;
              LogDebug("SiStripHitEfficiency:HitEff") << "Checking location of trajectory: abs(yloc) = " << abs(yloc)
                                                      << "  abs(xloc) = " << abs(xloc) << endl;
              LogDebug("SiStripHitEfficiency:HitEff")
                  << "Checking location of cluster hit: yloc = " << FinalCluster[4] << "+-" << FinalCluster[5]
                  << "  xloc = " << FinalCluster[2] << "+-" << FinalCluster[3] << endl;
              LogDebug("SiStripHitEfficiency:HitEff") << "Final cluster signal to noise = " << FinalCluster[6] << endl;

              float exclusionWidth = 0.4;
              float TOBexclusion = 0.0;
              float TECexRing5 = -0.89;
              float TECexRing6 = -0.56;
              float TECexRing7 = 0.60;
              //Added by Chris Edelmaier to do TEC bonding exclusion
              int subdetector = ((iidd >> 25) & 0x7);
              int ringnumber = ((iidd >> 5) & 0x7);

              //New TOB and TEC bonding region exclusion zone
              if ((TKlayers >= 5 && TKlayers < 11) ||
                  ((subdetector == 6) && ((ringnumber >= 5) && (ringnumber <= 7)))) {
                //There are only 2 cases that we need to exclude for
                float highzone = 0.0;
                float lowzone = 0.0;
                float higherr = yloc + 5.0 * yErr;
                float lowerr = yloc - 5.0 * yErr;
                if (TKlayers >= 5 && TKlayers < 11) {
                  //TOB zone
                  highzone = TOBexclusion + exclusionWidth;
                  lowzone = TOBexclusion - exclusionWidth;
                } else if (ringnumber == 5) {
                  //TEC ring 5
                  highzone = TECexRing5 + exclusionWidth;
                  lowzone = TECexRing5 - exclusionWidth;
                } else if (ringnumber == 6) {
                  //TEC ring 6
                  highzone = TECexRing6 + exclusionWidth;
                  lowzone = TECexRing6 - exclusionWidth;
                } else if (ringnumber == 7) {
                  //TEC ring 7
                  highzone = TECexRing7 + exclusionWidth;
                  lowzone = TECexRing7 - exclusionWidth;
                }
                //Now that we have our exclusion region, we just have to properly identify it
                if ((highzone <= higherr) && (highzone >= lowerr))
                  withinAcceptance = false;
                if ((lowzone >= lowerr) && (lowzone <= higherr))
                  withinAcceptance = false;
                if ((higherr <= highzone) && (higherr >= lowzone))
                  withinAcceptance = false;
                if ((lowerr >= lowzone) && (lowerr <= highzone))
                  withinAcceptance = false;
              }

              // fill ntuple varibles
              //get global position from module id number iidd
              TrajGlbX = xglob;
              TrajGlbY = yglob;
              TrajGlbZ = zglob;

              TrajLocErrX = xErr;
              TrajLocErrY = yErr;

              Id = iidd;
              run = run_nr;
              event = ev_nr;
              bunchx = bunch_nr;
              //if ( SiStripQuality_->IsModuleBad(iidd) ) {
              if (SiStripQuality_->getBadApvs(iidd) != 0) {
                SiStripQualBad = 1;
                LogDebug("SiStripHitEfficiency:HitEff") << "strip is bad from SiStripQuality" << endl;
              } else {
                SiStripQualBad = 0;
                LogDebug("SiStripHitEfficiency:HitEff") << "strip is good from SiStripQuality" << endl;
              }

              //check for FED-detected errors and include those in SiStripQualBad
              for (unsigned int ii = 0; ii < fedErrorIds.size(); ii++) {
                if (iidd == fedErrorIds[ii].rawId())
                  SiStripQualBad = 1;
              }

              TrajLocX = xloc;
              TrajLocY = yloc;
              TrajLocAngleX = angleX;
              TrajLocAngleY = angleY;
              ResX = FinalCluster[0];
              ResXSig = FinalResSig;
              if (FinalResSig != FinalCluster[1])
                LogDebug("SiStripHitEfficiency:HitEff")
                    << "Problem with best cluster selection because FinalResSig = " << FinalResSig
                    << " and FinalCluster[1] = " << FinalCluster[1] << endl;
              ClusterLocX = FinalCluster[2];
              ClusterLocY = FinalCluster[4];
              ClusterLocErrX = FinalCluster[3];
              ClusterLocErrY = FinalCluster[5];
              ClusterStoN = FinalCluster[6];

              // CM of APV crossed by traj
              if (addCommonMode_)
                if (commonModeDigis.isValid() && TrajStrip >= 0 && TrajStrip <= 768) {
                  edm::DetSetVector<SiStripRawDigi>::const_iterator digiframe = commonModeDigis->find(iidd);
                  if (digiframe != commonModeDigis->end())
                    if ((unsigned)TrajStrip / 128 < digiframe->data.size())
                      commonMode = digiframe->data.at(TrajStrip / 128).adc();
                }

              LogDebug("SiStripHitEfficiency:HitEff") << "before check good" << endl;

              if (FinalResSig < 999.0) {  //could make requirement on track/hit consistency, but for
                //now take anything with a hit on the module
                LogDebug("SiStripHitEfficiency:HitEff")
                    << "hit being counted as good " << FinalCluster[0] << " FinalRecHit " << iidd << "   TKlayers  "
                    << TKlayers << " xloc " << xloc << " yloc  " << yloc << " module " << iidd
                    << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/" << ((iidd & 0x3) == 1) << "/"
                    << ((iidd & 0x3) == 2) << endl;
                ModIsBad = 0;
                traj->Fill();
              } else {
                LogDebug("SiStripHitEfficiency:HitEff")
                    << "hit being counted as bad   ######### Invalid RPhi FinalResX " << FinalCluster[0]
                    << " FinalRecHit " << iidd << "   TKlayers  " << TKlayers << " xloc " << xloc << " yloc  " << yloc
                    << " module " << iidd << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/"
                    << ((iidd & 0x3) == 1) << "/" << ((iidd & 0x3) == 2) << endl;
                ModIsBad = 1;
                traj->Fill();

                LogDebug("SiStripHitEfficiency:HitEff") << " RPhi Error " << sqrt(xErr * xErr + yErr * yErr)
                                                        << " ErrorX " << xErr << " yErr " << yErr << endl;
              }
              LogDebug("SiStripHitEfficiency:HitEff") << "after good location check" << endl;
            }
            LogDebug("SiStripHitEfficiency:HitEff") << "after list of clusters" << endl;
          }
          LogDebug("SiStripHitEfficiency:HitEff") << "After layers=TKLayers if" << endl;
        }
        LogDebug("SiStripHitEfficiency:HitEff") << "After looping over TrajAtValidHit list" << endl;
      }
      LogDebug("SiStripHitEfficiency:HitEff") << "end TMeasurement loop" << endl;
    }
    LogDebug("SiStripHitEfficiency:HitEff") << "end of trajectories loop" << endl;
  }
}

beginJob()

void HitEff::beginJob ( )

overrideprivatevirtual

Reimplemented from edm::one::EDAnalyzerBase.

Definition at line 115 of file HitEff.cc.

References addCommonMode_, addLumi_, bunchx, nano_mu_local_reco_cff::chi2, ClusterLocErrX, ClusterLocErrY, ClusterLocX, ClusterLocY, ClusterStoN, commonMode, compSettings, events, EventTrackCKF, compareTotals::fs, highPurity, Id, instLumi, missHitPerLayer, ModIsBad, nHits, AlCaHLTBitMon_ParallelJobs::p, pT, PU, ResX, ResXSig, run, SiStripQualBad, totalNbHits, tquality, traj, TrajGlbX, TrajGlbY, TrajGlbZ, trajHitValid, TrajLocAngleX, TrajLocAngleY, TrajLocErrX, TrajLocErrY, TrajLocX, TrajLocY, whatlayer, and withinAcceptance.

                      {
  edm::Service<TFileService> fs;
  if (compSettings > 0) {
    edm::LogInfo("SiStripHitEfficiency:HitEff") << "the compressions settings are:" << compSettings << std::endl;
    fs->file().SetCompressionSettings(compSettings);
  }

  traj = fs->make<TTree>("traj", "tree of trajectory positions");
#ifdef ExtendedCALIBTree
  traj->Branch("timeDT", &timeDT, "timeDT/F");
  traj->Branch("timeDTErr", &timeDTErr, "timeDTErr/F");
  traj->Branch("timeDTDOF", &timeDTDOF, "timeDTDOF/I");
  traj->Branch("timeECAL", &timeECAL, "timeECAL/F");
  traj->Branch("dedx", &dedx, "dedx/F");
  traj->Branch("dedxNOM", &dedxNOM, "dedxNOM/I");
  traj->Branch("nLostHits", &nLostHits, "nLostHits/I");
  traj->Branch("chi2", &chi2, "chi2/F");
  traj->Branch("p", &p, "p/F");
#endif
  traj->Branch("TrajGlbX", &TrajGlbX, "TrajGlbX/F");
  traj->Branch("TrajGlbY", &TrajGlbY, "TrajGlbY/F");
  traj->Branch("TrajGlbZ", &TrajGlbZ, "TrajGlbZ/F");
  traj->Branch("TrajLocX", &TrajLocX, "TrajLocX/F");
  traj->Branch("TrajLocY", &TrajLocY, "TrajLocY/F");
  traj->Branch("TrajLocAngleX", &TrajLocAngleX, "TrajLocAngleX/F");
  traj->Branch("TrajLocAngleY", &TrajLocAngleY, "TrajLocAngleY/F");
  traj->Branch("TrajLocErrX", &TrajLocErrX, "TrajLocErrX/F");
  traj->Branch("TrajLocErrY", &TrajLocErrY, "TrajLocErrY/F");
  traj->Branch("ClusterLocX", &ClusterLocX, "ClusterLocX/F");
  traj->Branch("ClusterLocY", &ClusterLocY, "ClusterLocY/F");
  traj->Branch("ClusterLocErrX", &ClusterLocErrX, "ClusterLocErrX/F");
  traj->Branch("ClusterLocErrY", &ClusterLocErrY, "ClusterLocErrY/F");
  traj->Branch("ClusterStoN", &ClusterStoN, "ClusterStoN/F");
  traj->Branch("ResX", &ResX, "ResX/F");
  traj->Branch("ResXSig", &ResXSig, "ResXSig/F");
  traj->Branch("ModIsBad", &ModIsBad, "ModIsBad/i");
  traj->Branch("SiStripQualBad", &SiStripQualBad, "SiStripQualBad/i");
  traj->Branch("withinAcceptance", &withinAcceptance, "withinAcceptance/O");
  traj->Branch("nHits", &nHits, "nHits/I");
  traj->Branch("pT", &pT, "pT/F");
  traj->Branch("highPurity", &highPurity, "highPurity/O");
  traj->Branch("trajHitValid", &trajHitValid, "trajHitValid/i");
  traj->Branch("Id", &Id, "Id/i");
  traj->Branch("run", &run, "run/i");
  traj->Branch("event", &event, "event/i");
  traj->Branch("layer", &whatlayer, "layer/i");
  traj->Branch("tquality", &tquality, "tquality/I");
  traj->Branch("bunchx", &bunchx, "bunchx/I");
  if (addLumi_) {
    traj->Branch("instLumi", &instLumi, "instLumi/F");
    traj->Branch("PU", &PU, "PU/F");
  }
  if (addCommonMode_)
    traj->Branch("commonMode", &commonMode, "commonMode/F");

  events = 0;
  EventTrackCKF = 0;

  totalNbHits = 0;
  missHitPerLayer.resize(k_END_OF_LAYERS, 0);
}

endJob()

void HitEff::endJob ( void  )

overrideprivatevirtual

Reimplemented from edm::one::EDAnalyzerBase.

Definition at line 1180 of file HitEff.cc.

References doMissingHitsRecovery_, events, EventTrackCKF, hitRecoveryCounters, hitTotalCounters, mps_fire::i, LogDebug, missHitPerLayer, totalNbHits, and traj.

Referenced by o2olib.O2ORunMgr::executeJob().

                    {
  traj->GetDirectory()->cd();
  traj->Write();

  LogDebug("SiStripHitEfficiency:HitEff") << " Events Analysed             " << events << endl;
  LogDebug("SiStripHitEfficiency:HitEff") << " Number Of Tracked events    " << EventTrackCKF << endl;

  if (doMissingHitsRecovery_) {
    float totTIB = 0.0;
    float totTOB = 0.0;
    float totTID = 0.0;
    float totTEC = 0.0;

    float totTIBrepro = 0.0;
    float totTOBrepro = 0.0;
    float totTIDrepro = 0.0;
    float totTECrepro = 0.0;

    edm::LogInfo("SiStripHitEfficiency:HitEff") << "Within TIB :";
    for (int i = 0; i <= k_LayersAtTIBEnd; i++) {
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Layer " << i << " has : " << missHitPerLayer[i] << "/" << totalNbHits << " = "
          << (missHitPerLayer[i] * 1.0 / totalNbHits) * 100 << " % of missing hit";
      totTIB += missHitPerLayer[i];
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Removing recovered hits : layer " << i << " has : " << missHitPerLayer[i] - hitRecoveryCounters[i] << "/"
          << totalNbHits << " = " << ((missHitPerLayer[i] - hitRecoveryCounters[i]) * 1.0 / totalNbHits) * 100
          << " % of missing hit";
      totTIBrepro += (missHitPerLayer[i] - hitRecoveryCounters[i]);
    }
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "TOTAL % of missing hits within TIB :" << (totTIB * 1.0 / totalNbHits) * 100 << "%";
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "AFTER repropagation :" << (totTIBrepro * 1.0 / totalNbHits) * 100 << "%";

    edm::LogInfo("SiStripHitEfficiency:HitEff") << "Within TOB :";
    for (int i = k_LayersAtTIBEnd + 1; i <= k_LayersAtTOBEnd; i++) {
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Layer " << i << " has : " << missHitPerLayer[i] << "/" << totalNbHits << " = "
          << (missHitPerLayer[i] * 1.0 / totalNbHits) * 100 << " % of missing hit";
      totTOB += missHitPerLayer[i];
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Removing recovered hits : layer " << i << " has : " << missHitPerLayer[i] - hitRecoveryCounters[i] << "/"
          << totalNbHits << " = " << ((missHitPerLayer[i] - hitRecoveryCounters[i]) * 1.0 / totalNbHits) * 100
          << " % of missing hit";
      totTOBrepro += (missHitPerLayer[i] - hitRecoveryCounters[i]);
    }
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "TOTAL % of missing hits within TOB :" << (totTOB * 1.0 / totalNbHits) * 100 << "%";
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "AFTER repropagation :" << (totTOBrepro * 1.0 / totalNbHits) * 100 << "%";

    edm::LogInfo("SiStripHitEfficiency:HitEff") << "Within TID :";
    for (int i = k_LayersAtTOBEnd + 1; i <= k_LayersAtTIDEnd; i++) {
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Layer " << i << " has : " << missHitPerLayer[i] << "/" << totalNbHits << " = "
          << (missHitPerLayer[i] * 1.0 / totalNbHits) * 100 << " % of missing hit";
      totTID += missHitPerLayer[i];
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Removing recovered hits : layer " << i << " has : " << missHitPerLayer[i] - hitRecoveryCounters[i] << "/"
          << totalNbHits << " = " << ((missHitPerLayer[i] - hitRecoveryCounters[i]) * 1.0 / totalNbHits) * 100
          << " % of missing hit";
      totTIDrepro += (missHitPerLayer[i] - hitRecoveryCounters[i]);
    }
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "TOTAL % of missing hits within TID :" << (totTID * 1.0 / totalNbHits) * 100 << "%";
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "AFTER repropagation :" << (totTIDrepro * 1.0 / totalNbHits) * 100 << "%";

    edm::LogInfo("SiStripHitEfficiency:HitEff") << "Within TEC :";
    for (int i = k_LayersAtTIDEnd + 1; i < k_END_OF_LAYERS; i++) {
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Layer " << i << " has : " << missHitPerLayer[i] << "/" << totalNbHits << " = "
          << (missHitPerLayer[i] * 1.0 / totalNbHits) * 100 << " % of missing hit";
      totTEC += missHitPerLayer[i];
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << "Removing recovered hits : layer " << i << " has : " << missHitPerLayer[i] - hitRecoveryCounters[i] << "/"
          << totalNbHits << " = " << ((missHitPerLayer[i] - hitRecoveryCounters[i]) * 1.0 / totalNbHits) * 100
          << " % of missing hit";
      totTECrepro += (missHitPerLayer[i] - hitRecoveryCounters[i]);
    }
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "TOTAL % of missing hits within TEC :" << (totTEC * 1.0 / totalNbHits) * 100 << "%";
    edm::LogInfo("SiStripHitEfficiency:HitEff")
        << "AFTER repropagation :" << (totTECrepro * 1.0 / totalNbHits) * 100 << "%";

    edm::LogInfo("SiStripHitEfficiency:HitEff") << " Hit recovery summary:";

    for (int ilayer = 0; ilayer < k_END_OF_LAYERS; ilayer++) {
      edm::LogInfo("SiStripHitEfficiency:HitEff")
          << " layer " << ilayer << ": " << hitRecoveryCounters[ilayer] << " / " << hitTotalCounters[ilayer];
    }
  }
}

Member Data Documentation

addCommonMode_

bool HitEff::addCommonMode_

private

Definition at line 72 of file HitEff.h.

Referenced by analyze(), beginJob(), and HitEff().

addLumi_

bool HitEff::addLumi_

private

Definition at line 71 of file HitEff.h.

Referenced by analyze(), beginJob(), and HitEff().

bunchx

unsigned int HitEff::bunchx

private

Definition at line 135 of file HitEff.h.

Referenced by analyze(), and beginJob().

chi2MeasurementEstimatorToken_

const edm::ESGetToken<Chi2MeasurementEstimatorBase, TrackingComponentsRecord> HitEff::chi2MeasurementEstimatorToken_

private

Definition at line 96 of file HitEff.h.

Referenced by analyze().

ClusterLocErrX

float HitEff::ClusterLocErrX

private

Definition at line 126 of file HitEff.h.

Referenced by analyze(), and beginJob().

ClusterLocErrY

float HitEff::ClusterLocErrY

private

Definition at line 126 of file HitEff.h.

Referenced by analyze(), and beginJob().

ClusterLocX

float HitEff::ClusterLocX

private

Definition at line 126 of file HitEff.h.

Referenced by analyze(), and beginJob().

ClusterLocY

float HitEff::ClusterLocY

private

Definition at line 126 of file HitEff.h.

Referenced by analyze(), and beginJob().

clusters_token_

const edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster> > HitEff::clusters_token_

private

Definition at line 83 of file HitEff.h.

Referenced by analyze().

ClusterStoN

float HitEff::ClusterStoN

private

Definition at line 126 of file HitEff.h.

Referenced by analyze(), and beginJob().

combinatorialTracks_token_

const edm::EDGetTokenT<reco::TrackCollection> HitEff::combinatorialTracks_token_

private

Definition at line 80 of file HitEff.h.

Referenced by analyze().

commonMode

float HitEff::commonMode

private

Definition at line 138 of file HitEff.h.

Referenced by analyze(), and beginJob().

commonModeToken_

const edm::EDGetTokenT<edm::DetSetVector<SiStripRawDigi> > HitEff::commonModeToken_

private

Definition at line 67 of file HitEff.h.

Referenced by analyze().

compSettings

int HitEff::compSettings

private

Definition at line 104 of file HitEff.h.

Referenced by beginJob(), and HitEff().

conf_

edm::ParameterSet HitEff::conf_

private

Definition at line 99 of file HitEff.h.

Referenced by HitEff().

cpeToken_

const edm::ESGetToken<StripClusterParameterEstimator, TkStripCPERecord> HitEff::cpeToken_

private

Definition at line 92 of file HitEff.h.

Referenced by analyze().

cutOnTracks_

bool HitEff::cutOnTracks_

private

Definition at line 73 of file HitEff.h.

Referenced by analyze(), and HitEff().

DEBUG

bool HitEff::DEBUG

private

Definition at line 106 of file HitEff.h.

Referenced by HitEff().

digisCol_token_

const edm::EDGetTokenT<DetIdCollection> HitEff::digisCol_token_

private

Definition at line 84 of file HitEff.h.

Referenced by analyze().

digisVec_token_

const edm::EDGetTokenT<DetIdVector> HitEff::digisVec_token_

private

Definition at line 85 of file HitEff.h.

Referenced by analyze().

doMissingHitsRecovery_

bool HitEff::doMissingHitsRecovery_

private

Definition at line 78 of file HitEff.h.

Referenced by analyze(), endJob(), and HitEff().

event

unsigned int HitEff::event

private

Definition at line 135 of file HitEff.h.

Referenced by Types.EventID::cppID(), looper.Looper::process(), and core.AutoHandle.AutoHandle::product().

events

int HitEff::events

private

Definition at line 102 of file HitEff.h.

Referenced by eventsfwlite.Events::__getattr__(), eventsfwlite.Events::__getitem__(), eventsfwlite.Events::__len__(), analyze(), beginJob(), endJob(), looper.Looper::loop(), and looper.Looper::process().

EventTrackCKF

int HitEff::EventTrackCKF

private

Definition at line 102 of file HitEff.h.

Referenced by analyze(), beginJob(), and endJob().

geomToken_

const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> HitEff::geomToken_

private

Definition at line 91 of file HitEff.h.

Referenced by analyze().

highPurity

bool HitEff::highPurity

private

Definition at line 132 of file HitEff.h.

Referenced by analyze(), and beginJob().

hitRecoveryCounters

std::vector<unsigned int> HitEff::hitRecoveryCounters

private

Definition at line 109 of file HitEff.h.

Referenced by analyze(), endJob(), and HitEff().

hitTotalCounters

std::vector<unsigned int> HitEff::hitTotalCounters

private

Definition at line 110 of file HitEff.h.

Referenced by analyze(), endJob(), and HitEff().

Id

unsigned int HitEff::Id

private

Definition at line 129 of file HitEff.h.

Referenced by analyze(), and beginJob().

instLumi

float HitEff::instLumi

private

Definition at line 137 of file HitEff.h.

Referenced by analyze(), and beginJob().

layers

unsigned int HitEff::layers

private

Definition at line 105 of file HitEff.h.

Referenced by analyze(), and HitEff().

magFieldToken_

const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> HitEff::magFieldToken_

private

Definition at line 94 of file HitEff.h.

Referenced by analyze().

measurementTkToken_

const edm::ESGetToken<MeasurementTracker, CkfComponentsRecord> HitEff::measurementTkToken_

private

Definition at line 95 of file HitEff.h.

Referenced by analyze().

metaDataToken_

const edm::EDGetTokenT<OnlineLuminosityRecord> HitEff::metaDataToken_

private

Definition at line 66 of file HitEff.h.

Referenced by analyze().

missHitPerLayer

std::vector<int> HitEff::missHitPerLayer

private

Definition at line 122 of file HitEff.h.

Referenced by analyze(), beginJob(), and endJob().

ModIsBad

unsigned int HitEff::ModIsBad

private

Definition at line 128 of file HitEff.h.

Referenced by analyze(), and beginJob().

nHits

int HitEff::nHits

private

Definition at line 133 of file HitEff.h.

Referenced by analyze(), and beginJob().

propagatorToken_

const edm::ESGetToken<Propagator, TrackingComponentsRecord> HitEff::propagatorToken_

private

Definition at line 97 of file HitEff.h.

Referenced by analyze().

pT

float HitEff::pT

private

Definition at line 134 of file HitEff.h.

Referenced by analyze(), and beginJob().

PU

float HitEff::PU

private

Definition at line 137 of file HitEff.h.

Referenced by analyze(), beginJob(), upgradeWorkflowComponents.UpgradeWorkflow::init(), and upgradeWorkflowComponents.UpgradeWorkflow::setupPU().

ResX

float HitEff::ResX

private

Definition at line 127 of file HitEff.h.

Referenced by analyze(), and beginJob().

ResXSig

float HitEff::ResXSig

private

Definition at line 127 of file HitEff.h.

Referenced by analyze(), and beginJob().

run

unsigned int HitEff::run

private

Definition at line 135 of file HitEff.h.

Referenced by generateEDF.LumiInfo::__str__(), MatrixUtil.InputInfo::__str__(), DTWorkflow.DTWorkflow::all(), analyze(), beginJob(), Types.EventID::cppID(), Types.LuminosityBlockID::cppID(), MatrixUtil.InputInfo::das(), personalPlayback.Playback::do_create_lumi(), personalPlayback.Playback::do_exec(), o2olib.O2OTool::execute(), generateEDF.LumiInfo::fixXingInfo(), MatrixUtil.InputInfo::lumiRanges(), MatrixUtil.InputInfo::queries(), and dqm_interfaces.DirWalkerFile::walk().

scalerToken_

const edm::EDGetTokenT<LumiScalersCollection> HitEff::scalerToken_

private

Definition at line 65 of file HitEff.h.

Referenced by analyze().

siStripClusterInfo_

SiStripClusterInfo HitEff::siStripClusterInfo_

private

Definition at line 69 of file HitEff.h.

Referenced by analyze().

SiStripQualBad

unsigned int HitEff::SiStripQualBad

private

Definition at line 130 of file HitEff.h.

Referenced by analyze(), and beginJob().

siStripQualityToken_

const edm::ESGetToken<SiStripQuality, SiStripQualityRcd> HitEff::siStripQualityToken_

private

Definition at line 93 of file HitEff.h.

Referenced by analyze().

topoToken_

const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> HitEff::topoToken_

private

Definition at line 90 of file HitEff.h.

Referenced by analyze().

totalNbHits

int HitEff::totalNbHits

private

Definition at line 121 of file HitEff.h.

Referenced by analyze(), beginJob(), and endJob().

tquality

int HitEff::tquality

private

Definition at line 136 of file HitEff.h.

Referenced by beginJob().

trackerEvent_token_

const edm::EDGetTokenT<MeasurementTrackerEvent> HitEff::trackerEvent_token_

private

Definition at line 86 of file HitEff.h.

Referenced by analyze().

trackMultiplicityCut_

unsigned int HitEff::trackMultiplicityCut_

private

Definition at line 74 of file HitEff.h.

Referenced by analyze(), and HitEff().

traj

TTree* HitEff::traj

private

Definition at line 101 of file HitEff.h.

Referenced by analyze(), beginJob(), and endJob().

trajectories_token_

const edm::EDGetTokenT<std::vector<Trajectory> > HitEff::trajectories_token_

private

Definition at line 81 of file HitEff.h.

Referenced by analyze().

TrajGlbX

float HitEff::TrajGlbX

private

Definition at line 123 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajGlbY

float HitEff::TrajGlbY

private

Definition at line 123 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajGlbZ

float HitEff::TrajGlbZ

private

Definition at line 123 of file HitEff.h.

Referenced by analyze(), and beginJob().

trajHitValid

unsigned int HitEff::trajHitValid

private

Definition at line 135 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocAngleX

float HitEff::TrajLocAngleX

private

Definition at line 124 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocAngleY

float HitEff::TrajLocAngleY

private

Definition at line 124 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocErrX

float HitEff::TrajLocErrX

private

Definition at line 125 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocErrY

float HitEff::TrajLocErrY

private

Definition at line 125 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocX

float HitEff::TrajLocX

private

Definition at line 124 of file HitEff.h.

Referenced by analyze(), and beginJob().

TrajLocY

float HitEff::TrajLocY

private

Definition at line 124 of file HitEff.h.

Referenced by analyze(), and beginJob().

trajTrackAsso_token_

const edm::EDGetTokenT<TrajTrackAssociationCollection> HitEff::trajTrackAsso_token_

private

Definition at line 82 of file HitEff.h.

Referenced by analyze().

useAllHitsFromTracksWithMissingHits_

bool HitEff::useAllHitsFromTracksWithMissingHits_

private

Definition at line 77 of file HitEff.h.

Referenced by analyze(), and HitEff().

useFirstMeas_

bool HitEff::useFirstMeas_

private

Definition at line 75 of file HitEff.h.

Referenced by analyze(), and HitEff().

useLastMeas_

bool HitEff::useLastMeas_

private

Definition at line 76 of file HitEff.h.

Referenced by analyze(), and HitEff().

whatlayer

unsigned int HitEff::whatlayer

private

Definition at line 107 of file HitEff.h.

Referenced by analyze(), and beginJob().

withinAcceptance

bool HitEff::withinAcceptance

private

Definition at line 131 of file HitEff.h.

Referenced by analyze(), and beginJob().