#include <CSCEfficiency.h>

Inheritance diagram for CSCEfficiency:

Classes
struct	ChamberHistos

struct	StationHistos

Public Member Functions
	CSCEfficiency (const edm::ParameterSet &pset)
	Constructor. More...

virtual	~CSCEfficiency ()
	Destructor. More...

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

ModuleDescription const &	moduleDescription () const

virtual	~EDFilter ()

Public Member Functions inherited from edm::ProducerBase
void	callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)

std::vector< edm::ProductResolverIndex > const &	indiciesForPutProducts (BranchType iBranchType) const

	ProducerBase ()

void	registerProducts (ProducerBase , ProductRegistry , ModuleDescription const &)

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

void	resolvePutIndicies (BranchType iBranchType, std::unordered_multimap< std::string, edm::ProductResolverIndex > const &iIndicies, std::string const &moduleLabel)

virtual	~ProducerBase () noexcept(false)

Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo >	consumesInfo () const

	EDConsumerBase ()

	EDConsumerBase (EDConsumerBase const &)=delete

	EDConsumerBase (EDConsumerBase &&)=default

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::vector< ModuleDescription const * > &modules, 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

bool	registeredToConsumeMany (TypeID const &, BranchType) const

void	updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)

virtual	~EDConsumerBase () noexcept(false)

Private Member Functions
bool	applyTrigger (edm::Handle< edm::TriggerResults > &hltR, const edm::TriggerNames &triggerNames)

virtual void	beginJob ()

void	chamberCandidates (int station, int ring, float phi, std::vector< int > &coupleOfChambers)

bool	checkLocal (double yLocal, double yBoundary, int station, int ring)

void	chooseDirection (CLHEP::Hep3Vector &innerPosition, CLHEP::Hep3Vector &outerPosition)

bool	efficienciesPerChamber (CSCDetId &id, const CSCChamber *cscChamber, FreeTrajectoryState &ftsChamber)

virtual void	endJob ()

double	extrapolate1D (double initPosition, double initDirection, double parameterOfTheLine)

void	fillDigiInfo (edm::Handle< CSCALCTDigiCollection > &alcts, edm::Handle< CSCCLCTDigiCollection > &clcts, edm::Handle< CSCCorrelatedLCTDigiCollection > &correlatedlcts, edm::Handle< CSCWireDigiCollection > &wires, edm::Handle< CSCStripDigiCollection > &strips, edm::Handle< edm::PSimHitContainer > &simhits, edm::Handle< CSCRecHit2DCollection > &rechits, edm::Handle< CSCSegmentCollection > &segments, edm::ESHandle< CSCGeometry > &cscGeom)

void	fillLCT_info (edm::Handle< CSCALCTDigiCollection > &alcts, edm::Handle< CSCCLCTDigiCollection > &clcts, edm::Handle< CSCCorrelatedLCTDigiCollection > &correlatedlcts)

void	fillRechitsSegments_info (edm::Handle< CSCRecHit2DCollection > &rechits, edm::Handle< CSCSegmentCollection > &segments, edm::ESHandle< CSCGeometry > &cscGeom)

void	fillSimhit_info (edm::Handle< edm::PSimHitContainer > &simHits)

void	fillStrips_info (edm::Handle< CSCStripDigiCollection > &strips)

void	fillWG_info (edm::Handle< CSCWireDigiCollection > &wires, edm::ESHandle< CSCGeometry > &cscGeom)

virtual bool	filter (edm::Event &event, const edm::EventSetup &eventSetup)

FreeTrajectoryState	getFromCLHEP (const CLHEP::Hep3Vector &p3, const CLHEP::Hep3Vector &r3, int charge, const AlgebraicSymMatrix66 &cov, const MagneticField *field)

void	getFromFTS (const FreeTrajectoryState &fts, CLHEP::Hep3Vector &p3, CLHEP::Hep3Vector &r3, int &charge, AlgebraicSymMatrix66 &cov)

bool	inSensitiveLocalRegion (double xLocal, double yLocal, int station, int ring)

void	linearExtrapolation (GlobalPoint initialPosition, GlobalVector initialDirection, float zSurface, std::vector< float > &posZY)

double	lineParameter (double initZPosition, double destZPosition, double initZDirection)

TrajectoryStateOnSurface	propagate (FreeTrajectoryState &ftsStart, const BoundPlane &bp)

const Propagator *	propagator (std::string propagatorName) const

bool	recHitSegment_Efficiencies (CSCDetId &cscDetId, const CSCChamber *cscChamber, FreeTrajectoryState &ftsChamber)

bool	recSimHitEfficiency (CSCDetId &id, FreeTrajectoryState &ftsChamber)

void	returnTypes (CSCDetId &id, int &ec, int &st, int &rg, int &ch, int &secondRing)

void	ringCandidates (int station, float absEta, std::map< std::string, bool > &chamberTypes)

bool	stripWire_Efficiencies (CSCDetId &cscDetId, FreeTrajectoryState &ftsChamber)

Private Attributes
edm::EDGetTokenT< CSCALCTDigiCollection >	al_token

TH1F *	ALCTPerEvent

bool	allALCT [2][4][4][(36-1+1)]

bool	allCLCT [2][4][4][(36-1+1)]

bool	allCorrLCT [2][4][4][(36-1+1)]

std::vector< std::pair< LocalPoint, bool > >	allRechits [2][4][4][(36-1+1)][6]

std::vector< std::pair< LocalPoint, LocalVector > >	allSegments [2][4][4][(36-1+1)]

std::vector< std::pair< LocalPoint, int > >	allSimhits [2][4][4][(36-1+1)][6]

std::vector< std::pair< int, float > >	allStrips [2][4][4][(36-1+1)][6]

std::vector< std::pair< std::pair< int, float >, int > >	allWG [2][4][4][(36-1+1)][6]

bool	alongZ

bool	andOr

bool	applyIPangleCuts

struct CSCEfficiency::ChamberHistos	ChHist [2][4][3][36-1+1]

edm::EDGetTokenT< CSCCLCTDigiCollection >	cl_token

TH1F *	CLCTPerEvent

edm::EDGetTokenT< CSCCorrelatedLCTDigiCollection >	co_token

TH1F *	DataFlow

double	distanceFromDeadZone

bool	emptyChambers [2][4][4][(36-1+1)]

bool	getAbsoluteEfficiency

edm::EDGetTokenT< edm::TriggerResults >	ht_token

bool	isBeamdata

bool	isData

bool	isIPdata

double	local_DX_DZ_Max

double	local_DY_DZ_Max

double	local_DY_DZ_Min

bool	magField

double	maxNormChi2

double	maxP

double	minP

unsigned int	minTrackHits

std::vector< std::string >	myTriggers

int	nEventsAnalyzed

bool	passTheEvent

std::vector< int >	pointToTriggers

bool	printalot

unsigned int	printout_NEvents

TH1F *	recHitsPerEvent

edm::EDGetTokenT< CSCRecHit2DCollection >	rh_token

std::string	rootFileName

edm::EDGetTokenT< CSCStripDigiCollection >	sd_token

edm::EDGetTokenT< CSCSegmentCollection >	se_token

TH1F *	segmentsPerEvent

edm::EDGetTokenT< edm::PSimHitContainer >	sh_token

struct CSCEfficiency::StationHistos	StHist [2][4]

TFile *	theFile

MuonServiceProxy *	theService

edm::EDGetTokenT< edm::View< reco::Track > >	tk_token

TH1F *	TriggersFired

bool	useDigis

bool	useTrigger

edm::EDGetTokenT< CSCWireDigiCollection >	wd_token

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

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

Public Types inherited from edm::EDConsumerBase
typedef ProductLabels	Labels

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

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &)

Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	consumes (edm::InputTag const &tag)

EDGetToken	consumes (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	consumes (TypeToGet const &id, edm::InputTag const &tag)

ConsumesCollector	consumesCollector ()
	Use a ConsumesCollector to gather consumes information from helper functions. More...

template<typename ProductType , BranchType B = InEvent>
void	consumesMany ()

void	consumesMany (const TypeToGet &id)

template<BranchType B>
void	consumesMany (const TypeToGet &id)

template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	mayConsume (edm::InputTag const &tag)

EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

Detailed Description

Efficiency calculations Stoyan Stoynev, Northwestern University

Definition at line 114 of file CSCEfficiency.h.

Constructor & Destructor Documentation

CSCEfficiency::CSCEfficiency ( const edm::ParameterSet & pset )

Constructor.

Definition at line 1624 of file CSCEfficiency.cc.

References diMuonSkim_cfi::andOr, RecoEcal_EventContent_cff::ec, FirstCh, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), createfilelist::int, MuonMETValueMapProducer_cff::maxNormChi2, NewRecoMuonValidator_cfi::maxP, alignBH_cfg::minP, HPSPFTaus_cff::minTrackHits, MuonServiceProxy_cff::MuonServiceProxy, NumCh, mps_fire::Path, DTAnalyzerDetailed_cfi::rootFileName, AlCaHLTBitMon_QueryRunRegistry::string, and interactiveExample::theFile.

                                                        {
 
   // const float Xmin = -70;
   //const float Xmax = 70;
   //const int nXbins = int(4.*(Xmax - Xmin));
   const float Ymin = -165;
   const float Ymax = 165;
   const int nYbins = int((Ymax - Ymin)/2);
   const float Layer_min = -0.5;
   const float Layer_max = 9.5;
   const int nLayer_bins = int(Layer_max - Layer_min);
   //
 
   //---- Get the input parameters
   printout_NEvents  = pset.getUntrackedParameter<unsigned int>("printout_NEvents",0);
   rootFileName     = pset.getUntrackedParameter<string>("rootFileName","cscHists.root");
 
   isData  = pset.getUntrackedParameter<bool>("runOnData",true);// 
   isIPdata  = pset.getUntrackedParameter<bool>("IPdata",false);// 
   isBeamdata  = pset.getUntrackedParameter<bool>("Beamdata",false);//
   getAbsoluteEfficiency  = pset.getUntrackedParameter<bool>("getAbsoluteEfficiency",true);//
   useDigis = pset.getUntrackedParameter<bool>("useDigis", true);// 
   distanceFromDeadZone = pset.getUntrackedParameter<double>("distanceFromDeadZone", 10.);// 
   minP = pset.getUntrackedParameter<double>("minP",20.);//
   maxP = pset.getUntrackedParameter<double>("maxP",100.);//
   maxNormChi2 = pset.getUntrackedParameter<double>("maxNormChi2", 3.);//
   minTrackHits = pset.getUntrackedParameter<unsigned int>("minTrackHits",10);//
 
   applyIPangleCuts = pset.getUntrackedParameter<bool>("applyIPangleCuts", false);// 
     local_DY_DZ_Max = pset.getUntrackedParameter<double>("local_DY_DZ_Max",-0.1);//
       local_DY_DZ_Min = pset.getUntrackedParameter<double>("local_DY_DZ_Min",-0.8);//
         local_DX_DZ_Max = pset.getUntrackedParameter<double>("local_DX_DZ_Max",0.2);//
 
   sd_token = consumes<CSCStripDigiCollection>( pset.getParameter<edm::InputTag>("stripDigiTag") );
   wd_token = consumes<CSCWireDigiCollection>( pset.getParameter<edm::InputTag>("wireDigiTag") );
   al_token = consumes<CSCALCTDigiCollection>( pset.getParameter<edm::InputTag>("alctDigiTag") );
   cl_token = consumes<CSCCLCTDigiCollection>( pset.getParameter<edm::InputTag>("clctDigiTag") );
   co_token = consumes<CSCCorrelatedLCTDigiCollection>( pset.getParameter<edm::InputTag>("corrlctDigiTag") );
   rh_token = consumes<CSCRecHit2DCollection>( pset.getParameter<edm::InputTag>("rechitTag") );
   se_token = consumes<CSCSegmentCollection>( pset.getParameter<edm::InputTag>("segmentTag") );
   tk_token = consumes<edm::View<reco::Track> >( pset.getParameter<edm::InputTag>("tracksTag") );
   sh_token = consumes<edm::PSimHitContainer>( pset.getParameter<edm::InputTag>("simHitTag") );
 
 
   edm::ParameterSet serviceParameters = pset.getParameter<edm::ParameterSet>("ServiceParameters");
   // maybe use the service for getting magnetic field, propagators, etc. ...
   theService        = new MuonServiceProxy(serviceParameters);
 
   // Trigger
   useTrigger =  pset.getUntrackedParameter<bool>("useTrigger", false);
 
   ht_token = consumes<edm::TriggerResults>( pset.getParameter<edm::InputTag>("HLTriggerResults") );
 
   myTriggers = pset.getParameter<std::vector <std::string> >("myTriggers");
   andOr =  pset.getUntrackedParameter<bool>("andOr");
   pointToTriggers.clear();
 
 
   //---- set counter to zero
   nEventsAnalyzed = 0;
   //---- set presence of magnetic field
   magField = true;
   //
   std::string Path = "AllChambers/";
   std::string FullName;
   //---- File with output histograms 
   theFile = new TFile(rootFileName.c_str(), "RECREATE");
   theFile->cd();
   //---- Book histograms for the analysis
   char SpecName[50];
   
   sprintf(SpecName,"DataFlow"); 
   DataFlow =  
     new TH1F(SpecName,"Data flow;condition number;entries",40,-0.5,39.5);
   //
   sprintf(SpecName,"TriggersFired"); 
   TriggersFired =
     new TH1F(SpecName,"Triggers fired;trigger number;entries",140,-0.5,139.5);
   //
   int Chan = 50;
   float minChan = -0.5;
   float maxChan = 49.5;
   //
   sprintf(SpecName,"ALCTPerEvent");    
   ALCTPerEvent = new TH1F(SpecName,"ALCTs per event;N digis;entries",Chan,minChan,maxChan);
   //
   sprintf(SpecName,"CLCTPerEvent");    
   CLCTPerEvent = new TH1F(SpecName,"CLCTs per event;N digis;entries",Chan,minChan,maxChan);
   //
     sprintf(SpecName,"recHitsPerEvent");    
   recHitsPerEvent = new TH1F(SpecName,"RecHits per event;N digis;entries",150,-0.5,149.5);
   //
   sprintf(SpecName,"segmentsPerEvent");    
   segmentsPerEvent = new TH1F(SpecName,"segments per event;N digis;entries",Chan,minChan,maxChan);
   //
   //---- Book groups of histograms (for any chamber)
 
   map<std::string,bool>::iterator iter;  
   for(int ec = 0;ec<2;++ec){
     for(int st = 0;st<4;++st){
       theFile->cd();
       sprintf(SpecName,"Stations__E%d_S%d",ec+1, st+1);
       theFile->mkdir(SpecName);
       theFile->cd(SpecName);
 
       //
       sprintf(SpecName,"segmentChi2_ndf_St%d",st+1);
       StHist[ec][st].segmentChi2_ndf = 
     new TH1F(SpecName,"Chi2/ndf of a segment;chi2/ndf;entries",100,0.,20.);
       //
       sprintf(SpecName,"hitsInSegment_St%d",st+1);
       StHist[ec][st].hitsInSegment = 
     new TH1F(SpecName,"Number of hits in a segment;nHits;entries",7,-0.5,6.5);
       //
       Chan = 170;
       minChan = 0.85;
       maxChan = 2.55;
       //
       sprintf(SpecName,"AllSegments_eta_St%d",st+1);
       StHist[ec][st].AllSegments_eta = 
     new TH1F(SpecName,"All segments in eta;eta;entries",Chan,minChan,maxChan);
       //
       sprintf(SpecName,"EfficientSegments_eta_St%d",st+1);
       StHist[ec][st].EfficientSegments_eta = 
     new TH1F(SpecName,"Efficient segments in eta;eta;entries",Chan,minChan,maxChan);
       //
       sprintf(SpecName,"ResidualSegments_St%d",st+1);
       StHist[ec][st].ResidualSegments = 
     new TH1F(SpecName,"Residual (segments);residual,cm;entries",75,0.,15.);
       //
       Chan = 200;
       minChan = -800.;
       maxChan = 800.;
       int Chan2 = 200;
       float minChan2 = -800.;
       float maxChan2 = 800.;
       
       sprintf(SpecName,"EfficientSegments_XY_St%d",st+1);
       StHist[ec][st].EfficientSegments_XY = new TH2F(SpecName,"Efficient segments in XY;X;Y",
                              Chan,minChan,maxChan,Chan2,minChan2,maxChan2);
       sprintf(SpecName,"InefficientSegments_XY_St%d",st+1);
       StHist[ec][st].InefficientSegments_XY = new TH2F(SpecName,"Inefficient segments in XY;X;Y",
                                Chan,minChan,maxChan,Chan2,minChan2,maxChan2);
       //
       Chan = 80;
       minChan = 0;
       maxChan = 3.2;
       sprintf(SpecName,"EfficientALCT_momTheta_St%d",st+1);
       StHist[ec][st].EfficientALCT_momTheta = new TH1F(SpecName,"Efficient ALCT in theta (momentum);theta, rad;entries",
                                Chan,minChan,maxChan);
       //
       sprintf(SpecName,"InefficientALCT_momTheta_St%d",st+1);
       StHist[ec][st].InefficientALCT_momTheta = new TH1F(SpecName,"Inefficient ALCT in theta (momentum);theta, rad;entries",
                              Chan,minChan,maxChan);
       //
       Chan = 160;
       minChan = -3.2;
       maxChan = 3.2;
       sprintf(SpecName,"EfficientCLCT_momPhi_St%d",st+1);
       StHist[ec][st].EfficientCLCT_momPhi = new TH1F(SpecName,"Efficient CLCT in phi (momentum);phi, rad;entries",
                              Chan,minChan,maxChan);
       //
       sprintf(SpecName,"InefficientCLCT_momPhi_St%d",st+1);
       StHist[ec][st].InefficientCLCT_momPhi = new TH1F(SpecName,"Inefficient CLCT in phi (momentum);phi, rad;entries",
                                Chan,minChan,maxChan);
       //
       theFile->cd();
       for(int rg = 0;rg<3;++rg){
     if(0!=st && rg>1){
       continue;
     }
     else if(1==rg && 3==st){
       continue;
     }
     for(int iChamber=FirstCh;iChamber<FirstCh+NumCh;iChamber++){
       if(0!=st && 0==rg && iChamber >18){
         continue;
       }
       theFile->cd();
       sprintf(SpecName,"Chambers__E%d_S%d_R%d_Chamber_%d",ec+1, st+1, rg+1,iChamber);
       theFile->mkdir(SpecName);
       theFile->cd(SpecName);
       //
 
       sprintf(SpecName,"EfficientRechits_inSegment_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientRechits_inSegment = 
         new TH1F(SpecName,"Existing RecHit given a segment;layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"InefficientSingleHits_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientSingleHits = 
         new TH1F(SpecName,"Single RecHits not in the segment;layers (1-6);entries ",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"AllSingleHits_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].AllSingleHits = 
         new TH1F(SpecName,"Single RecHits given a segment; layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"digiAppearanceCount_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].digiAppearanceCount = 
         new TH1F(SpecName,"Digi appearance (no-yes): segment(0,1), ALCT(2,3), CLCT(4,5), CorrLCT(6,7); digi type;entries",
              8,-0.5,7.5);
       //
       Chan = 100;
       minChan = -1.1;
       maxChan = 0.9;
       sprintf(SpecName,"EfficientALCT_dydz_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientALCT_dydz = 
         new TH1F(SpecName,"Efficient ALCT; local dy/dz (ME 3 and 4 flipped);entries",
              Chan, minChan, maxChan);
       //
       sprintf(SpecName,"InefficientALCT_dydz_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientALCT_dydz = 
         new TH1F(SpecName,"Inefficient ALCT; local dy/dz (ME 3 and 4 flipped);entries",
              Chan, minChan, maxChan);
       //
       Chan = 100;
       minChan = -1.;
       maxChan = 1.0;
       sprintf(SpecName,"EfficientCLCT_dxdz_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientCLCT_dxdz = 
         new TH1F(SpecName,"Efficient CLCT; local dxdz;entries",
              Chan, minChan, maxChan);
       //
       sprintf(SpecName,"InefficientCLCT_dxdz_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientCLCT_dxdz = 
         new TH1F(SpecName,"Inefficient CLCT; local dxdz;entries",
              Chan, minChan, maxChan);
       //
       sprintf(SpecName,"EfficientRechits_good_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientRechits_good = 
         new TH1F(SpecName,"Existing RecHit - sensitive area only;layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"EfficientStrips_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientStrips = 
         new TH1F(SpecName,"Existing strip;layer (1-6); entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"EfficientWireGroups_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientWireGroups = 
         new TH1F(SpecName,"Existing WireGroups;layer (1-6); entries ",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"StripWiresCorrelations_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].StripWiresCorrelations = 
         new TH1F(SpecName,"StripWire correlations;; entries ",5,0.5,5.5);
       //
       Chan = 80;
       minChan = 0;
       maxChan = 3.2;
       sprintf(SpecName,"NoWires_momTheta_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].NoWires_momTheta = 
         new TH1F(SpecName,"No wires (all strips present) - in theta (momentum);theta, rad;entries",
              Chan,minChan,maxChan);
       //
       Chan = 160;
       minChan = -3.2;
       maxChan = 3.2;
       sprintf(SpecName,"NoStrips_momPhi_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].NoStrips_momPhi = 
         new TH1F(SpecName,"No strips (all wires present) - in phi (momentum);phi, rad;entries",
              Chan,minChan,maxChan);
       //
       for(int iLayer=0; iLayer<6;iLayer++){
         sprintf(SpecName,"Y_InefficientRecHits_inSegment_Ch%d_L%d",iChamber,iLayer);
         ChHist[ec][st][rg][iChamber-FirstCh].Y_InefficientRecHits_inSegment.push_back
           (new TH1F(SpecName,"Missing RecHit/layer in a segment (local system, whole chamber);Y, cm; entries",
             nYbins,Ymin, Ymax));
         //
         sprintf(SpecName,"Y_EfficientRecHits_inSegment_Ch%d_L%d",iChamber,iLayer);
         ChHist[ec][st][rg][iChamber-FirstCh].Y_EfficientRecHits_inSegment.push_back
           (new TH1F(SpecName,"Efficient (extrapolated from the segment) RecHit/layer in a segment (local system, whole chamber);Y, cm; entries",
             nYbins,Ymin, Ymax));
         //
             Chan = 200;
             minChan = -0.2;
             maxChan = 0.2;
             sprintf(SpecName,"Phi_InefficientRecHits_inSegment_Ch%d_L%d",iChamber,iLayer);
             ChHist[ec][st][rg][iChamber-FirstCh].Phi_InefficientRecHits_inSegment.push_back
               (new TH1F(SpecName,"Missing RecHit/layer in a segment (local system, whole chamber);Phi, rad; entries",
                         Chan, minChan, maxChan));
             //
             sprintf(SpecName,"Phi_EfficientRecHits_inSegment_Ch%d_L%d",iChamber,iLayer);
             ChHist[ec][st][rg][iChamber-FirstCh].Phi_EfficientRecHits_inSegment.push_back
               (new TH1F(SpecName,"Efficient (extrapolated from the segment) in a segment (local system, whole chamber);Phi, rad; entries",
                         Chan, minChan, maxChan));
         
       }
       //
       sprintf(SpecName,"Sim_Rechits_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].SimRechits = 
         new TH1F(SpecName,"Existing RecHit (Sim);layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"Sim_Simhits_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].SimSimhits = 
         new TH1F(SpecName,"Existing SimHit (Sim);layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       /*
       sprintf(SpecName,"Sim_Rechits_each_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].SimRechits_each = 
         new TH1F(SpecName,"Existing RecHit (Sim), each;layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       //
       sprintf(SpecName,"Sim_Simhits_each_Ch%d",iChamber);
       ChHist[ec][st][rg][iChamber-FirstCh].SimSimhits_each = 
         new TH1F(SpecName,"Existing SimHit (Sim), each;layers (1-6);entries",nLayer_bins,Layer_min,Layer_max);
       */
       theFile->cd();
     }
       }
     }
   }
 }

CSCEfficiency::~CSCEfficiency ( )

virtual

Destructor.

Definition at line 1934 of file CSCEfficiency.cc.

References create_public_pileup_plots::bins, gather_cfg::cout, RecoEcal_EventContent_cff::ec, FirstCh, NumCh, and interactiveExample::theFile.

                              {
   if (theService) delete theService; 
   // Write the histos to a file
   theFile->cd();
   //
   char SpecName[60];
   std::vector<float> bins, Efficiency, EffError;
   std::vector<float> eff(2);
 
   //---- loop over chambers
   std::map <std::string, bool> chamberTypes;
   chamberTypes["ME11"] = false;
   chamberTypes["ME12"] = false;
   chamberTypes["ME13"] = false;
   chamberTypes["ME21"] = false;
   chamberTypes["ME22"] = false;
   chamberTypes["ME31"] = false;
   chamberTypes["ME32"] = false;
   chamberTypes["ME41"] = false;
   
   map<std::string,bool>::iterator iter;  
   std::cout<<" Writing proper histogram structure (patience)..."<<std::endl;
   for(int ec = 0;ec<2;++ec){
     for(int st = 0;st<4;++st){
       snprintf(SpecName, sizeof(SpecName), "Stations__E%d_S%d",ec+1, st+1);
       theFile->cd(SpecName);
       StHist[ec][st].segmentChi2_ndf->Write();
       StHist[ec][st].hitsInSegment->Write();
       StHist[ec][st].AllSegments_eta->Write();
       StHist[ec][st].EfficientSegments_eta->Write();
       StHist[ec][st].ResidualSegments->Write();
       StHist[ec][st].EfficientSegments_XY->Write();
       StHist[ec][st].InefficientSegments_XY->Write();
       StHist[ec][st].EfficientALCT_momTheta->Write();
       StHist[ec][st].InefficientALCT_momTheta->Write();
       StHist[ec][st].EfficientCLCT_momPhi->Write();
       StHist[ec][st].InefficientCLCT_momPhi->Write();
       for(int rg = 0;rg<3;++rg){
     if(0!=st && rg>1){
       continue;
     }
     else if(1==rg && 3==st){
       continue;
     }
     for(int iChamber=FirstCh;iChamber<FirstCh+NumCh;iChamber++){
       if(0!=st && 0==rg && iChamber >18){
         continue;
       }
       snprintf(SpecName, sizeof(SpecName), "Chambers__E%d_S%d_R%d_Chamber_%d",ec+1, st+1, rg+1,iChamber);
       theFile->cd(SpecName);
       
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientRechits_inSegment->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].AllSingleHits->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].digiAppearanceCount->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientALCT_dydz->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientALCT_dydz->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientCLCT_dxdz->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientCLCT_dxdz->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].InefficientSingleHits->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientRechits_good->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientStrips->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].StripWiresCorrelations->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].NoWires_momTheta->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].NoStrips_momPhi->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].EfficientWireGroups->Write();
       for(unsigned int iLayer = 0; iLayer< 6; iLayer++){
         ChHist[ec][st][rg][iChamber-FirstCh].Y_InefficientRecHits_inSegment[iLayer]->Write();
         ChHist[ec][st][rg][iChamber-FirstCh].Y_EfficientRecHits_inSegment[iLayer]->Write();
             ChHist[ec][st][rg][iChamber-FirstCh].Phi_InefficientRecHits_inSegment[iLayer]->Write();
             ChHist[ec][st][rg][iChamber-FirstCh].Phi_EfficientRecHits_inSegment[iLayer]->Write();
       }
       ChHist[ec][st][rg][iChamber-FirstCh].SimRechits->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].SimSimhits->Write();
       /*
       ChHist[ec][st][rg][iChamber-FirstCh].SimRechits_each->Write();
       ChHist[ec][st][rg][iChamber-FirstCh].SimSimhits_each->Write();
       */
       //
       theFile->cd(SpecName);
       theFile->cd();
     }
       }
     }
   }
   //
   snprintf(SpecName, sizeof(SpecName), "AllChambers");
   theFile->mkdir(SpecName);
   theFile->cd(SpecName);
   DataFlow->Write(); 
   TriggersFired->Write();
   ALCTPerEvent->Write();
   CLCTPerEvent->Write();
   recHitsPerEvent->Write();
   segmentsPerEvent->Write();
   //
   theFile->cd(SpecName);
   //---- Close the file
   theFile->Close();
 }

Member Function Documentation

bool CSCEfficiency::applyTrigger	(	edm::Handle< edm::TriggerResults > &	hltR,
		const edm::TriggerNames &	triggerNames
	)

private

Definition at line 1542 of file CSCEfficiency.cc.

References edm::HLTGlobalStatus::accept(), diMuonSkim_cfi::andOr, gather_cfg::cout, edm::HLTGlobalStatus::error(), edm::HandleBase::isValid(), edm::TriggerNames::triggerNames(), and edm::HLTGlobalStatus::wasrun().

                                                                       {
   bool triggerPassed = true;
   std::vector<std::string>  hlNames=triggerNames.triggerNames();
   pointToTriggers.clear();
   for(size_t imyT = 0;imyT<myTriggers.size();++imyT){
     for (size_t iT=0; iT<hlNames.size(); ++iT) {
       //std::cout<<" iT = "<<iT<<" hlNames[iT] = "<<hlNames[iT]<<
       //" : wasrun = "<<hltR->wasrun(iT)<<" accept = "<<
       //     hltR->accept(iT)<<" !error = "<< 
       //    !hltR->error(iT)<<std::endl; 
       if(!imyT){       
         if(hltR->wasrun(iT) &&
            hltR->accept(iT) &&
            !hltR->error(iT) ){
            TriggersFired->Fill(iT); 
         }
       }
       if(hlNames[iT]==myTriggers[imyT]){
     pointToTriggers.push_back(iT);
     if(imyT){
       break;
     }
       }
     }
   }
   if(pointToTriggers.size()!=myTriggers.size()){
     pointToTriggers.clear();
     if(printalot){
       std::cout<<" Not all trigger names found - all trigger specifications will be ignored. Check your cfg file!"<<std::endl;
     }
   }
   else{
     if(pointToTriggers.size()){
       if(printalot){
         std::cout<<"The following triggers will be required in the event: "<<std::endl;
         for(size_t imyT =0; imyT <pointToTriggers.size();++imyT){
       std::cout<<"  "<<hlNames[pointToTriggers[imyT]];
         }
         std::cout<<std::endl;
         std::cout<<" in condition (AND/OR) : "<<!andOr<<"/"<<andOr<<std::endl;
       }
     }
   }
 
   if (hltR.isValid()) {
     if(!pointToTriggers.size()){
       if(printalot){
         std::cout<<" No triggers specified in the configuration or all ignored - no trigger information will be considered"<<std::endl;
       }
     }
     for(size_t imyT =0; imyT <pointToTriggers.size();++imyT){
       if(hltR->wasrun(pointToTriggers[imyT]) && 
      hltR->accept(pointToTriggers[imyT]) && 
      !hltR->error(pointToTriggers[imyT]) ){
     triggerPassed = true;
     if(andOr){
       break;
     }
       }
       else{
     triggerPassed = false;
     if(!andOr){
       triggerPassed = false;
       break;
     }
       }
     }
   }
   else{
     if(printalot){
       std::cout<<" TriggerResults handle returns invalid state?! No trigger information will be considered"<<std::endl;
     }
   }
   if(printalot){
     std::cout<<" Trigger passed: "<<triggerPassed<<std::endl;
   }
   return triggerPassed;
 }

void CSCEfficiency::beginJob ( void )

privatevirtual

Reimplemented from edm::EDFilter.

Definition at line 2036 of file CSCEfficiency.cc.

2037 {

2038 }

void CSCEfficiency::chamberCandidates	(	int	station,
		int	ring,
		float	phi,
		std::vector< int > &	coupleOfChambers
	)

private

Definition at line 1006 of file CSCEfficiency.cc.

References gather_cfg::cout, createfilelist::int, and M_PI.

                                                                                                         {
   coupleOfChambers.clear();
   // -pi< phi<+pi
   float phi_zero = 0.;// check! the phi at the "edge" of Ch 1
   float phi_const = 2.*M_PI/36.;
   int last_chamber = 36;
   int first_chamber = 1;
   if(1 != station && 1==ring){ // 18 chambers in the ring
     phi_const*=2;
     last_chamber /= 2;
   }
   if(phi<0.){
     if (printalot) std::cout<<" info: negative phi = "<<phi<<std::endl;
     phi += 2*M_PI;
   }
   float chamber_float = (phi - phi_zero)/phi_const;
   int chamber_int = int(chamber_float);
   if (chamber_float - float(chamber_int) -0.5 <0.){
     if(0!=chamber_int ){
       coupleOfChambers.push_back(chamber_int);
     }
     else{
       coupleOfChambers.push_back(last_chamber);
     }
     coupleOfChambers.push_back(chamber_int+1); 
     
   }
   else{
     coupleOfChambers.push_back(chamber_int+1);
     if(last_chamber!=chamber_int+1){
       coupleOfChambers.push_back(chamber_int+2);
     } 
     else{
       coupleOfChambers.push_back(first_chamber);
     }
   }
   if (printalot) std::cout<<" phi = "<<phi<<" phi_zero = "<<phi_zero<<" phi_const = "<<phi_const<<
     " candidate chambers: first ch = "<<coupleOfChambers[0]<<" second ch = "<<coupleOfChambers[1]<<std::endl;
 }

bool CSCEfficiency::checkLocal	(	double	yLocal,
		double	yBoundary,
		int	station,
		int	ring
	)

private

Definition at line 603 of file CSCEfficiency.cc.

                                                                                     {
 //---- check if it is in a good local region (sensitive area - geometrical and HV boundaries excluded) 
   bool pass = false;
   std::vector <float> deadZoneCenter(6);
   const float deadZoneHalf = 0.32*7/2;// wire spacing * (wires missing + 1)/2
   float cutZone = deadZoneHalf + distanceFromDeadZone;//cm
   //---- hardcoded... not good
   if(station>1 && station<5){
     if(2==ring){
       deadZoneCenter[0]= -162.48 ;
       deadZoneCenter[1] = -81.8744;
       deadZoneCenter[2] = -21.18165;
       deadZoneCenter[3] = 39.51105;
       deadZoneCenter[4] = 100.2939;
       deadZoneCenter[5] = 160.58;
       
       if(yLocal >yBoundary &&
      ((yLocal> deadZoneCenter[0] + cutZone && yLocal< deadZoneCenter[1] - cutZone) ||
       (yLocal> deadZoneCenter[1] + cutZone && yLocal< deadZoneCenter[2] - cutZone) ||
       (yLocal> deadZoneCenter[2] + cutZone && yLocal< deadZoneCenter[3] - cutZone) ||
       (yLocal> deadZoneCenter[3] + cutZone && yLocal< deadZoneCenter[4] - cutZone) ||
       (yLocal> deadZoneCenter[4] + cutZone && yLocal< deadZoneCenter[5] - cutZone))){
     pass = true;
       }
     }
     else if(1==ring){
       if(2==station){
     deadZoneCenter[0]= -95.94 ;
     deadZoneCenter[1] = -27.47;
     deadZoneCenter[2] = 33.67;
     deadZoneCenter[3] = 93.72;
       }
       else if(3==station){
     deadZoneCenter[0]= -85.97 ;
     deadZoneCenter[1] = -36.21;
     deadZoneCenter[2] = 23.68;
     deadZoneCenter[3] = 84.04;
       }
       else if(4==station){
     deadZoneCenter[0]= -75.82;
     deadZoneCenter[1] = -26.14;
     deadZoneCenter[2] = 23.85;
     deadZoneCenter[3] = 73.91;
       }
       if(yLocal >yBoundary &&
      ((yLocal> deadZoneCenter[0] + cutZone && yLocal< deadZoneCenter[1] - cutZone) ||
       (yLocal> deadZoneCenter[1] + cutZone && yLocal< deadZoneCenter[2] - cutZone) ||
       (yLocal> deadZoneCenter[2] + cutZone && yLocal< deadZoneCenter[3] - cutZone))){
     pass = true;
       }
     }
   }
   else if(1==station){
     if(3==ring){
       deadZoneCenter[0]= -83.155 ;
       deadZoneCenter[1] = -22.7401;
       deadZoneCenter[2] = 27.86665;
       deadZoneCenter[3] = 81.005;
       if(yLocal > yBoundary &&
      ((yLocal> deadZoneCenter[0] + cutZone && yLocal< deadZoneCenter[1] - cutZone) ||
       (yLocal> deadZoneCenter[1] + cutZone && yLocal< deadZoneCenter[2] - cutZone) ||
       (yLocal> deadZoneCenter[2] + cutZone && yLocal< deadZoneCenter[3] - cutZone))){
     pass = true;
       }
     }
     else if(2==ring){
       deadZoneCenter[0]= -86.285 ;
       deadZoneCenter[1] = -32.88305;
       deadZoneCenter[2] = 32.867423;
       deadZoneCenter[3] = 88.205;
       if(yLocal > (yBoundary) &&
      ((yLocal> deadZoneCenter[0] + cutZone && yLocal< deadZoneCenter[1] - cutZone) ||
       (yLocal> deadZoneCenter[1] + cutZone && yLocal< deadZoneCenter[2] - cutZone) ||
       (yLocal> deadZoneCenter[2] + cutZone && yLocal< deadZoneCenter[3] - cutZone))){
     pass = true;
       }
     }
     else{
       deadZoneCenter[0]= -81.0;
       deadZoneCenter[1] = 81.0;
       if(yLocal > (yBoundary) &&
      ((yLocal> deadZoneCenter[0] + cutZone && yLocal< deadZoneCenter[1] - cutZone) )){
     pass = true;
       }
     }
   }
   return pass;
 }

void CSCEfficiency::chooseDirection	(	CLHEP::Hep3Vector &	innerPosition,
		CLHEP::Hep3Vector &	outerPosition
	)

private

Definition at line 1484 of file CSCEfficiency.cc.

References allConversions_cfi::dz.

                                                                                                    {
 
   //---- Be careful with trigger conditions too
   if(!isIPdata){
     float dy = outerPosition.y() - innerPosition.y();
     float dz = outerPosition.z() - innerPosition.z();
     if(isBeamdata){
       if(dz>0){
     alongZ = true;
       }
       else{
     alongZ = false;
       }
     }
     else{//cosmics
       if(dy/dz>0){
     alongZ = false;
       }
       else{
     alongZ = true;
       }
     }
   }
 }

bool CSCEfficiency::efficienciesPerChamber	(	CSCDetId &	id,
		const CSCChamber *	cscChamber,
		FreeTrajectoryState &	ftsChamber
	)

private

Definition at line 1047 of file CSCEfficiency.cc.

References gather_cfg::cout, BeamSpotFakeParameters_cfi::dxdz, BeamSpotFakeParameters_cfi::dydz, RecoEcal_EventContent_cff::ec, FreeTrajectoryState::momentum(), MillePedeFileConverter_cfg::out, PV3DBase< T, PVType, FrameType >::phi(), PV3DBase< T, PVType, FrameType >::theta(), GeomDet::toLocal(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

                                                                                                                       {
   int ec, st, rg, ch, secondRing;
   returnTypes(id, ec, st, rg, ch, secondRing);
 
   LocalVector localDir = cscChamber->toLocal(ftsChamber.momentum());
   if(printalot){
     std::cout<<" global dir = "<<ftsChamber.momentum()<<std::endl;
     std::cout<<" local dir = "<<localDir<<std::endl;
     std::cout<<" local theta = "<<localDir.theta()<<std::endl;
   }
   float dxdz = localDir.x()/localDir.z();
   float dydz = localDir.y()/localDir.z();
   if(2==st || 3==st){
     if(printalot){
       std::cout<<"st 3 or 4 ... flip dy/dz"<<std::endl;
     }
     dydz = - dydz;
   }
   if(printalot){
     std::cout<<"dy/dz = "<<dydz<<std::endl;
   }
   // Apply angle cut
   bool out = true;
   if(applyIPangleCuts){
     if(dydz>local_DY_DZ_Max || dydz<local_DY_DZ_Min || fabs(dxdz)>local_DX_DZ_Max){ 
       out = false;
     }
   }
 
   // Segments
   bool firstCondition = allSegments[ec][st][rg][ch].size() ? true : false;
   bool secondCondition = false;
   //---- ME1 is special as usual - ME1a and ME1b are actually one chamber
   if(secondRing>-1){
     secondCondition = allSegments[ec][st][secondRing][ch].size() ? true : false;
   }
   if(firstCondition || secondCondition){
     if(out){
       ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(1);
     }  
   }
   else{
     if(out){
       ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(0);
     } 
   }
 
   if(useDigis){
     // ALCTs
     firstCondition = allALCT[ec][st][rg][ch];
     secondCondition = false;
     if(secondRing>-1){
       secondCondition = allALCT[ec][st][secondRing][ch];
     }
     if(firstCondition || secondCondition){
       if(out){
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(3);
       }
       // always apply partial angle cuts for this kind of histos 
       if(fabs(dxdz)<local_DX_DZ_Max){
     StHist[ec][st].EfficientALCT_momTheta->Fill(ftsChamber.momentum().theta());
     ChHist[ec][st][rg][ch].EfficientALCT_dydz->Fill(dydz);
       }
     }
     else{
       if(out){
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(2);
       }
       if(fabs(dxdz)<local_DX_DZ_Max){ 
     StHist[ec][st].InefficientALCT_momTheta->Fill(ftsChamber.momentum().theta());
     ChHist[ec][st][rg][ch].InefficientALCT_dydz->Fill(dydz);
       }
       if(printalot){
     std::cout<<" missing ALCT (dy/dz = "<<dydz<<")";
     printf("\t\tendcap/station/ring/chamber: %i/%i/%i/%i\n",ec+1,st+1,rg+1,ch+1);
       }
     }
     
     // CLCTs
     firstCondition = allCLCT[ec][st][rg][ch];
     secondCondition = false;
     if(secondRing>-1){
       secondCondition = allCLCT[ec][st][secondRing][ch];
     }
     if(firstCondition || secondCondition){
       if(out){
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(5);
       }
       if(dydz<local_DY_DZ_Max && dydz>local_DY_DZ_Min){
     StHist[ec][st].EfficientCLCT_momPhi->Fill(ftsChamber.momentum().phi() );// - phi chamber...
     ChHist[ec][st][rg][ch].EfficientCLCT_dxdz->Fill(dxdz);
       }
     }
     else{
       if(out){
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(4);
       }
       if(dydz<local_DY_DZ_Max && dydz>local_DY_DZ_Min){ 
     StHist[ec][st].InefficientCLCT_momPhi->Fill(ftsChamber.momentum().phi());// - phi chamber...
     ChHist[ec][st][rg][ch].InefficientCLCT_dxdz->Fill(dxdz);
       }
       if(printalot){
     std::cout<<" missing CLCT  (dx/dz = "<<dxdz<<")";
     printf("\t\tendcap/station/ring/chamber: %i/%i/%i/%i\n",ec+1,st+1,rg+1,ch+1);
       }
     }
     if(out){
       // CorrLCTs
       firstCondition = allCorrLCT[ec][st][rg][ch];
       secondCondition = false;
       if(secondRing>-1){
     secondCondition = allCorrLCT[ec][st][secondRing][ch];
       }
       if(firstCondition || secondCondition){
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(7);
       }
       else{
     ChHist[ec][st][rg][ch].digiAppearanceCount->Fill(6);
       }
     }
   }
   return out;
 }

void CSCEfficiency::endJob ( void )

privatevirtual

Reimplemented from edm::EDFilter.

Definition at line 2042 of file CSCEfficiency.cc.

References DEFINE_FWK_MODULE.

2042 {

2043 }

double CSCEfficiency::extrapolate1D	(	double	initPosition,
		double	initDirection,
		double	parameterOfTheLine
	)

private

Definition at line 1474 of file CSCEfficiency.cc.

                                                                                                        {
   double extrapolatedPosition = initPosition + initDirection*parameterOfTheLine;
   return extrapolatedPosition;
 }

void CSCEfficiency::fillDigiInfo	(	edm::Handle< CSCALCTDigiCollection > &	alcts,
		edm::Handle< CSCCLCTDigiCollection > &	clcts,
		edm::Handle< CSCCorrelatedLCTDigiCollection > &	correlatedlcts,
		edm::Handle< CSCWireDigiCollection > &	wires,
		edm::Handle< CSCStripDigiCollection > &	strips,
		edm::Handle< edm::PSimHitContainer > &	simhits,
		edm::Handle< CSCRecHit2DCollection > &	rechits,
		edm::Handle< CSCSegmentCollection > &	segments,
		edm::ESHandle< CSCGeometry > &	cscGeom
	)

private

Definition at line 692 of file CSCEfficiency.cc.

References NumCh.

                                                    {
   for(int iE=0;iE<2;iE++){
     for(int iS=0;iS<4;iS++){
       for(int iR=0;iR<4;iR++){
     for(int iC=0;iC<NumCh;iC++){
       allSegments[iE][iS][iR][iC].clear(); 
       allCLCT[iE][iS][iR][iC] = allALCT[iE][iS][iR][iC] = allCorrLCT[iE][iS][iR][iC] = false;
       for(int iL=0;iL<6;iL++){
         allStrips[iE][iS][iR][iC][iL].clear(); 
         allWG[iE][iS][iR][iC][iL].clear();
         allRechits[iE][iS][iR][iC][iL].clear(); 
         allSimhits[iE][iS][iR][iC][iL].clear(); 
       }
     }
       }
     }
   }
   //
   if(useDigis){
     fillLCT_info(alcts, clcts, correlatedlcts);
     fillWG_info(wires, cscGeom);
     fillStrips_info(strips);
   }
   fillRechitsSegments_info(rechits, segments, cscGeom);
   if(!isData){
     fillSimhit_info(simhits);
   }
 }

void CSCEfficiency::fillLCT_info	(	edm::Handle< CSCALCTDigiCollection > &	alcts,
		edm::Handle< CSCCLCTDigiCollection > &	clcts,
		edm::Handle< CSCCorrelatedLCTDigiCollection > &	correlatedlcts
	)

private

Definition at line 730 of file CSCEfficiency.cc.

References FirstCh, and plotBeamSpotDB::last.

                                                                             {
   //---- ALCTDigis
   int nSize = 0;
   for (CSCALCTDigiCollection::DigiRangeIterator j=alcts->begin(); j!=alcts->end(); j++) {
     ++nSize;
     const CSCDetId& id = (*j).first;
     const CSCALCTDigiCollection::Range& range =(*j).second;
     for (CSCALCTDigiCollection::const_iterator digiIt =
        range.first; digiIt!=range.second;
      ++digiIt){
       // Valid digi in the chamber (or in neighbouring chamber) 
       if((*digiIt).isValid()){
     allALCT[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh] = true;
       }
     }// for digis in layer
   }// end of for (j=...
   ALCTPerEvent->Fill(nSize);
   //---- CLCTDigis
   nSize = 0;
   for (CSCCLCTDigiCollection::DigiRangeIterator j=clcts->begin(); j!=clcts->end(); j++) {
     ++nSize;
     const CSCDetId& id = (*j).first;
     std::vector<CSCCLCTDigi>::const_iterator digiIt = (*j).second.first;
     std::vector<CSCCLCTDigi>::const_iterator last = (*j).second.second;
     for( ; digiIt != last; ++digiIt) {
       // Valid digi in the chamber (or in neighbouring chamber) 
       if((*digiIt).isValid()){
     allCLCT[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh] = true;
       }
     }
   }
   CLCTPerEvent->Fill(nSize);
   //---- CorrLCTDigis
   for (CSCCorrelatedLCTDigiCollection::DigiRangeIterator j=correlatedlcts->begin(); j!=correlatedlcts->end(); j++) {
     const CSCDetId& id = (*j).first;
     std::vector<CSCCorrelatedLCTDigi>::const_iterator digiIt = (*j).second.first;
     std::vector<CSCCorrelatedLCTDigi>::const_iterator last = (*j).second.second;
     for( ; digiIt != last; ++digiIt) {
       // Valid digi in the chamber (or in neighbouring chamber) 
       if((*digiIt).isValid()){
     allCorrLCT[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh] = true;
       }
     }
   }
 }

void CSCEfficiency::fillRechitsSegments_info	(	edm::Handle< CSCRecHit2DCollection > &	rechits,
		edm::Handle< CSCSegmentCollection > &	segments,
		edm::ESHandle< CSCGeometry > &	cscGeom
	)

private

Definition at line 853 of file CSCEfficiency.cc.

References relativeConstraints::chamber, CSCDetId::chamber(), gather_cfg::cout, CSCDetId::endcap(), makeMuonMisalignmentScenario::endcap, FirstCh, CSCDetId::layer(), CSCGeometry::layer(), NumCh, relativeConstraints::ring, CSCDetId::ring(), findQualityFiles::size, relativeConstraints::station, CSCDetId::station(), GeomDet::toGlobal(), PV3DBase< T, PVType, FrameType >::x(), LocalError::xx(), LocalError::xy(), PV3DBase< T, PVType, FrameType >::y(), LocalError::yy(), and PV3DBase< T, PVType, FrameType >::z().

                           {
   //---- RECHITS AND SEGMENTS
   //---- Loop over rechits 
   if (printalot){ 
     //printf("\tGet the recHits collection.\t ");
     printf("  The size of the rechit collection is %i\n",int(rechits->size()));
     //printf("\t...start loop over rechits...\n");
   }
   recHitsPerEvent->Fill(rechits->size());
   //---- Build iterator for rechits and loop :
   CSCRecHit2DCollection::const_iterator recIt;
   for (recIt = rechits->begin(); recIt != rechits->end(); recIt++) {
     //---- Find chamber with rechits in CSC 
     CSCDetId id = (CSCDetId)(*recIt).cscDetId();
     if (printalot){
       const CSCLayer* csclayer = cscGeom->layer( id);
       LocalPoint rhitlocal = (*recIt).localPosition();  
       LocalError rerrlocal = (*recIt).localPositionError();  
       GlobalPoint rhitglobal= csclayer->toGlobal(rhitlocal);
       printf("\t\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber(),id.layer());
       printf("\t\tx,y,z: %f, %f, %f\texx,eey,exy: %f, %f, %f\tglobal x,y,z: %f, %f, %f \n",
          rhitlocal.x(), rhitlocal.y(), rhitlocal.z(), rerrlocal.xx(), rerrlocal.yy(), rerrlocal.xy(),
          rhitglobal.x(), rhitglobal.y(), rhitglobal.z());
     }
     std::pair <LocalPoint, bool> recHitPos((*recIt).localPosition(), false);  
     allRechits[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh][id.layer()-1].push_back(recHitPos);
   }
   //---- "Empty" chambers
   for(int iE=0;iE<2;iE++){
     for(int iS=0;iS<4;iS++){
       for(int iR=0;iR<4;iR++){
     for(int iC=0;iC<NumCh;iC++){
       int numLayers = 0;
       for(int iL=0;iL<6;iL++){
         if(allRechits[iE][iS][iR][iC][iL].size()){
           ++numLayers; 
         }
       }
       if(numLayers>1){
         emptyChambers[iE][iS][iR][iC] = false;
       }
       else{
         emptyChambers[iE][iS][iR][iC] = true;
       }
     }
       }
     }
   }
 
   //
   if (printalot){
     printf("  The size of the segment collection is %i\n", int(segments->size()));
     //printf("\t...start loop over segments...\n");
   }
   segmentsPerEvent->Fill(segments->size());
   for(CSCSegmentCollection::const_iterator it = segments->begin(); it != segments->end(); it++) {
     CSCDetId id  = (CSCDetId)(*it).cscDetId();
     StHist[id.endcap()-1][id.station()-1].segmentChi2_ndf->Fill((*it).chi2()/(*it).degreesOfFreedom());
     StHist[id.endcap()-1][id.station()-1].hitsInSegment->Fill((*it).nRecHits());
     if (printalot){ 
       printf("\tendcap/station/ring/chamber: %i %i %i %i\n",
          id.endcap(),id.station(),id.ring(),id.chamber());
       std::cout<<"\tposition(loc) = "<<(*it).localPosition()<<" error(loc) = "<<(*it).localPositionError()<<std::endl;
       std::cout<<"\t chi2/ndf = "<<(*it).chi2()/(*it).degreesOfFreedom()<<" nhits = "<<(*it).nRecHits() <<std::endl;
 
     }
     allSegments[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh].push_back
       (make_pair((*it).localPosition(), (*it).localDirection()));
 
 
     //---- try to get the CSC recHits that contribute to this segment.
     //if (printalot) printf("\tGet the recHits for this segment.\t");
     std::vector<CSCRecHit2D> theseRecHits = (*it).specificRecHits();
     int nRH = (*it).nRecHits();
     if (printalot){
       printf("\tGet the recHits for this segment.\t");
       printf("    nRH = %i\n",nRH);
     }
     //---- Find which of the rechits in the chamber is in the segment
     int layerRH = 0;
     for ( vector<CSCRecHit2D>::const_iterator iRH = theseRecHits.begin(); iRH != theseRecHits.end(); iRH++) {
       ++layerRH;
       CSCDetId idRH = (CSCDetId)(*iRH).cscDetId();
       if(printalot){ 
     printf("\t%i RH\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",
            layerRH,idRH.endcap(),idRH.station(),idRH.ring(),idRH.chamber(),idRH.layer());
       }
       for(size_t jRH = 0; 
       jRH<allRechits[idRH.endcap()-1][idRH.station()-1][idRH.ring()-1][idRH.chamber()-FirstCh][idRH.layer()-1].size();
       ++jRH){
     float xDiff = iRH->localPosition().x() - 
       allRechits[idRH.endcap()-1][idRH.station()-1][idRH.ring()-1][idRH.chamber()-FirstCh][idRH.layer()-1][jRH].first.x();
     float yDiff = iRH->localPosition().y() - 
       allRechits[idRH.endcap()-1][idRH.station()-1][idRH.ring()-1][idRH.chamber()-FirstCh][idRH.layer()-1][jRH].first.y();
     if(fabs(xDiff)<0.0001 && fabs(yDiff)<0.0001){
       std::pair <LocalPoint, bool> 
         recHitPos(allRechits[idRH.endcap()-1][idRH.station()-1][idRH.ring()-1][idRH.chamber()-FirstCh][idRH.layer()-1][jRH].first, true);
       allRechits[idRH.endcap()-1][idRH.station()-1][idRH.ring()-1][idRH.chamber()-FirstCh][idRH.layer()-1][jRH] = recHitPos;
       if(printalot){ 
         std::cout<<" number of the rechit (from zero) in the segment = "<< jRH<<std::endl;
       }
     }
       }
     }
   }
 }

void CSCEfficiency::fillSimhit_info ( edm::Handle< edm::PSimHitContainer > & simHits )

private

Definition at line 842 of file CSCEfficiency.cc.

References CSCDetId::chamber(), CSCDetId::endcap(), FirstCh, CSCDetId::layer(), CSCDetId::ring(), and CSCDetId::station().

                                                                           {
   //---- SIMHITS
   edm::PSimHitContainer::const_iterator dSHsimIter;
   for (dSHsimIter = simhits->begin(); dSHsimIter != simhits->end(); dSHsimIter++){
     // Get DetID for this simHit:
     CSCDetId sId = (CSCDetId)(*dSHsimIter).detUnitId();
     std::pair <LocalPoint, int> simHitPos((*dSHsimIter).localPosition(), (*dSHsimIter).particleType());  
     allSimhits[sId.endcap()-1][sId.station()-1][sId.ring()-1][sId.chamber()-FirstCh][sId.layer()-1].push_back(simHitPos);
   }
 }

void CSCEfficiency::fillStrips_info ( edm::Handle< CSCStripDigiCollection > & strips )

private

Definition at line 804 of file CSCEfficiency.cc.

References mps_update::diff, objects.autophobj::float, plotBeamSpotDB::last, and ctppsDiamondLocalTracks_cfi::threshold.

                                                                             {
   //---- STRIPS
   for (CSCStripDigiCollection::DigiRangeIterator j=strips->begin(); j!=strips->end(); j++) {
     CSCDetId id = (CSCDetId)(*j).first;
     int largestADCValue = -1;
     std::vector<CSCStripDigi>::const_iterator digiItr = (*j).second.first;
     std::vector<CSCStripDigi>::const_iterator last = (*j).second.second;
     for( ; digiItr != last; ++digiItr) {
       int maxADC=largestADCValue;
       int myStrip = digiItr->getStrip();
       std::vector<int> myADCVals = digiItr->getADCCounts();
       float thisPedestal = 0.5*(float)(myADCVals[0]+myADCVals[1]);
       float threshold = 13.3 ;
       float diff = 0.;
       float peakADC  = -1000.;
       for (unsigned int iCount = 0; iCount < myADCVals.size(); iCount++) {
     diff = (float)myADCVals[iCount]-thisPedestal;
     if (diff > threshold) { 
       if (myADCVals[iCount] > largestADCValue) {
         largestADCValue = myADCVals[iCount];
       }
     }
     if (diff > threshold && diff > peakADC) {
       peakADC  = diff;
     }
       }
       if(largestADCValue>maxADC){// FIX IT!!!
     maxADC = largestADCValue;
     std::pair <int, float> LayerSignal (myStrip, peakADC);
     
         //---- AllStrips contains basic information about strips 
         //---- (strip number and peak signal for most significant strip in the layer) 
     allStrips[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-1][id.layer()-1].clear();
     allStrips[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-1][id.layer()-1].push_back(LayerSignal);
       }
     }
   }
 }

void CSCEfficiency::fillWG_info	(	edm::Handle< CSCWireDigiCollection > &	wires,
		edm::ESHandle< CSCGeometry > &	cscGeom
	)

private

Definition at line 778 of file CSCEfficiency.cc.

References FirstCh, plotBeamSpotDB::last, CSCGeometry::layer(), and CSCLayerGeometry::yOfWireGroup().

                                                                                                          {
   //---- WIRE GROUPS
   for (CSCWireDigiCollection::DigiRangeIterator j=wires->begin(); j!=wires->end(); j++) {
     CSCDetId id = (CSCDetId)(*j).first;
     const CSCLayer *layer_p = cscGeom->layer (id);
     const CSCLayerGeometry *layerGeom = layer_p->geometry ();
     //
     std::vector<CSCWireDigi>::const_iterator digiItr = (*j).second.first;
     std::vector<CSCWireDigi>::const_iterator last = (*j).second.second;
     //
     for( ; digiItr != last; ++digiItr) {
       std::pair < int, float > WG_pos(digiItr->getWireGroup(), layerGeom->yOfWireGroup(digiItr->getWireGroup())); 
       std::pair <std::pair < int, float >, int >  LayerSignal(WG_pos, digiItr->getTimeBin()); 
       
       //---- AllWG contains basic information about WG (WG number and Y-position, time bin)
       allWG[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh]
     [id.layer()-1].push_back(LayerSignal);
       if(printalot){
     //std::cout<<" WG check : "<<std::endl;
     //printf("\t\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber(),id.layer());
     //std::cout<<" WG size = "<<allWG[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh]
     //[id.layer()-1].size()<<std::endl;
       }
     }
   }
 }

bool CSCEfficiency::filter	(	edm::Event &	event,
		const edm::EventSetup &	eventSetup
	)

privatevirtual

Definition at line 22 of file CSCEfficiency.cc.

                                                                            {
   passTheEvent = false;
   DataFlow->Fill(0.);  
   MuonPatternRecoDumper debug;
  
   //---- increment counter
   nEventsAnalyzed++;
   // printalot debug output
   printalot = (nEventsAnalyzed < int(printout_NEvents)); // 
   edm::RunNumber_t const iRun = event.id().run();
   edm::EventNumber_t const iEvent = event.id().event();
   if(0==fmod(double (nEventsAnalyzed) ,double(1000) )){
     if(printalot){
       printf("\n==enter==CSCEfficiency===== run %u\tevent %llu\tn Analyzed %i\n",iRun,iEvent,nEventsAnalyzed);
     }
   }
   theService->update(eventSetup);  
   //---- These declarations create handles to the types of records that you want
   //---- to retrieve from event "e".
   if (printalot) printf("\tget handles for digi collections\n");
   
   //---- Pass the handle to the method "getByType", which is used to retrieve
   //---- one and only one instance of the type in question out of event "e". If
   //---- zero or more than one instance exists in the event an exception is thrown.
   if (printalot) printf("\tpass handles\n");
   edm::Handle<CSCALCTDigiCollection> alcts;
   edm::Handle<CSCCLCTDigiCollection> clcts;
   edm::Handle<CSCCorrelatedLCTDigiCollection> correlatedlcts;
   edm::Handle<CSCWireDigiCollection> wires;
   edm::Handle<CSCStripDigiCollection> strips;
   edm::Handle<CSCRecHit2DCollection> rechits; 
   edm::Handle<CSCSegmentCollection> segments;
   edm::Handle<edm::View<reco::Track> > trackCollectionH;
   edm::Handle<edm::PSimHitContainer> simhits;
 
   if(useDigis){
     event.getByToken( wd_token, wires );
     event.getByToken( sd_token, strips );
     event.getByToken( al_token, alcts );
     event.getByToken( cl_token, clcts );
     event.getByToken( co_token, correlatedlcts );
   }
   if(!isData){
     event.getByToken( sh_token, simhits );
   }
   event.getByToken( rh_token, rechits );
   event.getByToken( se_token, segments );
   event.getByToken( tk_token, trackCollectionH );
   const edm::View<reco::Track>  trackCollection = *(trackCollectionH.product());
 
   //---- Get the CSC Geometry :
   if (printalot) printf("\tget the CSC geometry.\n");
   edm::ESHandle<CSCGeometry> cscGeom;
   eventSetup.get<MuonGeometryRecord>().get(cscGeom);
 
   // use theTrackingGeometry instead of cscGeom?
   edm::ESHandle<GlobalTrackingGeometry> theTrackingGeometry;
   eventSetup.get<GlobalTrackingGeometryRecord>().get(theTrackingGeometry);
 
   bool triggerPassed = true;
   if(useTrigger){
     // access the trigger information
     // trigger names can be find in HLTrigger/Configuration/python/HLT_2E30_cff.py (or?)
    // get hold of TriggerResults
     edm::Handle<edm::TriggerResults> hltR;
     event.getByToken( ht_token, hltR );
     const edm::TriggerNames & triggerNames = event.triggerNames(*hltR);
     triggerPassed = applyTrigger(hltR, triggerNames);
   }
   if(!triggerPassed){
     return triggerPassed;
   }
   DataFlow->Fill(1.);   
   GlobalPoint gpZero(0.,0.,0.);
   if(theService->magneticField()->inTesla(gpZero).mag2()<0.1){ 
     magField = false;
   }
   else{
     magField = true;
   }
 
   //---- store info from digis
   fillDigiInfo(alcts, clcts, correlatedlcts, wires, strips, simhits, rechits, segments, cscGeom);
   //
   edm::Handle<reco::MuonCollection> muons;
   edm::InputTag muonTag_("muons");
   event.getByLabel(muonTag_,muons);
 
   edm::Handle<reco::BeamSpot> beamSpotHandle;
   event.getByLabel("offlineBeamSpot", beamSpotHandle);
   reco::BeamSpot vertexBeamSpot = *beamSpotHandle;
   //
   std::vector <reco::MuonCollection::const_iterator> goodMuons_it;
   unsigned int nPositiveZ = 0;
   unsigned int nNegativeZ = 0;
   float muonOuterZPosition = -99999.;
   if(isIPdata){
     if (printalot)std::cout<<" muons.size() = "<<muons->size() <<std::endl;
     for ( reco::MuonCollection::const_iterator muon = muons->begin(); muon != muons->end(); ++muon ) {
       DataFlow->Fill(31.);  
       if (printalot) {
         std::cout<<"  iMuon = "<<muon-muons->begin()<<" charge = "<<muon->charge()<<" p = "<<muon->p()<<" pt = "<<muon->pt()<<
           " eta = "<<muon->eta()<<" phi = "<<muon->phi()<<
           " matches = "<<
           muon->matches().size()<<" matched Seg = "<<muon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration)<<" GLB/TR/STA = "<<
           muon->isGlobalMuon()<<"/"<<muon->isTrackerMuon()<<"/"<<muon->isStandAloneMuon()<<std::endl;
       }
       if(!(muon->isTrackerMuon() && muon->isGlobalMuon())){
     continue;
       }
       DataFlow->Fill(32.);
       double relISO = ( muon->isolationR03().sumPt +
                         muon->isolationR03().emEt +
                         muon->isolationR03().hadEt)/muon->track()->pt();
       if (printalot) {
         std::cout<<" relISO = "<<relISO<<" emVetoEt = "<<muon->isolationR03().emVetoEt<<" caloComp = "<<
           muon::caloCompatibility(*(muon))<<" dxy = "<<fabs(muon->track()->dxy(vertexBeamSpot.position()))<<std::endl;
       }
       if(
      //relISO>0.1 || muon::caloCompatibility(*(muon))<.90 ||
          fabs(muon->track()->dxy(vertexBeamSpot.position()))>0.2 || muon->pt()<6.){
         continue;
       }
       DataFlow->Fill(33.);
       if(muon->track()->hitPattern().numberOfValidPixelHits()<1 ||
      muon->track()->hitPattern().numberOfValidTrackerHits()<11 ||
      muon->combinedMuon()->hitPattern().numberOfValidMuonHits()<1 ||
      muon->combinedMuon()->normalizedChi2()>10. ||
      muon->numberOfMatches()<2){
     continue;
       }
       DataFlow->Fill(34.);
       float zOuter = muon->combinedMuon()->outerPosition().z();
       float rhoOuter = muon->combinedMuon()->outerPosition().rho();
       bool passDepth = true;
       // barrel region 
       //if ( fabs(zOuter) < 660. && rhoOuter > 400. && rhoOuter < 480.){
       if ( fabs(zOuter) < 660. && rhoOuter > 400. && rhoOuter < 540.){ 
     passDepth = false;
       }
       // endcap region
       //else if( fabs(zOuter) > 550. && fabs(zOuter) < 650. && rhoOuter < 300.){
       else if( fabs(zOuter) > 550. && fabs(zOuter) < 650. && rhoOuter < 300.){
     passDepth = false;
       }
       // overlap region
       //else if ( fabs(zOuter) > 680. && fabs(zOuter) < 730. && rhoOuter < 480.){
       else if ( fabs(zOuter) > 680. && fabs(zOuter) < 880. && rhoOuter < 540.){
     passDepth = false;
       }
       if(!passDepth){
     continue;
       }
       DataFlow->Fill(35.);
       goodMuons_it.push_back(muon);
       if(muon->track()->momentum().z()>0.){
     ++nPositiveZ;
       }
       if(muon->track()->momentum().z()<0.){
     ++nNegativeZ;
       }
     }
   }
 
   //
 
 
   if (printalot) std::cout<<"Start track loop over "<<trackCollection.size()<<" tracks"<<std::endl;
   for(edm::View<reco::Track>::size_type i=0; i<trackCollection.size(); ++i) {
     DataFlow->Fill(2.);
     edm::RefToBase<reco::Track> track(trackCollectionH, i);
     //std::cout<<" iTR = "<<i<<" eta = "<<track->eta()<<" phi = "<<track->phi()<<std::cout<<" pt = "<<track->pt()<<std::endl;
     if(isIPdata){
       if (printalot){
     std::cout<<" nNegativeZ = "<<nNegativeZ<<" nPositiveZ = "<<nPositiveZ<<std::endl;
       }
       if(nNegativeZ>1 || nPositiveZ>1){
         break;
       }
       bool trackOK = false;
       if (printalot){
     std::cout<<" goodMuons_it.size() = "<<goodMuons_it.size()<<std::endl;
       }
       for(size_t iM=0;iM<goodMuons_it.size();++iM){
     //std::cout<<" iM = "<<iM<<" eta = "<<goodMuons_it[iM]->track()->eta()<<
     //" phi = "<<goodMuons_it[iM]->track()->phi()<<
     //" pt = "<<goodMuons_it[iM]->track()->pt()<<std::endl;
         float deltaR = pow(track->phi()-goodMuons_it[iM]->track()->phi(),2) +
           pow(track->eta()-goodMuons_it[iM]->track()->eta(),2);
         deltaR = sqrt(deltaR);
         if (printalot){
           std::cout<<" TR mu match to a tr: deltaR = "<<deltaR<<" dPt = "<<
         track->pt()-goodMuons_it[iM]->track()->pt()<<std::endl;
         }
         if(deltaR>0.01 || fabs(track->pt()-goodMuons_it[iM]->track()->pt())>0.1 ){
           continue;
         }
         else{
           trackOK = true;
       if (printalot){
         std::cout<<" trackOK "<<std::endl;
       }
       muonOuterZPosition = goodMuons_it[iM]->combinedMuon()->outerPosition().z();
           break;
           //++nChosenTracks;
         }
       }
       if(!trackOK){
     if (printalot){
       std::cout<<" failed: trackOK "<<std::endl;
     }
         continue;
       }
     }
     else{
       //---- Do we need a better "clean track" definition?
       if(trackCollection.size()>2){
     break;
       }
       DataFlow->Fill(3.);
       if(!i && 2==trackCollection.size()){
     edm::View<reco::Track>::size_type tType = 1;
     edm::RefToBase<reco::Track> trackTwo(trackCollectionH, tType);
     if(track->outerPosition().z()*trackTwo->outerPosition().z()>0){// in one and the same "endcap"
       break;
     }
       }
     }
     DataFlow->Fill(4.);
     if (printalot){
       std::cout<<"i track = "<<i<<" P = "<<track->p()<<" chi2/ndf = "<<track->normalizedChi2()<<" nSeg = "<<segments->size()<<std::endl;
       std::cout<<"quality undef/loose/tight/high/confirmed/goodIt/size "<<
     track->quality(reco::Track::undefQuality)<<"/"<<
     track->quality(reco::Track::loose)<<"/"<<
     track->quality(reco::Track::tight)<<"/"<<
     track->quality(reco::Track::highPurity)<<"/"<<
     track->quality(reco::Track::confirmed)<<"/"<<
     track->quality(reco::Track::goodIterative)<<"/"<<
     track->quality(reco::Track::qualitySize)<<
     std::endl;
       std::cout<<" pt = "<< track->pt()<<" +-"<<track->ptError()<<" q/pt = "<<track->qoverp()<<" +- "<<track->qoverpError()<<std::endl;
       //std::cout<<" const Pmin = "<<minTrackMomentum<<" pMax = "<<maxTrackMomentum<<" maxNormChi2 = "<<maxNormChi2<<std::endl;
       std::cout<<" track inner position = "<<track->innerPosition()<<" outer position = "<<track->outerPosition()<<std::endl;
       std::cout<<"track eta (outer) = "<<track->outerPosition().eta()<<" phi (outer) = "<<
     track->outerPosition().phi()<<std::endl;
       if(fabs(track->innerPosition().z())>500.){      
     DetId innerDetId(track->innerDetId());
     std::cout<<" dump inner state MUON detid  = "<<debug.dumpMuonId(innerDetId)<<std::endl;
       }
       if(fabs(track->outerPosition().z())>500.){
         DetId outerDetId(track->outerDetId());
     std::cout<<" dump outer state MUON detid  = "<<debug.dumpMuonId(outerDetId)<<std::endl;
       }
 
       std::cout<<" nHits = "<<track->found()<<std::endl;
       /*
       trackingRecHit_iterator rhbegin = track->recHitsBegin();
       trackingRecHit_iterator rhend = track->recHitsEnd();
       int iRH = 0;
       for(trackingRecHit_iterator recHit = rhbegin; recHit != rhend; ++recHit){
         const GeomDet* geomDet = theTrackingGeometry->idToDet((*recHit)->geographicalId());
     std::cout<<"hit "<<iRH<<" loc pos = " <<(*recHit)->localPosition()<<
       " glob pos = " <<geomDet->toGlobal((*recHit)->localPosition())<<std::endl;
         ++iRH;
       }
       */
     }
     float dpT_ov_pT = 0.;
     if(fabs(track->pt())>0.001){
       dpT_ov_pT =  track->ptError()/ track->pt();
     }
     //---- These define a "good" track
     if(track->normalizedChi2()>maxNormChi2){// quality
       break;
     }
     DataFlow->Fill(5.);
     if(track->found()<minTrackHits){// enough data points
       break;
     }
     DataFlow->Fill(6.);
     if(!segments->size()){// better have something in the CSC 
       break;
     }
     DataFlow->Fill(7.);
     if(magField && (track->p()<minP || track->p()>maxP)){// proper energy range 
       break;
     }
     DataFlow->Fill(8.);
     if(magField && (dpT_ov_pT >0.5) ){// not too crazy uncertainty
       break;
     }
     DataFlow->Fill(9.);
 
     passTheEvent = true; 
     if (printalot) std::cout<<"good Track"<<std::endl;
     CLHEP::Hep3Vector r3T_inner(track->innerPosition().x(),track->innerPosition().y(),track->innerPosition().z());
     CLHEP::Hep3Vector r3T(track->outerPosition().x(),track->outerPosition().y(),track->outerPosition().z());
     chooseDirection(r3T_inner, r3T);// for non-IP
 
     CLHEP::Hep3Vector p3T(track->outerMomentum().x(),track->outerMomentum().y(),track->outerMomentum().z());
     CLHEP::Hep3Vector p3_propagated, r3_propagated;
     AlgebraicSymMatrix66 cov_propagated, covT;
     covT *= 1e-20;
     cov_propagated *= 1e-20;
     int charge = track->charge(); 
     FreeTrajectoryState ftsStart = getFromCLHEP(p3T, r3T, charge, covT, &*(theService->magneticField()));
     if (printalot){
       std::cout<<" p = "<<track->p()<<" norm chi2 = "<<track->normalizedChi2()<<std::endl;
       std::cout<<" dump the very first FTS  = "<<debug.dumpFTS(ftsStart)<<std::endl;
     }
     TrajectoryStateOnSurface tSOSDest;
     int endcap = 0;
     //---- which endcap to look at
     if(track->outerPosition().z()>0){
       endcap = 1;
     }
     else{
       endcap = 2;
     }
     int chamber = 1;
     //---- a "refference" CSCDetId for each ring
     std::vector< CSCDetId > refME;
     for(int iS=1;iS<5;++iS){
       for(int iR=1;iR<4;++iR){
     if(1!=iS && iR>2){
       continue;
     }
     else if(4==iS && iR>1){
       continue;
     } 
     refME.push_back( CSCDetId(endcap, iS, iR, chamber));
       }
     }
     //---- loop over the "refference" CSCDetIds
     for(size_t iSt = 0; iSt<refME.size();++iSt){
       if (printalot){
     std::cout<<"loop iStatation = "<<iSt<<std::endl;
     std::cout<<"refME[iSt]: st = "<<refME[iSt].station()<<" rg = "<<refME[iSt].ring()<<std::endl;
       }
       std::map <std::string, bool> chamberTypes;
       chamberTypes["ME11"] = false;
       chamberTypes["ME12"] = false;
       chamberTypes["ME13"] = false;
       chamberTypes["ME21"] = false;
       chamberTypes["ME22"] = false;
       chamberTypes["ME31"] = false;
       chamberTypes["ME32"] = false;
       chamberTypes["ME41"] = false;
       const CSCChamber* cscChamber_base = cscGeom->chamber(refME[iSt].chamberId());
       DetId detId = cscChamber_base->geographicalId();
       if (printalot){
     std::cout<<" base iStation : eta = "<<cscGeom->idToDet(detId)->surface().position().eta()<<" phi = "<<
       cscGeom->idToDet(detId)->surface().position().phi() << " y = " <<cscGeom->idToDet(detId)->surface().position().y()<<std::endl;
     std::cout<<" dump base iStation detid  = "<<debug.dumpMuonId(detId)<<std::endl;
     std::cout<<" dump FTS start  = "<<debug.dumpFTS(ftsStart)<<std::endl;
       }
       //---- propagate to this ME
       tSOSDest = propagate(ftsStart, cscGeom->idToDet(detId)->surface());
       if(tSOSDest.isValid()){
     ftsStart = *tSOSDest.freeState();
     if (printalot) std::cout<<"  dump FTS end   = "<<debug.dumpFTS(ftsStart)<<std::endl;
     getFromFTS(ftsStart, p3_propagated, r3_propagated, charge, cov_propagated);
     float feta =  fabs(r3_propagated.eta());
     float phi =  r3_propagated.phi();
     //---- which rings are (possibly) penetrated
     ringCandidates(refME[iSt].station(), feta, chamberTypes);
 
     map<std::string,bool>::iterator iter;   
     int iterations = 0;
     //---- loop over ring candidates 
     for( iter = chamberTypes.begin(); iter != chamberTypes.end(); iter++ ) {
       ++iterations;
       //---- is this ME a machinig candidate station 
       if(iter->second && (iterations-1)==int(iSt)){
         if (printalot){
           std::cout<<" Chamber type "<< iter->first<<" is a candidate..."<<std::endl;
           std::cout<<" station() = "<< refME[iSt].station()<<" ring() = "<<refME[iSt].ring()<<" iSt = "<<iSt<<std::endl;
         }
         std::vector <int> coupleOfChambers;
         //---- which chamber (and its closes neighbor) is penetrated by the track - candidates
         chamberCandidates(refME[iSt].station(), refME[iSt].ring(), phi, coupleOfChambers);
         //---- loop over the two chamber candidates
         for(size_t iCh =0;iCh<coupleOfChambers.size();++iCh){
           DataFlow->Fill(11.);  
           if (printalot) std::cout<<" Check chamber N = "<<coupleOfChambers.at(iCh)<<std::endl;;
           if((!getAbsoluteEfficiency) && (true == emptyChambers
                           [refME[iSt].endcap()-1]
                           [refME[iSt].station()-1]
                           [refME[iSt].ring()-1]
                           [coupleOfChambers.at(iCh)-FirstCh])){
         continue;
           }
           CSCDetId theCSCId(refME[iSt].endcap(), refME[iSt].station(), refME[iSt].ring(), coupleOfChambers.at(iCh));
           const CSCChamber* cscChamber = cscGeom->chamber(theCSCId.chamberId());
           const BoundPlane bpCh = cscGeom->idToDet(cscChamber->geographicalId())->surface();
           float zFTS = ftsStart.position().z();
           float dz = fabs(bpCh.position().z() - zFTS);
           float zDistInner = track->innerPosition().z() - bpCh.position().z();
           float zDistOuter = track->outerPosition().z() - bpCh.position().z();
           //---- only detectors between the inner and outer points of the track are considered for non IP-data
           if(printalot){
         std::cout<<" zIn = "<<track->innerPosition().z()<<" zOut = "<<track->outerPosition().z()<<" zSurf = "<<bpCh.position().z()<<std::endl;
           }
           if(!isIPdata && (zDistInner*zDistOuter>0. || fabs(zDistInner)<15. || fabs(zDistOuter)<15.)){ // for non IP-data
         if(printalot){
           std::cout<<" Not an intermediate (as defined) point... Skip."<<std::endl;
         }
         continue;
           }
           if(isIPdata && fabs(track->eta())<1.8){
         if(fabs(muonOuterZPosition) - fabs(bpCh.position().z())<0 || 
            fabs(muonOuterZPosition-bpCh.position().z())<15.){
           continue;
         }
           }
               DataFlow->Fill(13.);
           //---- propagate to the chamber (from this ME) if it is a different surface (odd/even chambers)
           if(dz>0.1){// i.e. non-zero (float 0 check is bad) 
         //if(fabs(zChanmber - zFTS ) > 0.1){
         tSOSDest = propagate(ftsStart, cscGeom->idToDet(cscChamber->geographicalId())->surface());           
         if(tSOSDest.isValid()){
           ftsStart = *tSOSDest.freeState();
         }
         else{
           if(printalot) std::cout<<"TSOS not valid! Break."<<std::endl;
           break;
         }
           }
               else{
                 if(printalot) std::cout<<" info: dz<0.1"<<std::endl;
               }
           DataFlow->Fill(15.);  
           FreeTrajectoryState ftsInit = ftsStart;
           bool inDeadZone = false;
           //---- loop over the 6 layers
           for(int iLayer = 0;iLayer<6;++iLayer){
         bool extrapolationPassed = true;
         if (printalot){
           std::cout<<" iLayer = "<<iLayer<<"   dump FTS init  = "<<debug.dumpFTS(ftsInit)<<std::endl;
           std::cout<<" dump detid  = "<<debug.dumpMuonId(cscChamber->geographicalId())<<std::endl;
           std::cout<<"Surface to propagate to:  pos = "<<cscChamber->layer(iLayer+1)->surface().position()<<" eta = "
                <<cscChamber->layer(iLayer+1)->surface().position().eta()<<" phi = "
                <<cscChamber->layer(iLayer+1)->surface().position().phi()<<std::endl;
         }
         //---- propagate to this layer
         tSOSDest = propagate(ftsInit, cscChamber->layer(iLayer+1)->surface());
         if(tSOSDest.isValid()){
           ftsInit = *tSOSDest.freeState();
           if (printalot) std::cout<<" Propagation between layers successful:  dump FTS end  = "<<debug.dumpFTS(ftsInit)<<std::endl;
           getFromFTS(ftsInit, p3_propagated, r3_propagated, charge, cov_propagated);
         }
         else{
           if (printalot) std::cout<<"Propagation between layers not successful - notValid TSOS"<<std::endl;
           extrapolationPassed = false;
           inDeadZone = true;
         }
         //}
         //---- Extrapolation passed? For each layer?
         if(extrapolationPassed){
           GlobalPoint theExtrapolationPoint(r3_propagated.x(),r3_propagated.y(),r3_propagated.z());
           LocalPoint theLocalPoint = cscChamber->layer(iLayer+1)->toLocal(theExtrapolationPoint);
           //std::cout<<" Candidate chamber: extrapolated LocalPoint = "<<theLocalPoint<<std::endl;
           inDeadZone = ( inDeadZone ||
                  !inSensitiveLocalRegion(theLocalPoint.x(), theLocalPoint.y(), 
                              refME[iSt].station(), refME[iSt].ring()));
           if (printalot){
                     std::cout<<" Candidate chamber: extrapolated LocalPoint = "<<theLocalPoint<<"inDeadZone = "<<inDeadZone<<std::endl;
                   }
           //---- break if in dead zone for any layer ("clean" tracks)
           if(inDeadZone){
             break;
           }
         }
         else{
           break;
         }
           }
           DataFlow->Fill(17.);
           //---- Is a track in a sensitive area for each layer?  
           if(!inDeadZone){//---- for any layer
         DataFlow->Fill(19.);  
         if (printalot) std::cout<<"Do efficiencies..."<<std::endl;
         //---- Do efficiencies
         // angle cuts applied (if configured)
         bool angle_flag = true; angle_flag = efficienciesPerChamber(theCSCId, cscChamber, ftsStart);  
         if(useDigis && angle_flag){
           stripWire_Efficiencies(theCSCId, ftsStart);
         }
                 if(angle_flag){
           recHitSegment_Efficiencies(theCSCId, cscChamber, ftsStart);
           if(!isData){
             recSimHitEfficiency(theCSCId, ftsStart);
           }
         }
           }
               else{
                 if(printalot) std::cout<<" Not in active area for all layers"<<std::endl;
           }
         }
         if(tSOSDest.isValid()){
           ftsStart = *tSOSDest.freeState();
         }
       }
     }
       }
       else{
         if (printalot) std::cout<<" TSOS not valid..."<<std::endl;
       }
     }
   }
   //---- End
   if (printalot) printf("==exit===CSCEfficiency===== run %u\tevent %llu\n\n",iRun,iEvent);
   return passTheEvent;
 }

FreeTrajectoryState CSCEfficiency::getFromCLHEP	(	const CLHEP::Hep3Vector &	p3,
		const CLHEP::Hep3Vector &	r3,
		int	charge,
		const AlgebraicSymMatrix66 &	cov,
		const MagneticField *	field
	)

private

Definition at line 1451 of file CSCEfficiency.cc.

                                                                                   {
 
   GlobalVector p3GV(p3.x(), p3.y(), p3.z());
   GlobalPoint r3GP(r3.x(), r3.y(), r3.z());
   GlobalTrajectoryParameters tPars(r3GP, p3GV, charge, field);
 
   CartesianTrajectoryError tCov(cov);
 
   return cov.kRows == 6 ? FreeTrajectoryState(tPars, tCov) : FreeTrajectoryState(tPars) ;
 }

void CSCEfficiency::getFromFTS	(	const FreeTrajectoryState &	fts,
		CLHEP::Hep3Vector &	p3,
		CLHEP::Hep3Vector &	r3,
		int &	charge,
		AlgebraicSymMatrix66 &	cov
	)

private

Definition at line 1436 of file CSCEfficiency.cc.

References FreeTrajectoryState::cartesianError(), FreeTrajectoryState::charge(), FreeTrajectoryState::hasError(), CartesianTrajectoryError::matrix(), FreeTrajectoryState::momentum(), FreeTrajectoryState::position(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

                                                        {
   
   GlobalVector p3GV = fts.momentum();
   GlobalPoint r3GP = fts.position();
 
   p3.set(p3GV.x(), p3GV.y(), p3GV.z());
   r3.set(r3GP.x(), r3GP.y(), r3GP.z());
 
   charge = fts.charge();
   cov = fts.hasError() ? fts.cartesianError().matrix() : AlgebraicSymMatrix66();
 
 }

bool CSCEfficiency::inSensitiveLocalRegion	(	double	xLocal,
		double	yLocal,
		int	station,
		int	ring
	)

private

Definition at line 538 of file CSCEfficiency.cc.

References funct::abs().

                                                                                              {
   //---- Good region means sensitive area of a chamber. "Local" stands for the local system 
   bool pass = false;
   std::vector <double> chamberBounds(3);// the sensitive area
   float y_center = 99999.;
   //---- hardcoded... not good
   if(station>1 && station<5){
     if(2==ring){
       chamberBounds[0] = 66.46/2; // (+-)x1 shorter
       chamberBounds[1] = 127.15/2; // (+-)x2 longer 
       chamberBounds[2] = 323.06/2;
       y_center = -0.95;
     }
     else{
       if(2==station){
     chamberBounds[0] = 54.00/2; // (+-)x1 shorter
     chamberBounds[1] = 125.71/2; // (+-)x2 longer 
     chamberBounds[2] = 189.66/2;
     y_center = -0.955;
       }
       else if(3==station){
     chamberBounds[0] = 61.40/2; // (+-)x1 shorter
     chamberBounds[1] = 125.71/2; // (+-)x2 longer 
     chamberBounds[2] = 169.70/2;
     y_center = -0.97;
       }
       else if(4==station){
     chamberBounds[0] = 69.01/2; // (+-)x1 shorter
     chamberBounds[1] = 125.65/2; // (+-)x2 longer 
     chamberBounds[2] = 149.42/2;
     y_center = -0.94;
       }
     }
   }
   else if(1==station){
     if(3==ring){
       chamberBounds[0] = 63.40/2; // (+-)x1 shorter
       chamberBounds[1] = 92.10/2; // (+-)x2 longer 
       chamberBounds[2] = 164.16/2;
       y_center = -1.075;
     }
     else if(2==ring){
       chamberBounds[0] = 51.00/2; // (+-)x1 shorter
       chamberBounds[1] = 83.74/2; // (+-)x2 longer 
       chamberBounds[2] = 174.49/2;
       y_center = -0.96;
     }
     else{// to be investigated
       chamberBounds[0] = 30./2;//40./2; // (+-)x1 shorter
       chamberBounds[1] = 60./2;//100./2; // (+-)x2 longer 
       chamberBounds[2] = 160./2;//142./2;
       y_center = 0.;
     }
   }
   double yUp = chamberBounds[2] + y_center;
   double yDown = - chamberBounds[2] + y_center;
   double xBound1Shifted = chamberBounds[0]-distanceFromDeadZone;//
   double xBound2Shifted = chamberBounds[1]-distanceFromDeadZone;//
   double lineSlope = (yUp - yDown)/(xBound2Shifted-xBound1Shifted);
   double lineConst = yUp - lineSlope*xBound2Shifted;
   double yBoundary =  lineSlope*abs(xLocal) + lineConst;
   pass = checkLocal(yLocal, yBoundary, station, ring);
   return pass;
 }

void CSCEfficiency::linearExtrapolation	(	GlobalPoint	initialPosition,
		GlobalVector	initialDirection,
		float	zSurface,
		std::vector< float > &	posZY
	)

private

Definition at line 1464 of file CSCEfficiency.cc.

References PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), ctppsDiamondLocalTracks_cfi::yPosition, and PV3DBase< T, PVType, FrameType >::z().

                                                              {
   double paramLine =  lineParameter(initialPosition.z(), zSurface, initialDirection.z());
   double xPosition = extrapolate1D(initialPosition.x(), initialDirection.x(),paramLine);
   double yPosition = extrapolate1D(initialPosition.y(), initialDirection.y(),paramLine);
   posZY.clear();
   posZY.push_back(xPosition);
   posZY.push_back(yPosition);
 }

double CSCEfficiency::lineParameter	(	double	initZPosition,
		double	destZPosition,
		double	initZDirection
	)

private

Definition at line 1479 of file CSCEfficiency.cc.

                                                                                                     {
   double paramLine = (destZPosition-initZPosition)/initZDirection;
   return paramLine;
 }

TrajectoryStateOnSurface CSCEfficiency::propagate	(	FreeTrajectoryState &	ftsStart,
		const BoundPlane &	bp
	)

private

Definition at line 1514 of file CSCEfficiency.cc.

References PhotonConversionTrajectorySeedProducerFromQuadruplets_cfi::propagator, TSGForRoadSearch_cfi::propagatorName, and AlCaHLTBitMon_QueryRunRegistry::string.

                                                                                                          {
   TrajectoryStateOnSurface tSOSDest;
   std::string propagatorName;
 /*
 // it would work if cosmic muons had properly assigned direction...
   bool dzPositive = bpDest.position().z() - ftsStart.position().z() > 0 ? true : false;
  //---- Be careful with trigger conditions too
   if(!isIPdata){
     bool rightDirection = !(alongZ^dzPositive);
     if(rightDirection){
       if(printalot) std::cout<<" propagate along momentum"<<std::endl;
       propagatorName = "SteppingHelixPropagatorAlong";
     }
     else{
       if(printalot) std::cout<<" propagate opposite momentum"<<std::endl;
       propagatorName = "SteppingHelixPropagatorOpposite";
     }
   }
   else{
     if(printalot) std::cout<<" propagate any (momentum)"<<std::endl;
     propagatorName = "SteppingHelixPropagatorAny";
   }
 */
   propagatorName = "SteppingHelixPropagatorAny";
   tSOSDest = propagator(propagatorName)->propagate(ftsStart, bpDest);   
   return tSOSDest;
 }

const Propagator * CSCEfficiency::propagator ( std::string propagatorName ) const

private

Definition at line 1509 of file CSCEfficiency.cc.

                                                                           {
   return &*theService->propagator(propagatorName);
 }

bool CSCEfficiency::recHitSegment_Efficiencies	(	CSCDetId &	cscDetId,
		const CSCChamber *	cscChamber,
		FreeTrajectoryState &	ftsChamber
	)

private

Definition at line 1295 of file CSCEfficiency.cc.

References relativeConstraints::chamber, gather_cfg::cout, RecoEcal_EventContent_cff::ec, makeMuonMisalignmentScenario::endcap, PV3DBase< T, PVType, FrameType >::eta(), plotBeamSpotDB::first, FirstCh, CSCLayer::geometry(), CSCChamber::layer(), M_PI, CSCLayerGeometry::nearestStrip(), FreeTrajectoryState::position(), funct::pow(), relativeConstraints::ring, edm::second(), findQualityFiles::size, mathSSE::sqrt(), relativeConstraints::station, CSCLayerGeometry::stripAngle(), GeomDet::surface(), GeomDet::toGlobal(), GeomDet::toLocal(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

                                                                                                                           {
   int ec, st, rg, ch, secondRing;
   returnTypes(id, ec, st, rg, ch, secondRing);
   bool firstCondition, secondCondition;
 
   std::vector <bool> missingLayers_rh(6);
   std::vector <int> usedInSegment(6);
   // Rechits
   if(printalot) std::cout<<"RecHits eff"<<std::endl;
   for(int iLayer=0;iLayer<6;++iLayer){
     firstCondition = allRechits[ec][st][rg][ch][iLayer].size() ? true : false;
     secondCondition = false;
     int thisRing = rg;
     if(secondRing>-1){
       secondCondition = allRechits[ec][st][secondRing][ch][iLayer].size() ? true : false;
       if(secondCondition){
     thisRing = secondRing;
       }
     }
     if(firstCondition || secondCondition){
       ChHist[ec][st][rg][ch].EfficientRechits_good->Fill(iLayer+1);
       for(size_t iR=0;
       iR<allRechits[ec][st][thisRing][ch][iLayer].size();
       ++iR){
     if(allRechits[ec][st][thisRing][ch][iLayer][iR].second){
       usedInSegment[iLayer] = 1;
       break;
     }
     else{
       usedInSegment[iLayer] = -1;
     }
       }
     }
     else{
       missingLayers_rh[iLayer] = true;
       if(printalot){
     std::cout<<"missing rechits ";
     printf("\t\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber(),iLayer+1);
       }
     }
   }
   GlobalVector globalDir;
   GlobalPoint globalPos;
  // Segments
   firstCondition = allSegments[ec][st][rg][ch].size() ? true : false;
   secondCondition = false;
   int secondSize = 0;
   int thisRing = rg;
   if(secondRing>-1){
     secondCondition = allSegments[ec][st][secondRing][ch].size() ? true : false;
     secondSize = allSegments[ec][st][secondRing][ch].size();
     if(secondCondition){
       thisRing = secondRing;
     }
   }
   if(firstCondition || secondCondition){
     if (printalot) std::cout<<"segments - start ec = "<<ec<<" st = "<<st<<" rg = "<<rg<<" ch = "<<ch<<std::endl;
     StHist[ec][st].EfficientSegments_XY->Fill(ftsChamber.position().x(),ftsChamber.position().y());
     if(1==allSegments[ec][st][rg][ch].size() + secondSize){
       globalDir = cscChamber->toGlobal(allSegments[ec][st][thisRing][ch][0].second);
       globalPos = cscChamber->toGlobal(allSegments[ec][st][thisRing][ch][0].first);
       StHist[ec][st].EfficientSegments_eta->Fill(fabs(ftsChamber.position().eta()));
       double residual = sqrt(pow(ftsChamber.position().x() - globalPos.x(),2)+
                  pow(ftsChamber.position().y() - globalPos.y(),2)+
                  pow(ftsChamber.position().z() - globalPos.z(),2));
       if (printalot) std::cout<<" fts.position() = "<<ftsChamber.position()<<" segPos = "<<globalPos<<" res = "<<residual<< std::endl;
       StHist[ec][st].ResidualSegments->Fill(residual);
     }
     for(int iLayer=0;iLayer<6;++iLayer){
       if(printalot) std::cout<<" iLayer = "<<iLayer<<" usedInSegment = "<<usedInSegment[iLayer]<<std::endl;
       if(0!=usedInSegment[iLayer]){
     if(-1==usedInSegment[iLayer]){
       ChHist[ec][st][rg][ch].InefficientSingleHits->Fill(iLayer+1);
     }
     ChHist[ec][st][rg][ch].AllSingleHits->Fill(iLayer+1);
       }
       firstCondition = allRechits[ec][st][rg][ch][iLayer].size() ? true : false;
       secondCondition = false;
       if(secondRing>-1){
     secondCondition = allRechits[ec][st][secondRing][ch][iLayer].size() ? true : false;
       }
       float stripAngle = 99999.;
       std::vector<float> posXY(2);
       bool oneSegment = false;
       if(1==allSegments[ec][st][rg][ch].size() + secondSize){
     oneSegment = true;
     const BoundPlane bp = cscChamber->layer(iLayer+1)->surface();
     linearExtrapolation(globalPos,globalDir, bp.position().z(), posXY);
         GlobalPoint gp_extrapol( posXY.at(0), posXY.at(1),bp.position().z());
         const LocalPoint lp_extrapol = cscChamber->layer(iLayer+1)->toLocal(gp_extrapol);
         posXY.at(0) = lp_extrapol.x();
         posXY.at(1) = lp_extrapol.y();
         int nearestStrip = cscChamber->layer(iLayer+1)->geometry()->nearestStrip(lp_extrapol);
         stripAngle = cscChamber->layer(iLayer+1)->geometry()->stripAngle(nearestStrip) - M_PI/2. ;
       }
       if(firstCondition || secondCondition){
     ChHist[ec][st][rg][ch].EfficientRechits_inSegment->Fill(iLayer+1);
     if(oneSegment){
       ChHist[ec][st][rg][ch].Y_EfficientRecHits_inSegment[iLayer]->Fill(posXY.at(1));
       ChHist[ec][st][rg][ch].Phi_EfficientRecHits_inSegment[iLayer]->Fill(stripAngle);
     }
       }
       else{
     if(oneSegment){
       ChHist[ec][st][rg][ch].Y_InefficientRecHits_inSegment[iLayer]->Fill(posXY.at(1));
       ChHist[ec][st][rg][ch].Phi_InefficientRecHits_inSegment[iLayer]->Fill(stripAngle);
     }
       }
     }
   }
   else{
     StHist[ec][st].InefficientSegments_XY->Fill(ftsChamber.position().x(),ftsChamber.position().y());
     if(printalot){
       std::cout<<"missing segment "<<std::endl;
       printf("\t\tendcap/station/ring/chamber: %i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber());
       std::cout<<" fts.position() = "<<ftsChamber.position()<<std::endl;
     }
   }
   // Normalization
   ChHist[ec][st][rg][ch].EfficientRechits_good->Fill(8);
   if(allSegments[ec][st][rg][ch].size()+secondSize<2){
     StHist[ec][st].AllSegments_eta->Fill(fabs(ftsChamber.position().eta()));
   }
   ChHist[ec][st][rg][id.chamber()-FirstCh].EfficientRechits_inSegment->Fill(9);
 
   return true;
 }

bool CSCEfficiency::recSimHitEfficiency	(	CSCDetId &	id,
		FreeTrajectoryState &	ftsChamber
	)

private

Definition at line 1248 of file CSCEfficiency.cc.

References RecoEcal_EventContent_cff::ec, edm::second(), and findQualityFiles::size.

                                                                                      {
   int ec, st, rg, ch, secondRing;
   returnTypes(id, ec, st, rg, ch, secondRing);
   bool firstCondition, secondCondition;
   for(int iLayer=0; iLayer<6;iLayer++){
     firstCondition = allSimhits[ec][st][rg][ch][iLayer].size() ? true : false;
     secondCondition = false;
     int thisRing = rg;
     if(secondRing>-1){
       secondCondition = allSimhits[ec][st][secondRing][ch][iLayer].size() ? true : false;
       if(secondCondition){
     thisRing = secondRing;
       }
     }
     if(firstCondition || secondCondition){
       for(size_t iSH=0;
       iSH<allSimhits[ec][st][thisRing][ch][iLayer].size();
       ++iSH){
     if(13 ==
        fabs(allSimhits[ec][st][thisRing][ch][iLayer][iSH].second)){
       ChHist[ec][st][rg][ch].SimSimhits->Fill(iLayer+1);
       if(allRechits[ec][st][thisRing][ch][iLayer].size()){
         ChHist[ec][st][rg][ch].SimRechits->Fill(iLayer+1);
       }
       break;
     }
       }
       //---- Next is not too usefull... 
       /*
       for(unsigned int iSimHits=0;
       iSimHits<allSimhits[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh][iLayer].size();
       iSimHits++){
     ChHist[ec][st][rg][id.chamber()-FirstCh].SimSimhits_each->Fill(iLayer+1);
       }
       for(unsigned int iRecHits=0;
       iRecHits<allRechits[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh][iLayer].size();
       iRecHits++){
     ChHist[ec][st][rg][id.chamber()-FirstCh].SimRechits_each->Fill(iLayer+1);
       }
       */
       //
     }   
   }
   return true;
 }

void CSCEfficiency::returnTypes	(	CSCDetId &	id,
		int &	ec,
		int &	st,
		int &	rg,
		int &	ch,
		int &	secondRing
	)

private

Definition at line 1423 of file CSCEfficiency.cc.

References FirstCh, relativeConstraints::ring, and relativeConstraints::station.

                                                                                                  {
   ec = id.endcap()-1;
   st = id.station()-1;
   rg = id.ring()-1;
   secondRing = -1;
   if(1==id.station() && (4==id.ring() || 1==id.ring()) ){
     rg = 0;
     secondRing = 3;
   }
   ch = id.chamber()-FirstCh;
 }

void CSCEfficiency::ringCandidates	(	int	station,
		float	absEta,
		std::map< std::string, bool > &	chamberTypes
	)

private

Definition at line 964 of file CSCEfficiency.cc.

                                                                                                   {
   // yeah, hardcoded again...
   switch (station){
   case 1:
     if(feta>0.85 && feta<1.18){//ME13
       chamberTypes["ME13"] = true;
     }      
     if(feta>1.18 && feta<1.7){//ME12
       chamberTypes["ME12"] = true;
     }
     if(feta>1.5 && feta<2.45){//ME11
       chamberTypes["ME11"] = true;
     }
     break;
   case 2:
     if(feta>0.95 && feta<1.6){//ME22
       chamberTypes["ME22"] = true;
       
     }  
     if(feta>1.55 && feta<2.45){//ME21
       chamberTypes["ME21"] = true;
     }
     break;
   case 3:
     if(feta>1.08 && feta<1.72){//ME32
       chamberTypes["ME32"] = true;
       
     }  
     if(feta>1.69 && feta<2.45){//ME31
       chamberTypes["ME31"] = true;
     }
     break;
   case 4:
     if(feta>1.78 && feta<2.45){//ME41
       chamberTypes["ME41"] = true;
     }
     break;
   default:
     break;
   }
 }

bool CSCEfficiency::stripWire_Efficiencies	(	CSCDetId &	cscDetId,
		FreeTrajectoryState &	ftsChamber
	)

private

Definition at line 1172 of file CSCEfficiency.cc.

References relativeConstraints::chamber, gather_cfg::cout, RecoEcal_EventContent_cff::ec, makeMuonMisalignmentScenario::endcap, FirstCh, FreeTrajectoryState::momentum(), relativeConstraints::ring, relativeConstraints::station, and PV3DBase< T, PVType, FrameType >::theta().

                                                                                         {
   int ec, st, rg, ch, secondRing;
   returnTypes(id, ec, st, rg, ch, secondRing);
 
   bool firstCondition, secondCondition;
   int missingLayers_s = 0;
   int missingLayers_wg = 0;
   for(int iLayer=0;iLayer<6;iLayer++){
  //----Strips
     if(printalot){ 
       printf("\t%i swEff: \tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",
          iLayer + 1,id.endcap(),id.station(),id.ring(),id.chamber(),iLayer+1);
       std::cout<<" size S = "<<allStrips[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh][iLayer].size()<<
     "size W = "<<allWG[id.endcap()-1][id.station()-1][id.ring()-1][id.chamber()-FirstCh][iLayer].size()<<std::endl;
 
     }
     firstCondition = allStrips[ec][st][rg][ch][iLayer].size() ? true : false;
     //allSegments[ec][st][rg][ch].size() ? true : false;
     secondCondition = false;
     if(secondRing>-1){
       secondCondition = allStrips[ec][st][secondRing][ch][iLayer].size() ? true : false;
     }
     if(firstCondition || secondCondition){
       ChHist[ec][st][rg][ch].EfficientStrips->Fill(iLayer+1);
     }
     else{
       if(printalot){ 
     std::cout<<"missing strips ";
     printf("\t\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber(),iLayer+1);
       }
     }
     // Wires
     firstCondition = allWG[ec][st][rg][ch][iLayer].size() ? true : false;
     secondCondition = false;
     if(secondRing>-1){
       secondCondition = allWG[ec][st][secondRing][ch][iLayer].size() ? true : false;
     }
     if(firstCondition || secondCondition){
       ChHist[ec][st][rg][ch].EfficientWireGroups->Fill(iLayer+1);
     }
     else{
       if(printalot){ 
     std::cout<<"missing wires ";
     printf("\t\tendcap/station/ring/chamber/layer: %i/%i/%i/%i/%i\n",id.endcap(),id.station(),id.ring(),id.chamber(),iLayer+1);
       }
     }
   }
   // Normalization
   if(6!=missingLayers_s){
     ChHist[ec][st][rg][ch].EfficientStrips->Fill(8);
   }
   if(6!=missingLayers_wg){
     ChHist[ec][st][rg][ch].EfficientWireGroups->Fill(8);
   }
   ChHist[ec][st][rg][ch].EfficientStrips->Fill(9);
   ChHist[ec][st][rg][ch].EfficientWireGroups->Fill(9);
 //
   ChHist[ec][st][rg][ch].StripWiresCorrelations->Fill(1);
   if(missingLayers_s!=missingLayers_wg){
     ChHist[ec][st][rg][ch].StripWiresCorrelations->Fill(2);
     if(6==missingLayers_wg){
       ChHist[ec][st][rg][ch].StripWiresCorrelations->Fill(3);
       ChHist[ec][st][rg][ch].NoWires_momTheta->Fill(ftsChamber.momentum().theta());
     }
     if(6==missingLayers_s){
       ChHist[ec][st][rg][ch].StripWiresCorrelations->Fill(4);
       ChHist[ec][st][rg][ch].NoStrips_momPhi->Fill(ftsChamber.momentum().theta());
     }
   }
   else if(6==missingLayers_s){
     ChHist[ec][st][rg][ch].StripWiresCorrelations->Fill(5);
   }
 
   return true;  
 }

Member Data Documentation

edm::EDGetTokenT<CSCALCTDigiCollection> CSCEfficiency::al_token

private

Definition at line 137 of file CSCEfficiency.h.

TH1F* CSCEfficiency::ALCTPerEvent

private

Definition at line 269 of file CSCEfficiency.h.

bool CSCEfficiency::allALCT[2][4][4][(36-1+1)]

private

Definition at line 191 of file CSCEfficiency.h.

bool CSCEfficiency::allCLCT[2][4][4][(36-1+1)]

private

Definition at line 190 of file CSCEfficiency.h.

bool CSCEfficiency::allCorrLCT[2][4][4][(36-1+1)]

private

Definition at line 192 of file CSCEfficiency.h.

std::vector<std::pair <LocalPoint, bool> > CSCEfficiency::allRechits[2][4][4][(36-1+1)][6]

private

Definition at line 206 of file CSCEfficiency.h.

std::vector<std::pair <LocalPoint, LocalVector> > CSCEfficiency::allSegments[2][4][4][(36-1+1)]

private

Definition at line 209 of file CSCEfficiency.h.

std::vector<std::pair <LocalPoint, int> > CSCEfficiency::allSimhits[2][4][4][(36-1+1)][6]

private

Definition at line 202 of file CSCEfficiency.h.

std::vector<std::pair <int, float> > CSCEfficiency::allStrips[2][4][4][(36-1+1)][6]

private

Definition at line 195 of file CSCEfficiency.h.

std::vector<std::pair <std::pair <int, float>, int> > CSCEfficiency::allWG[2][4][4][(36-1+1)][6]

private

Definition at line 198 of file CSCEfficiency.h.

bool CSCEfficiency::alongZ

private

Definition at line 184 of file CSCEfficiency.h.

bool CSCEfficiency::andOr

private

Definition at line 170 of file CSCEfficiency.h.

bool CSCEfficiency::applyIPangleCuts

private

Definition at line 161 of file CSCEfficiency.h.

struct CSCEfficiency::ChamberHistos CSCEfficiency::ChHist[2][4][3][ 36 - 1 +1]

private

edm::EDGetTokenT<CSCCLCTDigiCollection> CSCEfficiency::cl_token

private

Definition at line 138 of file CSCEfficiency.h.

TH1F* CSCEfficiency::CLCTPerEvent

private

Definition at line 270 of file CSCEfficiency.h.

edm::EDGetTokenT<CSCCorrelatedLCTDigiCollection> CSCEfficiency::co_token

private

Definition at line 139 of file CSCEfficiency.h.

TH1F* CSCEfficiency::DataFlow

private

Definition at line 266 of file CSCEfficiency.h.

double CSCEfficiency::distanceFromDeadZone

private

Definition at line 155 of file CSCEfficiency.h.

bool CSCEfficiency::emptyChambers[2][4][4][(36-1+1)]

private

Definition at line 212 of file CSCEfficiency.h.

bool CSCEfficiency::getAbsoluteEfficiency

private

Definition at line 153 of file CSCEfficiency.h.

edm::EDGetTokenT<edm::TriggerResults> CSCEfficiency::ht_token

private

Definition at line 146 of file CSCEfficiency.h.

bool CSCEfficiency::isBeamdata

private

Definition at line 152 of file CSCEfficiency.h.

bool CSCEfficiency::isData

private

Definition at line 150 of file CSCEfficiency.h.

bool CSCEfficiency::isIPdata

private

Definition at line 151 of file CSCEfficiency.h.

double CSCEfficiency::local_DX_DZ_Max

private

Definition at line 164 of file CSCEfficiency.h.

double CSCEfficiency::local_DY_DZ_Max

private

Definition at line 162 of file CSCEfficiency.h.

double CSCEfficiency::local_DY_DZ_Min

private

Definition at line 163 of file CSCEfficiency.h.

bool CSCEfficiency::magField

private

Definition at line 182 of file CSCEfficiency.h.

double CSCEfficiency::maxNormChi2

private

Definition at line 158 of file CSCEfficiency.h.

double CSCEfficiency::maxP

private

Definition at line 157 of file CSCEfficiency.h.

double CSCEfficiency::minP

private

Definition at line 156 of file CSCEfficiency.h.

unsigned int CSCEfficiency::minTrackHits

private

Definition at line 159 of file CSCEfficiency.h.

std::vector<std::string> CSCEfficiency::myTriggers

private

Definition at line 168 of file CSCEfficiency.h.

int CSCEfficiency::nEventsAnalyzed

private

Definition at line 180 of file CSCEfficiency.h.

bool CSCEfficiency::passTheEvent

private

Definition at line 186 of file CSCEfficiency.h.

std::vector<int> CSCEfficiency::pointToTriggers

private

Definition at line 169 of file CSCEfficiency.h.

bool CSCEfficiency::printalot

private

Definition at line 178 of file CSCEfficiency.h.

unsigned int CSCEfficiency::printout_NEvents

private

Definition at line 149 of file CSCEfficiency.h.

TH1F* CSCEfficiency::recHitsPerEvent

private

Definition at line 271 of file CSCEfficiency.h.

edm::EDGetTokenT<CSCRecHit2DCollection> CSCEfficiency::rh_token

private

Definition at line 140 of file CSCEfficiency.h.

std::string CSCEfficiency::rootFileName

private

Definition at line 131 of file CSCEfficiency.h.

edm::EDGetTokenT<CSCStripDigiCollection> CSCEfficiency::sd_token

private

Definition at line 136 of file CSCEfficiency.h.

edm::EDGetTokenT<CSCSegmentCollection> CSCEfficiency::se_token

private

Definition at line 141 of file CSCEfficiency.h.

TH1F* CSCEfficiency::segmentsPerEvent

private

Definition at line 272 of file CSCEfficiency.h.

edm::EDGetTokenT<edm::PSimHitContainer> CSCEfficiency::sh_token

private

Definition at line 144 of file CSCEfficiency.h.

struct CSCEfficiency::StationHistos CSCEfficiency::StHist[2][4]

private

TFile* CSCEfficiency::theFile

private

Definition at line 176 of file CSCEfficiency.h.

MuonServiceProxy* CSCEfficiency::theService

private

Definition at line 174 of file CSCEfficiency.h.

edm::EDGetTokenT<edm::View<reco::Track> > CSCEfficiency::tk_token

private

Definition at line 143 of file CSCEfficiency.h.

TH1F* CSCEfficiency::TriggersFired

private

Definition at line 267 of file CSCEfficiency.h.

bool CSCEfficiency::useDigis

private

Definition at line 154 of file CSCEfficiency.h.

bool CSCEfficiency::useTrigger

private

Definition at line 167 of file CSCEfficiency.h.

edm::EDGetTokenT<CSCWireDigiCollection> CSCEfficiency::wd_token

private

Definition at line 135 of file CSCEfficiency.h.

Classes

Public Member Functions

Private Member Functions

Private Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation