ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/CSCPriEff/src/TPTrackMuonSys.cc

Comparing UserCode/CSCPriEff/src/TPTrackMuonSys.cc (file contents):
Revision 1.3 by zhangjin, Thu Mar 14 21:09:48 2013 UTC vs.
Revision 1.7 by zhangjin, Tue May 7 21:12:35 2013 UTC

#	Line 65 | Line 65 | TPTrackMuonSys::TPTrackMuonSys(const edm
65
66
67		//////////////////// for LUTs
68	<
68	>	/*
69		bzero(srLUTs_, sizeof(srLUTs_));
70		//Int_t endcap =1, sector =1;
71		Bool_t TMB07=true;
#	Line 86 | Line 86 | TPTrackMuonSys::TPTrackMuonSys(const edm
86		} // stationItr loop
87		} // sectorItr loop
88		} // endcapItr loop
89	<
89	>	*/
90		// ------------------------------
91		// --- Summary ntuple booking ---
92		// ------------------------------
#	Line 253 | Line 253 | TPTrackMuonSys::TPTrackMuonSys(const edm
253		RunInfo->Branch("HLTDiMuName",&HLTDiMuName);
254		RunInfo->Branch("HLTDiMuObjModuleName",&HLTDiMuObjModuleName);
255		RunInfo->Branch("BadCSCChamberIndexList",&badChambersIndices);
256	–	fillChamberPosition();
256		}
257
258		TPTrackMuonSys::~TPTrackMuonSys(){
#	Line 298 | Line 297 | TPTrackMuonSys::analyze(const edm::Event
297		// setup.get<MuonGeometryRecord>().get(dtGeom);
298
299		setup.get<MuonGeometryRecord>().get(cscGeom);
300	<
300	>	#ifdef m_debug
301	>	//describe the chamber active region outline
302	>	if (nEventsAnalyzed==1) {
303	>	Short_t rings[]={11,14,12};//,12,13,21,31,41,22};
304	>	for (Byte_t iring=0;iring<3;iring++) {
305	>	CSCDetId id=CSCDetId(1, rings[iring]/10, rings[iring]%10, 1, 1);
306	>	printf("============ ME%d chamber outline=============\n",rings[iring]);
307	>	const CSCLayerGeometry *layerGeom = cscGeom->chamber(id)->layer(1)->geometry ();
308	>	const Byte_t NStrips=layerGeom->numberOfStrips(),NWires=layerGeom->numberOfWireGroups();
309	>	LocalPoint interSect_ = layerGeom->stripWireGroupIntersection(1, 1);
310	>	printf("(strip 1, wire group 1) at (%.3f,%.3f)\n",interSect_.x(),interSect_.y());
311	>	interSect_ = layerGeom->stripWireGroupIntersection(NStrips, 1);
312	>	printf("(strip %d, wire group 1) at (%.3f,%.3f)\n",NStrips,interSect_.x(),interSect_.y());
313	>	interSect_ = layerGeom->stripWireGroupIntersection(NStrips, NWires);
314	>	printf("(strip %d, wire group %d) at (%.3f,%.3f)\n",NStrips, NWires,interSect_.x(),interSect_.y());
315	>	interSect_ = layerGeom->stripWireGroupIntersection(1, NWires);
316	>	printf("(strip 1, wire group %d) at (%.3f,%.3f)\n", NWires,interSect_.x(),interSect_.y());
317	>	printf("   ======== middle x of wires(check)==========\n");
318	>	const CSCWireTopology* wireTopology = layerGeom->wireTopology();
319	>	for (Short_t wireGroup_id=10;wireGroup_id<=15;wireGroup_id++) {
320	>	Float_t wideWidth      = wireTopology->wideWidthOfPlane();
321	>	Float_t narrowWidth    = wireTopology->narrowWidthOfPlane();
322	>	Float_t length         = wireTopology->lengthOfPlane();
323	>	Float_t tangent = (wideWidth-narrowWidth)/(2.*length);
324	>	Float_t wireangle= wireTopology->wireAngle();
325	>	std::vector<float> wirecenters=wireTopology->wireValues(wireTopology->middleWireOfGroup(wireGroup_id));
326	>	printf("x center of wire %d is at %.3f=%.3f\n",wireGroup_id,wireTopology->wireValues(wireTopology->middleWireOfGroup(wireGroup_id))[0],wirecenters[2]*sin(wireangle)*tangent/2);
327	>	}
328	>	printf("y center of the wires is at %.3f\n",wireTopology->yOfWire(1)+0.5*wireTopology->lengthOfPlane());
329	>	}
330	>	}
331	>	#endif
332		// Get the propagators
333		setup.get<TrackingComponentsRecord>().get("SmartPropagatorAnyRK", propagatorAlong   );
334		setup.get<TrackingComponentsRecord>().get("SmartPropagatorAnyOpposite", propagatorOpposite);
#	Line 320 | Line 350 | TPTrackMuonSys::analyze(const edm::Event
350		// get CSC segment collection
351		event.getByLabel("cscSegments", cscSegments);
352
323	–	setup.get<MuonGeometryRecord>().get(rpcGeo);
324	–
325	–	edm::Handle<RPCRecHitCollection> rpcRecHits;
326	–	event.getByLabel("rpcRecHits",rpcRecHits);
327	–
328	–	//edm::Handle<DTRecSegment4DCollection> dtSegments;
329	–	//event.getByLabel( "dt4DSegments", dtSegments );
330	–
353		edm::Handle<CSCCorrelatedLCTDigiCollection> mpclcts;
354		try{
355		event.getByLabel("csctfunpacker","", mpclcts);
#	Line 335 | Line 357 | TPTrackMuonSys::analyze(const edm::Event
357		edm::LogError("")<< "Error! Can't get m_gTracksTag by label. ";
358		}
359
360	+	/*
361	+	// get RPC and DT
362	+	setup.get<MuonGeometryRecord>().get(rpcGeo);
363	+
364	+	edm::Handle<RPCRecHitCollection> rpcRecHits;
365	+	event.getByLabel("rpcRecHits",rpcRecHits);
366	+
367	+	edm::Handle<DTRecSegment4DCollection> dtSegments;
368	+	event.getByLabel( "dt4DSegments", dtSegments );
369	+	*/
370	+
371		//  muon trigger scales
372		edm::ESHandle<L1MuTriggerScales> trigscales_h;
373		setup.get<L1MuTriggerScalesRcd> ().get(trigscales_h);
#	Line 959 | Line 992 | TPTrackMuonSys::analyze(const edm::Event
992		Int_t st = j>2?j-2:0, ec=CSCEndCapPlus?1:2; // determine the station number...
993		Float_t trkEta = tsos.globalPosition().eta(), trkPhi = tsos.globalPosition().phi();
994
995	<	Int_t rg = ringCandidate(st+1, trkEta, trkPhi);
996	<	if (rg==-9999) continue;
995	>	Short_t rg = ringCandidate(st+1, trkEta, trkPhi);
996	>	if ( rg==-9999) continue;
997	>	if (!(
998	>	(j==0&&(rg==1||rg==4)) ||
999	>	(j==1&&rg==2) ||
1000	>	(j==2&&rg==3) ||
1001	>	(j>2)
1002	>	)) continue;//check if it is runing on the right ring for a certain Z
1003	>	else CSCRg[st]=rg;
1004
1005		//for chamber overlap region
1006		/* This part of code is useless. Just directly extrapolate to chamber
#	Line 971 | Line 1011 | TPTrackMuonSys::analyze(const edm::Event
1011		if(!CSCEndCapPlus)zzPlaneME = -zzPlaneME;
1012
1013		tsos = surfExtrapTrkSam(trackRef, zzPlaneME);
1014	<	if (!tsos.isValid()) continue;*/
975	<
1014	>	if (!tsos.isValid()) continue;
1015		trkEta = tsos.globalPosition().eta(); trkPhi =  tsos.globalPosition().phi();
1016	+	rg = ringCandidate(st+1, trkEta, trkPhi);
1017	+	if ( rg==-9999 ) continue;
1018	+	*/
1019	+
1020	+	//Calculate the solid angle between the muon and this extrapolated track
1021		if(trkPhi  < 0) trkPhi += 2*M_PI;
1022		/*dR between muon and track*/
1023		Float_t MPhi = phiMuVec[st];
1024	<	if(MPhi < 0 ) MPhi +=  2*M_PI;
1024	>	if (MPhi < 0 ) MPhi +=  2*M_PI;
1025		Float_t TPhi = tsos.globalPosition().phi();
1026	<	if(TPhi < 0 ) TPhi += 2*M_PI;
1026	>	if (TPhi < 0 ) TPhi += 2*M_PI;
1027		dRTkMu[st] = deltaR( etaMuVec[st], MPhi, tsos.globalPosition().eta() , TPhi);
1028		//      printf("Mu(%f,%f),Trk(%f,%f)-->dR(%f)",etaMuVec[st], MPhi, tsos.globalPosition().eta() , TPhi,dRTkMu[st]);
1029	<
1030	<	CSCRg[st] = ringCandidate(st+1, trkEta, trkPhi);
987	<	if ( CSCRg[st]==0 ) continue;
988	<
989	<	CSCChCand[st] = thisChamberCandidate(st+1, CSCRg[st], trkPhi);
1029	>	//Find the chamber candidate. If two chamber overlaps, simply divide them from middle. e.g. if chamber1 is from phi=1-3, chamber 2 is from phi=2-4. Then phi=1-2.5-->chamber 1; phi=2.5-4-->chamber 2. Since the distribution over phi is uniform, it won't bias the result.
1030	>	CSCChCand[st] = thisChamberCandidate(st+1, rg, trkPhi);
1031		CSCDetId Layer0Id=CSCDetId(ec, st+1, rg,  CSCChCand[st], 0);//layer 0 is the mid point of the chamber. It is not a real layer.
1032		CSCChBad[st] = badChambers_->isInBadChamber( Layer0Id );
1033	+	//skip not-existing ME42 chambers
1034	+	if (CSCChBad[st]&&st==3&&rg==2) continue;
1035		#ifdef jz_debug
1036	<	if (CSCChBad[st]) cerr<<(CSCEndCapPlus?"ME+":"ME-")<<st+1<<"/"<<CSCRg[st]<<"/"<<CSCChBad[st]<<" is a dead chamber."<<endl;
1036	>	if (CSCChBad[st]) cerr<<(CSCEndCapPlus?"ME+":"ME-")<<st+1<<"/"<<rg<<"/"<<CSCChBad[st]<<" is a dead chamber."<<endl;
1037		#endif
1038	+	/*Not necessary. Chamber search using phi is accurate enough.
1039	+	CSCDetId Layerid  = CSCDetId( ec, st+1, rg,  CSCChCand[st], 1 );
1040	+	vector<Float_t> EdgeAndDistToGap( GetEdgeAndDistToGap(trackRef,Layerid) );//values: 1-edge;2-err of edge;3-disttogap;4-err of dist to gap
1041	+	if (EdgeAndDistToGap[0]>0.||EdgeAndDistToGap[0]<-9990.) {
1042	+	cerr<<"ch"<<CSCChCand[st]<<", oldedge="<<EdgeAndDistToGap[0]<<"cm --> ";
1043	+	//try neighborhood chambers
1044	+	CSCDetId Layerid_plusone  = CSCDetId( ec, st+1, rg,  CSCChCand[st]+1, 1 );
1045	+	vector<Float_t> EdgeAndDistToGap_plusone( GetEdgeAndDistToGap(trackRef,Layerid_plusone) );
1046	+	CSCDetId Layerid_minusone  = CSCDetId( ec, st+1, rg,  CSCChCand[st]-1, 1 );
1047	+	vector<Float_t> EdgeAndDistToGap_minusone( GetEdgeAndDistToGap(trackRef,Layerid_minusone) );
1048	+	if (EdgeAndDistToGap_plusone[0]<EdgeAndDistToGap[0]&&EdgeAndDistToGap_plusone[0]>-9990.) {
1049	+	if (EdgeAndDistToGap_minusone[0]>EdgeAndDistToGap_plusone[0]) CSCChCand[st]+=1;
1050	+	else { if (EdgeAndDistToGap_minusone[0]>-9990.) CSCChCand[st]-=1;}
1051	+	}
1052	+	else { if (EdgeAndDistToGap_minusone[0]<EdgeAndDistToGap[0]&&EdgeAndDistToGap_minusone[0]>-9990.) CSCChCand[st]-=1;}
1053	+	cerr<<"ch"<<CSCChCand[st]<<", plusoneedge="<<EdgeAndDistToGap_plusone[0]<<"cm,"<<", minusoneedge="<<EdgeAndDistToGap_minusone[0]<<"cm."<<endl;
1054	+	}*/
1055		for (Int_t ly=1;ly<7;ly++) {
1056	<	CSCDetId ID  = CSCDetId( ec, st+1, rg,  CSCChCand[st], ly );
1057	<	vector<Float_t> EdgeAndDistToGap( GetEdgeAndDistToGap(trackRef,ID) );//values: 1-edge;2-err of edge;3-disttogap;4-err of dist to gap
1056	>	CSCDetId Layerid = CSCDetId( ec, st+1, rg,  CSCChCand[st], ly );
1057	>	vector<Float_t> EdgeAndDistToGap( GetEdgeAndDistToGap(trackRef,Layerid) );//values: 1-edge;2-err of edge;3-disttogap;4-err of dist to gap
1058		if (EdgeAndDistToGap[0]>CSCProjEdgeDist[st]) {
1059		CSCProjEdgeDist[st]=EdgeAndDistToGap[0];
1060		CSCProjEdgeDistErr[st]=EdgeAndDistToGap[1];
#	Line 1058 | Line 1118 | TPTrackMuonSys::analyze(const edm::Event
1118		CSCdYdZTTSeg[st] = dydz_trk - dydz_seg;
1119		}
1120		nTrkCountCSCSeg++;
1121	+	delete TrajToSeg;
1122		}
1123
1124		////// Loop over MPC infromation to look for LCTs....
1125	<	CSCDetId Layer3id  = CSCDetId( ec, st+1, rg,  CSCChCand[st], 3 );//go to layer 3 where corresponds to the LCTPos
1125	>	CSCDetId Layer3id  = CSCDetId( ec, st+1, rg,  CSCChCand[st], 3 );//go to layer 3 that corresponds to the LCTPos
1126		const GeomDet* Layer3gdet=cscGeom->idToDet(Layer3id);
1127		tsos=surfExtrapTrkSam(trackRef, Layer3gdet->surface().position().z());
1128		if (tsos.isValid()) {
#	Line 1076 | Line 1137 | TPTrackMuonSys::analyze(const edm::Event
1137		LocalError localTTErr =tsos.localError().positionError();
1138		CSCDxErrTTLCT[st] = sqrt(localTTErr.xx()); CSCDyErrTTLCT[st] = sqrt(localTTErr.yy());
1139		CSCDxyErrTTLCT[st] =sqrt(pow(CSCDxTTLCT[st],2)*localTTErr.xx() + pow(CSCDyTTLCT[st],2)*localTTErr.yy())/CSCDxyTTLCT[st];
1140	<	//cerr<<"ME"<<st+1<<CSCRg[st]<<":TTLCTx,TTLCTy,TTLCTxy:"<<CSCDxTTLCT[st]<<","<<CSCDyTTLCT[st]<<","<<CSCDxyTTLCT[st]<<endl;
1141	<	//cerr<<"segZ_TTy,LCTZ_TTy:"<<CSCSegyLc[st]-CSCDyTTSeg[st]<<","<<localL3pCSC.y()<<";"<<endl;
1142	<	}
1143	<	}
1140	>	delete LCTPos;
1141	>	//cerr<<"segZ_TTy,LCTZ_TTy:"<<CSCSegyLc[st]-CSCDyTTSeg[st]<<","<<localL3pCSC.y()<<";"<<endl;
1142	>	}//end of if found LCT
1143	>	}//end of if layer 3 extrapolation is available
1144		} // for loop for the stations -- j
1145		/*
1146		for (Byte_t st=0;st<4;st++) {
1147		lctSt[st] = ( CSCDxyTTLCT[st] >0. && CSCDxyTTLCT[st] < 40 )?1:0;
1148		segSt[st] = ( CSCDxyTTSeg[st] >0. && CSCDxyTTSeg[st] < 40)?1:0;
1149		}*/
1150	+	//cerr<<"eta="<<tracks_eta<<":ME1"<<CSCRg[0]<<":TTLCTx,TTLCTy,TTLCTxy:"<<CSCDxTTLCT[0]<<","<<CSCDyTTLCT[0]<<","<<CSCDxyTTLCT[0]<<endl;
1151		fractNtuple->Fill();
1152		firsttrackmatchingtoMuTag=false;
1153		} // loop over tracks...
#	Line 1272 | Line 1334 | vector<Float_t> TPTrackMuonSys::GetEdgeA
1334		TrajectoryStateOnSurface tsos=surfExtrapTrkSam(trackRef, gdet->surface().position().z());
1335		if (!tsos.isValid()) return result;
1336		LocalPoint localTTPos = gdet->surface().toLocal(tsos.freeState()->position());
1337	<
1276	<	const CSCChamberSpecs* chamberSpecs = cscGeom->chamber(detid)->specs();
1277	<	const CSCLayerGeometry* layerGeometry = (detid.layer()%2==0)?chamberSpecs->evenLayerGeometry(detid.endcap()):chamberSpecs->oddLayerGeometry(detid.endcap());
1278	<	const CSCWireTopology* wireTopology = layerGeometry->wireTopology();
1337	>	const CSCWireTopology* wireTopology = cscGeom->layer(detid)->geometry()->wireTopology();
1338		Float_t wideWidth      = wireTopology->wideWidthOfPlane();
1339		Float_t narrowWidth    = wireTopology->narrowWidthOfPlane();
1340		Float_t length         = wireTopology->lengthOfPlane();
#	Line 1289 | Line 1348 | vector<Float_t> TPTrackMuonSys::GetEdgeA
1348		Float_t yPrime  = localTTPos.y()+fabs(yOfFirstWire);
1349		// half trapezoid width at y' is 0.5 * narrowWidth + x side of triangle with the above tangent and side y'
1350		Float_t halfWidthAtYPrime = 0.5*narrowWidth+yPrime*tangent;
1351	+	// x offset between local origin and symmetry center of wire plane is zero
1352	+	// x offset of ME11s is also zero. x center of wire groups is not at zero, because it is not parallel to x. The wire groups of ME11s have a complex geometry, see the code in m_debug.
1353		Float_t edgex = fabs(localTTPos.x()) - halfWidthAtYPrime;
1354		Float_t edgey = fabs(localTTPos.y()-yCOWPOffset) - 0.5*length;
1355		LocalError localTTErr = tsos.localError().positionError();
#	Line 1394 | Line 1455 | Bool_t TPTrackMuonSys::matchTTwithCSCRec
1455		}
1456
1457		////// Match TT with RPCEndCapHit
1458	+	/*
1459		Bool_t TPTrackMuonSys::matchTTwithRPCEChit(Bool_t trackDir,
1460		Int_t j,
1461		reco::TrackRef trackRef,
#	Line 1448 | Line 1510 | Bool_t TPTrackMuonSys::matchTTwithRPCECh
1510
1511		return aMatch;
1512		}
1513	<
1513	>	*/
1514		//////////////  Get the matching with LCTs...
1515		LocalPoint * TPTrackMuonSys::matchTTwithLCTs(Float_t xPos, Float_t yPos, Short_t ec, Short_t st, Short_t &rg, Short_t cham,
1516		edm::Handle<CSCCorrelatedLCTDigiCollection> mpclcts, Float_t &dRTrkLCT, Int_t &lctBX ) {
#	Line 1456 | Line 1518 | LocalPoint * TPTrackMuonSys::matchTTwith
1518
1519		for (CSCCorrelatedLCTDigiCollection::DigiRangeIterator detMPCUnitIt = mpclcts->begin();
1520		detMPCUnitIt != mpclcts->end(); detMPCUnitIt++) {
1521	<	const CSCDetId& id = (*detMPCUnitIt).first;
1521	>	CSCDetId id = (*detMPCUnitIt).first;
1522
1523		if(ec != id.endcap())continue;
1524		if(st != id.station())continue;
#	Line 1465 | Line 1527 | LocalPoint * TPTrackMuonSys::matchTTwith
1527		Bool_t ed1 = (st == 1) && ((rg == 1 || rg == 4) && (id.ring() == 1 || id.ring() == 4));
1528		Bool_t ed2 = (st == 1) && ((rg == 2 && id.ring() == 2) || (rg == 3 && id.ring() == 3));
1529		Bool_t ed3 = (st != 1) && (rg == id.ring());
1530	<	Bool_t TMCSCMatch = (ed1 || ed2 || ed3);
1469	<	if(! TMCSCMatch)continue;
1470	<
1530	>	if ( !(ed1 || ed2 || ed3) ) continue;
1531		const CSCCorrelatedLCTDigiCollection::Range& MPCrange = (*detMPCUnitIt).second;
1532		for (CSCCorrelatedLCTDigiCollection::const_iterator mpcIt = MPCrange.first; mpcIt != MPCrange.second; mpcIt++) {
1533		Bool_t lct_valid = (*mpcIt).isValid();
1534		if(!lct_valid)continue;
1535		//In CSC offline/hlt software in general, is COUNT FROM ONE, such as CSCGeometry.
1536	<	//However, the LCT software counts from zero. strip_id is from 0 to 159. wireGroup_id is from 0 to 79. So here we need to plus one.
1537	<	Int_t wireGroup_id = (*mpcIt).getKeyWG()+1;
1538	<	Float_t strip_id = (*mpcIt).getStrip()/2.+1;
1536	>	//However, the LCT software counts from zero. strip_id is from 0 to 159. wireGroup_id is from 0 to 31(ME13),47(ME11),63(ME1234/2),95(ME31,41),111(ME21). So here we need to plus one.
1537	>	Byte_t wireGroup_id = (*mpcIt).getKeyWG()+1;
1538	>	Byte_t strip_id=(*mpcIt).getStrip()/2+1;
1539	>	//if ( id.ring()==4 ) cerr<<"Alert id.ring()==4"<<endl;
1540	>	Bool_t me11=(st == 1) && (id.ring() == 1 || id.ring() == 4),
1541	>	me11a = me11 && strip_id>64;
1542	>	//me11b = me11 && strip_id<=64;
1543	>	//http://cmssdt.cern.ch/SDT/lxr/source/CondFormats/CSCObjects/src/CSCChannelTranslator.cc
1544	>	// Translate a raw strip channel in range 1-80, iraw,  into
1545	>	// corresponding geometry-oriented channel in which increasing
1546	>	// channel number <-> strip number increasing with +ve local x.
1547	>	//geomStripChannel( CSCDetId(ec, st+1, rg,  CSCChCand[st], 0) )
1548	>	if ( me11a ) {
1549	>	strip_id-=64;
1550	>	// The CSCCorrelatedLCTDigi DetId does NOT distinguish ME11A and B. All of the DetIDs are labelled as ME11B (all ME11, none ME14)
1551	>	// However, stripWireGroupIntersection must know that since ME11A has 48 strips and ME11B has 64.
1552	>	id=CSCDetId(ec, 1, 4, cham, 3);
1553	>	//if ( id.endcap()==1 ) strip_id=17-strip_id;
1554	>	}
1555	>	//if ( me11b && id.endcap()!=1 ) strip_id=65-strip_id;
1556		const CSCLayerGeometry *layerGeom = cscGeom->chamber(id)->layer (3)->geometry ();
1557	<	const Int_t Nstrips=layerGeom->numberOfStrips();
1558	<	Bool_t me11a = ( (st == 1) && (id.ring() == 1 || id.ring() == 4) && strip_id>Nstrips );
1559	<	if ( me11a ) strip_id-=Nstrips;
1483	<	for(Int_t ii = 0; ii < 3; ii++){
1484	<	if ( strip_id>Nstrips ) LogWarning("Strip_id") << "Got "<<strip_id<<", but there are "<< Nstrips <<" strips in total." <<m_hlt.process();
1485	<	LocalPoint interSect_ = layerGeom->stripWireGroupIntersection(Int_t(strip_id), wireGroup_id);
1557	>	for(Byte_t ii = 0; ii < 3; ii++){
1558	>	// if ( strip_id>64 ) LogWarning("Strip_id") << "Got "<<strip_id<<", but there are "<< Nstrips <<" strips in total." <<m_hlt.process();
1559	>	LocalPoint interSect_ = layerGeom->stripWireGroupIntersection(strip_id, wireGroup_id);
1560		//      printf( "ME%d/%d: %.1f/%d, %d/%d: xLCT-xTT=%.2f-%.2f; yLCT-yTT=%.2f-%.2f \n",st,id.ring(),strip_id,Nstrips,wireGroup_id,layerGeom->numberOfWireGroups(),interSect_.x(),xPos,interSect_.y(),yPos);
1561		Float_t DeltaR_ = sqrt(pow((interSect_.x()-xPos),2) + pow((interSect_.y()-yPos),2));
1562		if( DeltaR_ < fabs(dRTrkLCT) ) {
1563	+	delete interSect;
1564		interSect=new LocalPoint(interSect_);
1565		dRTrkLCT =  DeltaR_ ;
1566		lctBX = (*mpcIt).getBX();
#	Line 1493 | Line 1568 | LocalPoint * TPTrackMuonSys::matchTTwith
1568		else rg=id.ring();
1569		//cout << "1: BX = " << (*mpcIt).getBX() << " BX0 = " << (*mpcIt).getBX0() << std::endl;
1570		} // for the matching if statement...
1571	<	if (me11a) strip_id+=16.;
1571	>	if (me11a) strip_id+=16;
1572		else break;
1573		}// end iteration over of ME11A triplet
1574		}// end iteration over lcts_mpc_data in one chamber
#	Line 1514 | Line 1589 | inline Float_t TPTrackMuonSys::Trajector
1589		TrajectoryStateOnSurface *TPTrackMuonSys::matchTTwithCSCSeg( reco::TrackRef trackRef, edm::Handle<CSCSegmentCollection> cscSegments,
1590		CSCSegmentCollection::const_iterator &cscSegOut, CSCDetId & idCSC ) {
1591		TrajectoryStateOnSurface *TrajSuf=NULL;
1592	<	Float_t deltaCSCR = 50.0;
1592	>	Float_t deltaCSCR = 9999.;
1593		for(CSCSegmentCollection::const_iterator segIt=cscSegments->begin(); segIt != cscSegments->end(); segIt++) {
1594		CSCDetId id  = (CSCDetId)(*segIt).cscDetId();
1595		if(idCSC.endcap() != id.endcap())continue;
#	Line 1533 | Line 1608 | TrajectoryStateOnSurface *TPTrackMuonSys
1608
1609		Float_t dR_= fabs( TrajectoryDistToSeg( &TrajSuf_, segIt ) );
1610		if ( dR_ < deltaCSCR ){
1611	+	delete TrajSuf;
1612		TrajSuf=new TrajectoryStateOnSurface(TrajSuf_);
1613		deltaCSCR = dR_;
1614		cscSegOut = segIt;
#	Line 1733 | Line 1809 | void TPTrackMuonSys::getCSCSegWkeyHalfSt
1809		#endif
1810		}
1811		#endif// end of GetCSCHitsBefore500
1812	<
1812	>	/*
1813		Bool_t TPTrackMuonSys::matchCSCSegWithLCT(edm::Handle<CSCCorrelatedLCTDigiCollection> mpclcts,
1814		CSCDetId & idCSC,
1815		Int_t TT,
#	Line 1853 | Line 1929 | Bool_t TPTrackMuonSys::matchCSCSegWithLC
1929		}// end iteration over lcts_mpc_data
1930		return (xMatch[0] || xMatch[1]);
1931		}
1932	<
1932	>	*/
1933		void TPTrackMuonSys::chamberCandidates(Int_t station, Float_t feta, Float_t phi, std::vector <int> &coupleOfChambers){
1934		// yeah, hardcoded again...
1935		coupleOfChambers.clear();
#	Line 1887 | Line 1963 | void TPTrackMuonSys::chamberCandidates(I
1963
1964
1965		Int_t TPTrackMuonSys::ringCandidate(Int_t station, Float_t feta, Float_t phi){
1890	–	// yeah, hardcoded again...
1891	–
1892	–	Int_t ring = -9999;
1893	–
1966		switch (station){
1967		case 1:
1968		if(fabs(feta)>=0.85 && fabs(feta)<1.18){//ME13
1969	<	ring = 3;
1969	>	return 3;
1970		}
1971		if(fabs(feta)>=1.18 && fabs(feta)<=1.5){//ME12 if(fabs(feta)>1.18 && fabs(feta)<1.7){//ME12
1972	<	ring = 2;
1972	>	return 2;
1973		}
1974	<	if(fabs(feta)>1.5 && fabs(feta)<2.45){//ME11
1975	<	ring = 1; //or 4;
1974	>	if(fabs(feta)>1.5 && fabs(feta)<2.1){//ME11
1975	>	return 1;
1976	>	}
1977	>	if(fabs(feta)>=2.1 && fabs(feta)<2.45){//ME11
1978	>	return 4;
1979		}
1980		break;
1981		case 2:
1982		if(fabs(feta)>0.95 && fabs(feta)<1.6){//ME22
1983	<	ring = 2;
1983	>	return 2;
1984		}
1985		if(fabs(feta)>1.55 && fabs(feta)<2.45){//ME21
1986	<	ring = 1;
1986	>	return 1;
1987		}
1988		break;
1989		case 3:
1990		if(fabs(feta)>1.08 && fabs(feta)<1.72){//ME32
1991	<	ring = 2;
1991	>	return 2;
1992		}
1993		if(fabs(feta)>1.69 && fabs(feta)<2.45){//ME31
1994	<	ring = 1;
1994	>	return 1;
1995		}
1996		break;
1997		case 4:
1998		if(fabs(feta)>1.78 && fabs(feta) <2.45){//ME41
1999	<	ring = 1;
1999	>	return 1;
2000		}
2001		if(fabs(feta)>1.15 && fabs(feta) <=1.78){//ME42
2002	<	ring = 2;
2002	>	return 2;
2003		}
2004		break;
2005		default:
2006		LogError("")<<"Invalid station: "<<station<<endl;
2007		break;
2008		}
2009	<	return ring;
2009	>	return -9999;
2010		}
2011
2012		Short_t TPTrackMuonSys::thisChamberCandidate(Short_t station, Short_t ring, Float_t phi){
2013	<
2014	<	//  double minDelPhi = 0.17; Int_t station, Int_t ring, Float_t phi
2015	<
2016	<	Float_t *MyPhiValues=NULL;
2017	<	//36 chambers
2018	<	if (station == 1)  MyPhiValues = StationOnePhi;
2019	<	if (station == 2 && ring == 2) MyPhiValues = StationTwoTwoPhi;
2020	<	if (station == 3 && ring == 2) MyPhiValues = StationThreeTwoPhi;
2021	<	if (station == 4 && ring == 2) MyPhiValues = StationFourTwoPhi;
2022	<	//18 chambers
1948	<	if(station == 2 && ring == 1) MyPhiValues = StationTwoOnePhi;
1949	<	if(station == 3 && ring == 1) MyPhiValues = StationThreeOnePhi;
1950	<	if(station == 4 && ring == 1) MyPhiValues = StationFourOnePhi;
1951	<
1952	<	if ( MyPhiValues==NULL ) {
1953	<	LogError("")<<"Invalid station and ring: "<<station<<"/"<<ring<<endl;
1954	<	return 0;
1955	<	}
1956	<
1957	<	Int_t iCh=0;
1958	<	if(phi < 0.) phi = phi + 2*M_PI;
1959	<
1960	<	Short_t nVal = 36;
1961	<	if(station != 1 && ring == 1) nVal = 18;
1962	<
1963	<	double minDelPhi = 100000.;
1964	<	for (Short_t i = 0; i < nVal; i++){
1965	<	double dphi = deltaPhi(MyPhiValues[i], phi);
1966	<	if(fabs(dphi) < minDelPhi){
1967	<	minDelPhi = fabs(dphi);
1968	<	iCh = i+2;
1969	<	}
1970	<	}
1971	<	if(iCh > nVal)iCh = 1;
1972	<	return iCh;
1973	<	}
1974	<
1975	<
1976	<	void TPTrackMuonSys::fillChamberPosition(){
1977	<	TVector3 x(0,0,1);
1978	<	for (Int_t i = 0; i < 36; i++){
1979	<	StationOnePhi[i] = 0.0;
1980	<	StationTwoTwoPhi[i] = 0.0;
1981	<	StationThreeTwoPhi[i] = 0.0;
1982	<	StationFourTwoPhi[i] = 0.0;
1983	<	if(i < 18){
1984	<	StationTwoOnePhi[i] = 0.0;
1985	<	StationThreeOnePhi[i] = 0.0;
1986	<	StationFourOnePhi[i] = 0.0;
1987	<	}
1988	<	}
1989	<
1990	<	/////////////////////////////////
1991	<
1992	<	TVector3 pp(180.5,0,0);
1993	<	for (Int_t i = 0; i < 36; i++){
1994	<	pp.Rotate(2*M_PI/36,x);
1995	<	StationOnePhi[i]  = pp.Phi();
1996	<	if(StationOnePhi[i] < 0)StationOnePhi[i] = StationOnePhi[i] + 2*M_PI;
1997	<	#ifdef m_debug
1998	<	std::cout << " PHI(1-" << i << ")" << StationOnePhi[i];
1999	<	#endif
2000	<	}
2001	<	std::cout << std::endl;
2002	<
2003	<	TVector3 pp1(241.73,0,0);
2004	<	pp1.Rotate(M_PI/36 ,x);
2005	<	for (Int_t i = 0; i < 18; i++){
2006	<	pp1.Rotate(2*M_PI/18,x);
2007	<	StationTwoOnePhi[i]  = pp1.Phi();
2008	<	if(StationTwoOnePhi[i] < 0)StationTwoOnePhi[i] = StationTwoOnePhi[i] + 2*M_PI;
2009	<	#ifdef m_debug
2010	<	std::cout << " PHI(2/1-" << i << ")" << StationTwoOnePhi[i];
2011	<	#endif
2012	<	}
2013	<	std::cout << std::endl;
2014	<
2015	<	TVector3 pp2(707.56,0,0);
2016	<	for (Int_t i = 0; i < 36; i++){
2017	<	pp2.Rotate(2*M_PI/36,x);
2018	<	StationTwoTwoPhi[i]  = pp2.Phi();
2019	<	if(StationTwoTwoPhi[i] < 0)StationTwoTwoPhi[i] = StationTwoTwoPhi[i] + 2*M_PI;
2020	<	#ifdef m_debug
2021	<	std::cout << " PHI(2/2-" << i << ")" << StationTwoTwoPhi[i];
2022	<	#endif
2023	<	}
2024	<	std::cout << std::endl;
2025	<
2026	<	TVector3 pp3(251.74, 0.0, 0.0);
2027	<	pp3.Rotate(M_PI/36,x);
2028	<	//pp3.Rotate(0.0511,x);
2029	<	for (Int_t i = 0; i < 18; i++){
2030	<	pp3.Rotate(2*M_PI/18,x);
2031	<	StationThreeOnePhi[i]  = pp3.Phi();
2032	<	if(StationThreeOnePhi[i] < 0)StationThreeOnePhi[i] = StationThreeOnePhi[i] + 2*M_PI;
2033	<	#ifdef m_debug
2034	<	std::cout << " PHI(3/1-" << i << ")" << StationThreeOnePhi[i];
2035	<	#endif
2036	<	}
2037	<	std::cout << std::endl;
2038	<
2039	<	TVector3 pp4(525.55,0.0,0);
2040	<	for (Int_t i = 0; i < 36; i++){
2041	<	pp4.Rotate(2*M_PI/36,x);
2042	<	StationThreeTwoPhi[i]  = pp4.Phi();
2043	<	if(StationThreeTwoPhi[i] < 0)StationThreeTwoPhi[i] = StationThreeTwoPhi[i] + 2*M_PI;
2044	<	#ifdef m_debug
2045	<	std::cout << " PHI(3/2-" << i << ")" << StationThreeTwoPhi[i];
2046	<	#endif
2047	<	}
2048	<	std::cout << std::endl;
2049	<
2050	<	TVector3 pp5(261.7,0.,0.);
2051	<	pp5.Rotate(M_PI/36,x);
2052	<	for (Int_t i = 0; i < 18; i++){
2053	<	pp5.Rotate(2*M_PI/18,x);
2054	<	StationFourOnePhi[i]  = pp5.Phi();
2055	<	if(StationFourOnePhi[i] < 0)StationFourOnePhi[i] = StationFourOnePhi[i] + 2*M_PI;
2056	<	#ifdef m_debug
2057	<	std::cout << " PHI(4/1-" << i << ")" << StationFourOnePhi[i];
2058	<	#endif
2059	<	}
2060	<	std::cout << std::endl;
2013	>	//search for chamber candidate based on CMS IN-2007/024
2014	>	//10 deg chambers are ME1,ME22,ME32,ME42 chambers; 20 deg chambers are ME21,31,41 chambers
2015	>	//Chambers one always starts from approx -5 deg.
2016	>	const Short_t nVal = (station>1 && ring == 1)?18:36;
2017	>	const Float_t ChamberSpan=2*M_PI/nVal;
2018	>	Float_t dphi = phi + M_PI/36;
2019	>	while (dphi >= 2*M_PI) dphi -= 2*M_PI;
2020	>	while (dphi < 0) dphi += 2*M_PI;
2021	>	Short_t ChCand=floor(dphi/ChamberSpan)+1;
2022	>	return ChCand>nVal?nVal:ChCand;
2023		}
2024
2025		///////////////////////
#	Line 2396 | Line 2358 | TPTrackMuonSys::MCParticleInfo TPTrackMu
2358		return TBA;
2359		}
2360
2361	<	//deadzone center is according to http://cmslxr.fnal.gov/lxr/source/RecoLocalMuon/CSCEfficiency/src/CSCEfficiency.cc#605
2361	>	//deadzone center is according to http://cmssdt.cern.ch/SDT/lxr/source/RecoLocalMuon/CSCEfficiency/src/CSCEfficiency.cc#605
2362		//wire spacing is according to CSCTDR
2363		Float_t TPTrackMuonSys::YDistToHVDeadZone(Float_t yLocal, Int_t StationAndRing){
2364	<	const Float_t deadZoneCenterME1_14[2] = {-81.0,81.0};
2365	<	const Float_t deadZoneCenterME1_2[4] = {-86.285,-32.88305,32.867423,88.205};
2366	<	const Float_t deadZoneCenterME1_3[4] = {-83.155,-22.7401,27.86665,81.005};
2367	<	const Float_t deadZoneCenterME2_1[4] = {-95.94,-27.47,33.67,93.72};
2368	<	const Float_t deadZoneCenterME3_1[4] = {-85.97,-36.21,23.68,84.04};
2369	<	const Float_t deadZoneCenterME4_1[4] = {-75.82,-26.14,23.85,73.91};
2370	<	const Float_t deadZoneCenterME234_2[6] = {-162.48,-81.8744,-21.18165,39.51105,100.2939,160.58};
2364	>	//the ME11 wires are not parallel to x, but no gap
2365	>	//chamber edges are not included.
2366	>	const Float_t deadZoneCenterME1_2[2] = {-32.88305,32.867423};
2367	>	const Float_t deadZoneCenterME1_3[2] = {-22.7401,27.86665};
2368	>	const Float_t deadZoneCenterME2_1[2] = {-27.47,33.67};
2369	>	const Float_t deadZoneCenterME3_1[2] = {-36.21,23.68};
2370	>	const Float_t deadZoneCenterME4_1[2] = {-26.14,23.85};
2371	>	const Float_t deadZoneCenterME234_2[4] = {-81.8744,-21.18165,39.51105,100.2939};
2372		const Float_t *deadZoneCenter;
2373	<	Float_t deadZoneHeightHalf=2.53/2.,minY=999999.;
2374	<	Byte_t nGaps=4;
2373	>	Float_t deadZoneHeightHalf=0.32*7/2;// wire spacing * (wires missing + 1)/2
2374	>	Float_t minY=999999.;
2375	>	Byte_t nGaps=2;
2376		switch (abs(StationAndRing)) {
2377		case 11:
2378		case 14:
2379	<	deadZoneCenter=deadZoneCenterME1_14;
2416	<	deadZoneHeightHalf=2/2.;
2417	<	nGaps=2;
2379	>	return 162;//the height of ME11
2380		break;
2381		case 12:
2382		deadZoneCenter=deadZoneCenterME1_2;
#	Line 2424 | Line 2386 | Float_t TPTrackMuonSys::YDistToHVDeadZon
2386		break;
2387		case 21:
2388		deadZoneCenter=deadZoneCenterME2_1;
2427	–	deadZoneHeightHalf=2.5/2.;
2389		break;
2390		case 31:
2391		deadZoneCenter=deadZoneCenterME3_1;
2431	–	deadZoneHeightHalf=2.5/2.;
2392		break;
2393		case 41:
2394		deadZoneCenter=deadZoneCenterME4_1;
2435	–	deadZoneHeightHalf=2.5/2.;
2395		break;
2396		default:
2397		deadZoneCenter=deadZoneCenterME234_2;
2398	<	nGaps=6;
2398	>	nGaps=4;
2399		}
2400		for ( Byte_t iGap=0;iGap<nGaps;iGap++ ) {
2401		Float_t newMinY=yLocal<deadZoneCenter[iGap]?deadZoneCenter[iGap]-deadZoneHeightHalf-yLocal:yLocal-(deadZoneCenter[iGap]+deadZoneHeightHalf);

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines