ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/Jeng/PVStudy/plugins/PVStudy.cc

Comparing UserCode/Jeng/PVStudy/plugins/PVStudy.cc (file contents):
Revision 1.9 by yygao, Thu Oct 29 02:05:54 2009 UTC vs.
Revision 1.11 by yygao, Wed Nov 11 19:28:50 2009 UTC

#	Line 72 | Line 72 | PVStudy::PVStudy(const edm::ParameterSet
72		splitTrackCollection2Tag_     = iConfig.getUntrackedParameter<edm::InputTag>("splitTrackCollection2");
73		vertexCollectionTag_       = iConfig.getUntrackedParameter<edm::InputTag>("vertexCollection");
74		splitVertexCollection1Tag_       = iConfig.getUntrackedParameter<edm::InputTag>("splitVertexCollection1");
75	<	splitVertexCollection2Tag_       = iConfig.getUntrackedParameter<edm::InputTag>("splitVertexCollection2");
75	>	splitVertexCollection2Tag_       = iConfig.getUntrackedParameter<edm::InputTag>("splitVertexCollection2");
76	>	pixelVertexCollectionTag_      = iConfig.getParameter<edm::InputTag>("pixelVertexCollectionTag");
77		verbose_         = iConfig.getUntrackedParameter<bool>("verbose",false);
78		realData_        = iConfig.getUntrackedParameter<bool>("realData",false);
79		analyze_         = iConfig.getUntrackedParameter<bool>("analyzeOnTheFly",false);
#	Line 110 | Line 111 | PVStudy::PVStudy(const edm::ParameterSet
111		//pvtx using all tracks
112		pvtxtree_->Branch("nrecPV",&nrecPV_,"nrecPV/I");
113		pvtxtree_->Branch("nTrkPV",&nTrkPV_,"nTrkPV[nrecPV]/I");
114	+	pvtxtree_->Branch("sumptsq",&sumptsq_,"sumptsq[nrecPV]/D");
115	+	pvtxtree_->Branch("isValid",&isValid_,"isValid/I");
116	+	pvtxtree_->Branch("isFake",&isFake_,"isFake/I");
117		pvtxtree_->Branch("recx",&recx_,"recx[nrecPV]/D");
118		pvtxtree_->Branch("recy",&recy_,"recy[nrecPV]/D");
119		pvtxtree_->Branch("recz",&recz_,"recz[nrecPV]/D");
120		pvtxtree_->Branch("recx_err",&recx_err_,"recx_err[nrecPV]/D");
121		pvtxtree_->Branch("recy_err",&recy_err_,"recy_err[nrecPV]/D");
122		pvtxtree_->Branch("recz_err",&recz_err_,"recz_err[nrecPV]/D");
123	+
124	+	pvtxtree_->Branch("min_zsep",&min_zsep_,"min_zsep/D");
125	+	pvtxtree_->Branch("min_ntrksep",&min_ntrksep_,"min_ntrksep/D");
126	+	pvtxtree_->Branch("min_sumpt2sep",&min_sumpt2sep_,"min_sumpt2sep/D");
127	+
128		//pvtx using splittracks1
129		pvtxtree_->Branch("nrecPV1_spl",&nrecPV1_spl_,"nrecPV1_spl/I");
130	<	pvtxtree_->Branch("nTrkPV1_spl",&nTrkPV1_spl_,"nTrkPV1_spl[nrecPV1_spl]/I");
130	>	pvtxtree_->Branch("nTrkPV1_spl",&nTrkPV1_spl_,"nTrkPV1_spl[nrecPV1_spl]/I");
131	>	pvtxtree_->Branch("sumptsq1_spl",&sumptsq1_spl_,"sumptsq1_spl[nrecPV1_spl]/D");
132	>	pvtxtree_->Branch("isValid1_spl",&isValid1_spl_,"isValid1_spl/I");
133	>	pvtxtree_->Branch("isFake1_spl",&isFake1_spl_,"isFake1_spl/I");
134		pvtxtree_->Branch("recx1_spl",&recx1_spl_,"recx1_spl[nrecPV1_spl]/D");
135		pvtxtree_->Branch("recy1_spl",&recy1_spl_,"recy1_spl[nrecPV1_spl]/D");
136		pvtxtree_->Branch("recz1_spl",&recz1_spl_,"recz1_spl[nrecPV1_spl]/D");
137		pvtxtree_->Branch("recx1_err_spl",&recx1_err_spl_,"recx1_err_spl[nrecPV1_spl]/D");
138		pvtxtree_->Branch("recy1_err_spl",&recy1_err_spl_,"recy1_err_spl[nrecPV1_spl]/D");
139	<	pvtxtree_->Branch("recz1_err_spl",&recz1_err_spl_,"recz1_err_spl[nrecPV1_spl]/D");
139	>	pvtxtree_->Branch("recz1_err_spl",&recz1_err_spl_,"recz1_err_spl[nrecPV1_spl]/D");
140	>
141		//pvtx using splittracks2
142		pvtxtree_->Branch("nrecPV2_spl",&nrecPV2_spl_,"nrecPV2_spl/I");
143	<	pvtxtree_->Branch("nTrkPV2_spl",&nTrkPV2_spl_,"nTrkPV2_spl[nrecPV2_spl]/I");
143	>	pvtxtree_->Branch("nTrkPV2_spl",&nTrkPV2_spl_,"nTrkPV2_spl[nrecPV2_spl]/I");
144	>	pvtxtree_->Branch("sumptsq2_spl",&sumptsq2_spl_,"sumptsq2_spl[nrecPV2_spl]/D");
145	>	pvtxtree_->Branch("isValid2_spl",&isValid2_spl_,"isValid2_spl/I");
146	>	pvtxtree_->Branch("isFake2_spl",&isFake2_spl_,"isFake2_spl/I");
147		pvtxtree_->Branch("recx2_spl",&recx2_spl_,"recx2_spl[nrecPV2_spl]/D");
148		pvtxtree_->Branch("recy2_spl",&recy2_spl_,"recy2_spl[nrecPV2_spl]/D");
149		pvtxtree_->Branch("recz2_spl",&recz2_spl_,"recz2_spl[nrecPV2_spl]/D");
150		pvtxtree_->Branch("recx2_err_spl",&recx2_err_spl_,"recx2_err_spl[nrecPV2_spl]/D");
151		pvtxtree_->Branch("recy2_err_spl",&recy2_err_spl_,"recy2_err_spl[nrecPV2_spl]/D");
152	<	pvtxtree_->Branch("recz2_err_spl",&recz2_err_spl_,"recz2_err_spl[nrecPV2_spl]/D");
152	>	pvtxtree_->Branch("recz2_err_spl",&recz2_err_spl_,"recz2_err_spl[nrecPV2_spl]/D");
153	>
154	>	//pixeVertices
155	>	pvtxtree_->Branch("nrecPV_pxlpvtx",&nrecPV_pxlpvtx_,"nrecPV_pxlpvtx/I");
156	>	pvtxtree_->Branch("nTrkPV_pxlpvtx",&nTrkPV_pxlpvtx_,"nTrkPV_pxlpvtx[nrecPV_pxlpvtx]/I");
157	>	pvtxtree_->Branch("sumptsq_pxlpvtx",&sumptsq_pxlpvtx_,"sumptsq_pxlpvtx[nrecPV_pxlpvtx]/D");
158	>	pvtxtree_->Branch("isValid_pxlpvtx",&isValid_pxlpvtx_,"isValid_pxlpvtx/I");
159	>	pvtxtree_->Branch("isFake_pxlpvtx",&isFake_pxlpvtx_,"isFake_pxlpvtx/I");
160	>	pvtxtree_->Branch("recx_pxlpvtx",&recx_pxlpvtx_,"recx_pxlpvtx[nrecPV_pxlpvtx]/D");
161	>	pvtxtree_->Branch("recy_pxlpvtx",&recy_pxlpvtx_,"recy_pxlpvtx[nrecPV_pxlpvtx]/D");
162	>	pvtxtree_->Branch("recz_pxlpvtx",&recz_pxlpvtx_,"recz_pxlpvtx[nrecPV_pxlpvtx]/D");
163	>	pvtxtree_->Branch("recx_err_pxlpvtx",&recx_err_pxlpvtx_,"recx_err_pxlpvtx[nrecPV_pxlpvtx]/D");
164	>	pvtxtree_->Branch("recy_err_pxlpvtx",&recy_err_pxlpvtx_,"recy_err_pxlpvtx[nrecPV_pxlpvtx]/D");
165	>	pvtxtree_->Branch("recz_err_pxlpvtx",&recz_err_pxlpvtx_,"recz_err_pxlpvtx[nrecPV_pxlpvtx]/D");
166	>
167		//Fill the variables in the twovtx pair (recvtx1, recvtx2)
168		pvtxtree_->Branch("nrecPV_twovtx",&nrecPV_twovtx_,"nrecPV_twovtx/I");
169		pvtxtree_->Branch("nTrkPV1_spl_twovtx",&nTrkPV1_spl_twovtx_,"nTrkPV1_spl_twovtx[nrecPV_twovtx]/I");
#	Line 198 | Line 228 | PVStudy::PVStudy(const edm::ParameterSet
228		//i=0 : x (as in unsplit track collection)
229		//i=1 : x1_spl_
230		//i=2 : x2_spl_
231	<
231	>
232	>
233		for(int i=0;i<3;i++)
234		{
235		string suffix;
#	Line 212 | Line 243 | PVStudy::PVStudy(const edm::ParameterSet
243		h["trkPt"+suffix]  = subDir.make<TH1D>(TString("trkPt"+suffix), TString("Pt of rec tracks "+suffix),100,0,100);
244		h["trkEta"+suffix] = subDir.make<TH1D>(TString("trkEta"+suffix), TString("#eta of rec tracks "+suffix),100,-3,3);
245		h["trkPhi"+suffix] = subDir.make<TH1D>(TString("trkPhi"+suffix), TString("#Phi of rec tracks "+suffix),100,-3.2,3.2);
246	<	h["trkDxy"+suffix] = subDir.make<TH1D>(TString("trkDxy"+suffix), TString("Dxy of rec tracks "+suffix),100,-5,5);
246	>	h["trkDxy"+suffix] = subDir.make<TH1D>(TString("trkDxy"+suffix), TString("Dxy of rec tracks "+suffix),100,-0.5,0.5);
247		h["trkDz"+suffix] = subDir.make<TH1D>(TString("trkDz"+suffix), TString("Dz of rec tracks "+suffix),100,-50,50);
248
249		h["nTrkPV"+suffix]   = subDir.make<TH1D>(TString("nTrkPV"+suffix), TString("Num of rec tracks in PV"+suffix),300,0,300);
#	Line 227 | Line 258 | PVStudy::PVStudy(const edm::ParameterSet
258		h["nrecPVDiff"]     = subDir.make<TH1D>("nrecPVDiff","nrecPV1-nRecPV2",21,-10.5,10.5);
259		h["nTrkPVDiff"]     = subDir.make<TH1D>("nTrkPVDiff","nTrkPV1-nTrkPV2",41,-20.5,20.5);
260		h["nTrkPVRelDiff"]  = subDir.make<TH1D>("nTrkPVRelDiff","(nTrkPV1-nTrkPV2)/(nTrkPV1+nTrkPV2)",100,-1,1);
261	<
261	>
262	>	// Histograms on comparing the multi-vertices
263	>	// Difference in reconstructed vtx position
264	>	h["min_xsep"]   = subDir.make<TH1D>("min_xsep", "min x diff of primary and secondary pvtx",300,0,0.1);
265	>	h["min_xsep_sign"]   = subDir.make<TH1D>("min_xsep_sign", "min x diff in signf of primary and secondary pvtx",300,0,5);
266	>	h["min_ysep"]   = subDir.make<TH1D>("min_ysep", "min y diff of primary and secondary pvtx",300,0,0.1);
267	>	h["min_ysep_sign"]   = subDir.make<TH1D>("min_ysep_sign", "min y diff in signf of primary and secondary pvtx",300,0,5);
268	>	h["min_zsep"]   = subDir.make<TH1D>("min_zsep", "min z diff of primary and secondary pvtx",300,0,5);
269	>	h["min_zsep_sign"]   = subDir.make<TH1D>("min_zsep_sign", "min z diff in signf of primary and secondary pvtx",300,0,200);
270	>	// Difference in reconstructed vtx position
271	>	h["min_ntrksep"]   = subDir.make<TH1D>("min_ntrksep", "min nTrk diff of primary and secondary pvtx",201,-50.5,150.5);
272	>	h["min_sumpt2sep"]   = subDir.make<TH1D>("min_sumpt2sep", "min sumpt2 diff of primary and secondary pvtx",300,0,10000);
273	>
274		//Book histograms sensitive to data/mc
275		for (vector<int>::const_iterator it= datamode.begin(); it != datamode.end() ; ++it) {
276		string suffix;
#	Line 306 | Line 349 | PVStudy::PVStudy(const edm::ParameterSet
349		h2["pullz"+suffix+"_nTrk"]->GetXaxis()->SetTitle("Num of tracks");
350		} // End of   if(analyze_) {
351		suffix.clear();
352	<	} // End of Book histograms sensitive to data/mc
352	>	} // End of Book histograms sensitive to data/mc
353	>
354	>	// Book histograms about pixelVertices
355	>	h["trkdz_pxlpvtxdz"]   = subDir.make<TH1D>("trkdz_pxlpvtxdz", "(Track dz - pixelpvtx dz) in cm",300,-0.5,0.5);
356	>	h["trkdz_pxlpvtxdz_pxlpvtxdzerr"]   = subDir.make<TH1D>("trkdz_pxlpvtxdz_pxlpvtxdzerr", "|Track dz - pixelpvtx dz| / pxlpvtxdzErr",300,0,100);
357	>	h["trkdz_pxlpvtxdz_trkdzerr"]   = subDir.make<TH1D>("trkdz_pxlpvtxdz_trkdzerr", "|Track dz - pixelpvtx dz| / trkdzErr",300,0,50);
358	>	h["trkdzErr_pxlpvtx"]   = subDir.make<TH1D>("trkdzErr_pxlpvtxdz", "Track dzErr of leading pixelpvtx ",300,0,0.5);
359	>	h["trkdzErr_pvtx"]   = subDir.make<TH1D>("trkdzErr_pvtx", "Track dzErr of the leading pvtx ",300,0,0.5);
360	>	h["dzErr_pxlpvtx"]   = subDir.make<TH1D>("dzErr_pxlpvtx", "zError of the leading pvtx ",300,0,0.5);
361	>	// Compare offlinePrimaryVertices with pixelVertices
362	>	h["nrecPV_minus_nrecPxlPV"] =  subDir.make<TH1D>("nrecPV_minus_nrecPxlPV", "nrecPV_minus_nrecPxlPV",21,-10.5,10.5);
363	>	h["recxPV_minus_recxPxlPV"] =  subDir.make<TH1D>("recxPV_minus_recxPxlPV", "recxPV_minus_recxPxlPV",300,-0.02,0.02);
364	>	h["recyPV_minus_recyPxlPV"] =  subDir.make<TH1D>("recyPV_minus_recyPxlPV", "recyPV_minus_recyPxlPV",300,-0.02,0.02);
365	>	h["reczPV_minus_reczPxlPV"] =  subDir.make<TH1D>("reczPV_minus_reczPxlPV", "reczPV_minus_reczPxlPV",300,-0.1,0.1);
366	>
367		// Book MC only plots
368		if (!realData_) {
369		h["genPart_T"]      = subDir.make<TH1D>("genPart_T","t component of gen particles",300,-0.5,0.5);
#	Line 555 | Line 612 | PVStudy::analyze(const edm::Event& iEven
612		} else {
613		cout << "trackCollection cannot be found -> using empty collection of same type." <<endl;
614		}
615	+
616		//splitTrackCollection1
617		static const reco::TrackCollection s_empty_splitTrackColl1;
618		const reco::TrackCollection *splitTrackColl1 = &s_empty_splitTrackColl1;
#	Line 565 | Line 623 | PVStudy::analyze(const edm::Event& iEven
623		} else {
624		cout << "splitTrackCollection1 cannot be found -> using empty collection of same type." <<endl;
625		}
626	+
627		//splitTrackCollection2
628		static const reco::TrackCollection s_empty_splitTrackColl2;
629		const reco::TrackCollection *splitTrackColl2 = &s_empty_splitTrackColl2;
#	Line 575 | Line 634 | PVStudy::analyze(const edm::Event& iEven
634		} else {
635		cout << "splitTrackCollection2 cannot be found -> using empty collection of same type." <<endl;
636		}
637	<
637	>
638	>
639	>	//=======================================================
640	>	// Fill trackparameters of the input tracks to pvtx fitter
641	>	//=======================================================
642	>	if(verbose_)
643	>	cout<<"Start filling track parameters of the input tracks to pvtx fitter."<<endl;
644	>	//fillTrackHisto(const reco::TrackCollection *trackColl, int datatype)
645	>	// datatype: unsplittracks (0); splittracks1 (1);  splittracks2 (2);
646	>	fillTrackHisto(trackColl, 0);
647	>	fillTrackHisto(splitTrackColl1, 1);
648	>	fillTrackHisto(splitTrackColl2, 2);
649	>	if(verbose_)
650	>	cout<<"End filling track parameters of the input tracks to pvtx fitter."<<endl;
651	>
652		//=======================================================
653		// PVTX accessors
654		//=======================================================
#	Line 589 | Line 662 | PVStudy::analyze(const edm::Event& iEven
662		} else {
663		cout << "vertexCollection cannot be found -> using empty collection of same type." <<endl;
664		}
665	+
666		//splitVertexCollection1
667		static const reco::VertexCollection s_empty_splitVertexColl1;
668		const reco::VertexCollection *splitVertexColl1 = &s_empty_splitVertexColl1;
#	Line 598 | Line 672 | PVStudy::analyze(const edm::Event& iEven
672		splitVertexColl1 = splitVertexCollection1Handle.product();
673		} else {
674		cout << "splitVertexCollection1 cannot be found -> using empty collection of same type." <<endl;
675	<	}
675	>	}
676	>
677		//splitVertexCollection2
678		static const reco::VertexCollection s_empty_splitVertexColl2;
679		const reco::VertexCollection *splitVertexColl2 = &s_empty_splitVertexColl2;
#	Line 609 | Line 684 | PVStudy::analyze(const edm::Event& iEven
684		} else {
685		cout << "splitVertexCollection2 cannot be found -> using empty collection of same type." <<endl;
686		}
687	<
687	>
688		if(verbose_) cout<<"End accessing the track and vertex collections"<<endl;
689
690		//=======================================================
#	Line 631 | Line 706 | PVStudy::analyze(const edm::Event& iEven
706		}catch(const Exception&){
707		cout << "Some problem occurred with the particle data table. This may not work !." <<endl;
708		}
709	+
710	+	//=======================================================
711	+	// GET pixelVertices
712	+	//=======================================================
713	+	static const reco::VertexCollection s_empty_pixelVertexColl;
714	+	const reco::VertexCollection *pixelVertexColl = &s_empty_pixelVertexColl;
715	+	edm::Handle<reco::VertexCollection>  pixelVertexCollectionHandle;
716	+	iEvent.getByLabel(pixelVertexCollectionTag_, pixelVertexCollectionHandle);
717	+	if( iEvent.getByLabel(pixelVertexCollectionTag_, pixelVertexCollectionHandle)) {
718	+	pixelVertexColl = pixelVertexCollectionHandle.product();
719	+	} else {
720	+	cout << "pixelVertexCollection cannot be found. -> using empty collection of same type." <<endl;
721	+	}
722
723		//=======================================================
724	<	// Fill trackparameters of the input tracks to pvtx fitter
724	>	// Fill pixelVertices related histograms
725		//=======================================================
726	<	if(verbose_)
727	<	cout<<"Start filling track parameters of the input tracks to pvtx fitter."<<endl;
728	<	//fillTrackHisto(const reco::TrackCollection *trackColl, int datatype)
729	<	// datatype: unsplittracks (0); splittracks1 (1);  splittracks2 (2);
730	<	fillTrackHisto(trackColl, 0);
731	<	fillTrackHisto(splitTrackColl1, 1);
732	<	fillTrackHisto(splitTrackColl2, 2);
733	<	if(verbose_)
734	<	cout<<"End filling track parameters of the input tracks to pvtx fitter."<<endl;
726	>	nrecPV_pxlpvtx_ = int (pixelVertexColl->size());
727	>	if(pixelVertexColl->size()>0 && pixelVertexColl->begin()->isValid() && !(pixelVertexColl->begin()->isFake())) {
728	>	//fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple) {
729	>	fillTrackHistoInPV(pixelVertexColl, 4, false, true);
730	>	h["dzErr_pxlpvtx"]->Fill( pixelVertexColl->begin()->zError());
731	>	// Get the dZ error of the tracks in the leading pixelVertexColl
732	>	for(reco::Vertex::trackRef_iterator t = (pixelVertexColl->begin())->tracks_begin();
733	>	t!= (pixelVertexColl->begin())->tracks_end(); t++) {
734	>
735	>	if ( (**t).charge() < -1 || (**t).charge() > 1 ) {
736	>	//        h["trkPtPV"]->Fill(0.);
737	>	}
738	>	else
739	>	h["trkdzErr_pxlpvtx"]->Fill((**t).dzError());
740	>	}
741	>	// Fill the track->dz() in the PV difference from first pixelpvtx
742	>	for(reco::Vertex::trackRef_iterator t = (vertexColl->begin())->tracks_begin();
743	>	t!= (vertexColl->begin())->tracks_end(); t++) {
744	>	// illegal charge
745	>	if ( (**t).charge() < -1 || (**t).charge() > 1 ) {
746	>	//        h["trkPtPV"]->Fill(0.);
747	>	}
748	>	else {
749	>	if(pixelVertexColl->size()>0) {
750	>	h["trkdz_pxlpvtxdz"]->Fill((**t).dz() -  pixelVertexColl->begin()->z());
751	>	h["trkdz_pxlpvtxdz_pxlpvtxdzerr"]->Fill(fabs((**t).dz() -  pixelVertexColl->begin()->z())/pixelVertexColl->begin()->zError());
752	>	h["trkdz_pxlpvtxdz_trkdzerr"]->Fill(fabs((**t).dz() -  pixelVertexColl->begin()->z())/(**t).dzError());
753	>	h["trkdzErr_pvtx"]->Fill((**t).dzError());
754	>	}
755	>	}
756	>	}
757	>	}
758
759		//=======================================================
760	<	// Fill reconstructed vertices ootb
760	>	// Fill number of reconstructed vertices
761		//=======================================================
762	+
763		if(verbose_)  {
764		cout<<"PVStudy::analyze() Printing vertexCollection: "<<endl;
765		printRecVtxs(vertexCollectionHandle);
#	Line 657 | Line 769 | PVStudy::analyze(const edm::Event& iEven
769		printRecVtxs(splitVertexCollection2Handle);
770		cout<<"Start filling pvtx parameters."<<endl;
771		}
772	<
772	>
773		nrecPV_ = int (vertexColl->size());
774		nrecPV1_spl_ = int (splitVertexColl1->size());
775		nrecPV2_spl_ = int (splitVertexColl2->size());
#	Line 666 | Line 778 | PVStudy::analyze(const edm::Event& iEven
778		h["nrecPV1_spl"]->Fill(nrecPV1_spl_);
779		h["nrecPV2_spl"]->Fill(nrecPV2_spl_);
780		h["nrecPVDiff"]->Fill(double(nrecPV1_spl_)-double(nrecPV2_spl_));
781	<
781	>
782	>	//=================================================================
783	>	// Fill track parameter ntuple/hist for tracks used in recoVertices
784	>	//=================================================================
785	>
786		//fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple) {
787	<	fillTrackHistoInPV(vertexColl, 0, true, true);
788	<	fillTrackHistoInPV(splitVertexColl1, 1, false, true);
789	<	fillTrackHistoInPV(splitVertexColl2, 2, false, true);
787	>	if(vertexColl->size() > 0 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake()))
788	>	fillTrackHistoInPV(vertexColl, 0, true, true);
789	>
790	>	if(splitVertexColl1->size() > 0 && splitVertexColl1->begin()->isValid() && !(splitVertexColl1->begin()->isFake()))
791	>	fillTrackHistoInPV(splitVertexColl1, 1, false, true);
792	>
793	>	if(splitVertexColl1->size() > 0 && splitVertexColl2->begin()->isValid() && !(splitVertexColl2->begin()->isFake()))
794	>	fillTrackHistoInPV(splitVertexColl2, 2, false, true);
795	>
796	>
797	>	//=======================================================
798	>	// Compare offlinePrimaryVertices with pixelVertices
799	>	//=======================================================
800	>	if( (pixelVertexColl->size()>0 && pixelVertexColl->begin()->isValid() && !(pixelVertexColl->begin()->isFake()))
801	>	&& (vertexColl->size()>0 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake())) ) {
802	>	h["nrecPV_minus_nrecPxlPV"]->Fill( double (nrecPV_ - nrecPV_pxlpvtx_));
803	>	// difference in reconstructed position of the leading pvtx
804	>	//cout<<"recx_[0] = "<< recx_[0] << "recx_pxlpvtx_[0] = "<< recx_pxlpvtx_[0]<<endl;
805	>	//cout<<"recy_[0] = "<< recy_[0] << "recy_pxlpvtx_[0] = "<< recy_pxlpvtx_[0]<<endl;
806	>	h["recxPV_minus_recxPxlPV"]->Fill (recx_[0] - recx_pxlpvtx_[0]);
807	>	h["recyPV_minus_recyPxlPV"]->Fill (recy_[0] - recy_pxlpvtx_[0]);
808	>	h["reczPV_minus_reczPxlPV"]->Fill (recz_[0] - recz_pxlpvtx_[0]);
809	>	}
810	>
811	>
812	>
813	>	//==========================================================
814	>	// Fill secondary/primary min separations for multi-vertices
815	>	//==========================================================
816	>
817	>	if(vertexColl->size()>1 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake()) ) {
818	>
819	>	double min_xsep = 9999.0;
820	>	double min_ysep = 9999.0;
821	>	double min_zsep = 9999.0;
822	>	double min_xsep_sign = 9999.0;
823	>	double min_ysep_sign = 9999.0;
824	>	double min_zsep_sign = 9999.0;
825	>	double min_ntrksep = 9999.0;
826	>	double min_sumpt2sep = 9999999.0;
827	>
828	>	if(verbose_) cout<<"leading pvtx z = "<< vertexColl->begin()->z() <<endl;
829	>
830	>	// Looping through the secondary vertices to find the mininum separation to the primary
831	>	for(reco::VertexCollection::const_iterator v=vertexColl->begin() + 1 ;
832	>	v!=vertexColl->end(); ++v) {
833	>	if(v->isValid() && ! v->isFake()) {
834	>	// xsep
835	>	if(fabs(v->x()- vertexColl->begin()->x())<min_xsep)
836	>	min_xsep = fabs(v->x()- vertexColl->begin()->x());
837	>	// ysep
838	>	if(fabs(v->y()- vertexColl->begin()->y())<min_ysep)
839	>	min_ysep = fabs(v->y()- vertexColl->begin()->y());
840	>	// zsep
841	>	if(fabs(v->z()- vertexColl->begin()->z())<min_zsep)
842	>	min_zsep = fabs(v->z()- vertexColl->begin()->z());
843	>	// xsep_sign
844	>	double xsep_sign = fabs(v->x()- vertexColl->begin()->x())/max(v->xError(), vertexColl->begin()->xError());
845	>	if(xsep_sign < min_xsep_sign)
846	>	min_xsep_sign = xsep_sign;
847	>	// ysep_sign
848	>	double ysep_sign = fabs(v->y()- vertexColl->begin()->y())/max(v->yError(), vertexColl->begin()->yError());
849	>	if(ysep_sign < min_ysep_sign)
850	>	min_ysep_sign = ysep_sign;
851	>	// zsep_sign
852	>	double zsep_sign = fabs(v->z()- vertexColl->begin()->z())/max(v->zError(), vertexColl->begin()->zError());
853	>	if(zsep_sign < min_zsep_sign)
854	>	min_zsep_sign = zsep_sign;
855	>	// ntrksep
856	>	if( (double(vertexColl->begin()->tracksSize()) - double(v->tracksSize())) < min_ntrksep)
857	>	min_ntrksep = double(vertexColl->begin()->tracksSize()) - double(v->tracksSize());
858	>	// sumpt2sep
859	>	if(fabs(sumPtSquared(*v) - sumPtSquared(*(vertexColl->begin()))) < min_sumpt2sep)
860	>	min_sumpt2sep = fabs(sumPtSquared(*v) - sumPtSquared(*(vertexColl->begin())));
861	>	}
862	>	}
863	>	h["min_xsep"]->Fill(min_xsep);
864	>	h["min_ysep"]->Fill(min_ysep);
865	>	h["min_zsep"]->Fill(min_zsep);
866	>	h["min_xsep_sign"]->Fill(min_xsep_sign);
867	>	h["min_ysep_sign"]->Fill(min_ysep_sign);
868	>	h["min_zsep_sign"]->Fill(min_zsep_sign);
869	>	h["min_ntrksep"]->Fill(min_ntrksep);
870	>	h["min_sumpt2sep"]->Fill(min_sumpt2sep);
871	>
872	>	min_zsep_ = min_zsep;
873	>	min_ntrksep_ = min_ntrksep;
874	>	min_sumpt2sep_ = min_sumpt2sep;
875	>
876	>
877	>	} // End of  if(vertexColl->size()>1) {
878	>
879		if(verbose_)
880		cout<<"End filling pvtx parameters."<<endl;
881
#	Line 690 | Line 895 | PVStudy::analyze(const edm::Event& iEven
895		reco::VertexCollection *matchedVertexColl2 = &s_empty_matchedVertexColl2;
896		matchedVertexColl2->reserve(splitVertexColl2->size());
897
898	<	for (size_t ivtx = 0; ivtx<splitVertexCollection1Handle->size(); ++ivtx) {
898	>	bool stopmatching = false;
899	>	for (size_t ivtx = 0; ivtx<splitVertexCollection1Handle->size() && !stopmatching; ++ivtx) {
900		reco::VertexRef recvtx1(splitVertexCollection1Handle, ivtx);
901	<	for (size_t jvtx = ivtx; jvtx < splitVertexCollection2Handle->size(); ++jvtx) {
901	>	if( !(recvtx1->isValid()) || recvtx1->isFake()) break;
902	>	for (size_t jvtx = ivtx; jvtx < splitVertexCollection2Handle->size(); ++jvtx) {
903	>	//------------------------------------------------------------------------
904	>	//== If only considering matching the first vertex of the splitVertexColl
905	>	//------------------------------------------------------------------------
906	>	if (ivtx!=0 || jvtx!=0) { stopmatching = true; break; }
907		reco::VertexRef recvtx2(splitVertexCollection2Handle, jvtx);
908	+	if( !(recvtx2->isValid()) || recvtx2->isFake()) break;
909	+	if(verbose_) {
910	+	cout<<"Matching splitVertexColl1: #  "<< ivtx<< " and splitVertexColl1: # "<<jvtx<<endl;
911	+	cout<<"recvtx1->z() = " << recvtx1->z() << "+/- "<< recvtx1->zError() << "." << endl;
912	+	cout<<"recvtx2->z() = " << recvtx2->z() << "+/- "<< recvtx2->zError() << "." << endl;
913	+	}
914		if(matchVertex(recvtx1, recvtx2, zsigncut_)) {
915		if(verbose_) {
916	<	cout<<"The two splitted vertices match in Z: "<<endl;
700	<	cout<<"recvtx1->z() = " << recvtx1->z() << "+/- "<< recvtx1->zError() << "." << endl;
701	<	cout<<"recvtx2->z() = " << recvtx2->z() << "+/- "<< recvtx2->zError() << "." << endl;
916	>	cout<<"The two splitted vertices match in Z. "<<endl;
917		}
918		int nTrkPV1 = recvtx1->tracksSize();
919		int nTrkPV2 = recvtx2->tracksSize();
#	Line 711 | Line 926 | PVStudy::analyze(const edm::Event& iEven
926		}
927		matchedVertexColl1->push_back(*recvtx1);
928		matchedVertexColl2->push_back(*recvtx2);
929	<	continue;
929	>	stopmatching = true; // stop the matching when the first match is found!
930	>	break;
931		}
932		else
933		cout<<"WARNING: (nTrkPV1-nTrkPV2)/(nTrkPV1+nTrkPV2) = " << ntrkreldiff<<", exceeding cut "<<ntrkdiffcut_<<endl;
#	Line 731 | Line 947 | PVStudy::analyze(const edm::Event& iEven
947		if(verbose_) {
948		cout<<"Begin analyzing the matchedVertexColl with size = "<< matchedVertexColl1->size()<< endl;
949		}
734	–
950
951		// Compare the reconstructed vertex position and calculate resolution/pulls
952		if(matchedVertexColl1->size() != 0 && matchedVertexColl2->size() != 0) {
#	Line 1135 | Line 1350 | void PVStudy::fillTrackHistoInPV(const r
1350		v!=vertexColl->end(); ++v, ++ivtx) {
1351		if(fillNtuple) {
1352		if(datatype == 0) {
1353	<	nTrkPV_[ivtx] = v->tracksSize();
1353	>	nTrkPV_[ivtx] = v->tracksSize();
1354	>	sumptsq_[ivtx] = sumPtSquared(*v);
1355	>	isValid_[ivtx] = v->isValid();
1356	>	isFake_[ivtx] = v->isFake();
1357		recx_[ivtx] = v->x();
1358		recy_[ivtx] = v->y();
1359		recz_[ivtx] = v->z();
1360		recx_err_[ivtx] = v->xError();
1361		recy_err_[ivtx] = v->yError();
1362		recz_err_[ivtx] = v->zError();
1363	+
1364	+
1365		}
1366		if(datatype == 1) {
1367	<	nTrkPV1_spl_[ivtx] = v->tracksSize();
1367	>	nTrkPV1_spl_[ivtx] = v->tracksSize();
1368	>	sumptsq1_spl_[ivtx] = sumPtSquared(*v);
1369	>	isValid1_spl_[ivtx] = v->isValid();
1370	>	isFake1_spl_[ivtx] = v->isFake();
1371		recx1_spl_[ivtx] = v->x();
1372		recy1_spl_[ivtx] = v->y();
1373		recz1_spl_[ivtx] = v->z();
1374		recx1_err_spl_[ivtx] = v->xError();
1375		recy1_err_spl_[ivtx] = v->yError();
1376		recz1_err_spl_[ivtx] = v->zError();
1377	+
1378		}
1379		if(datatype == 2) {
1380	<	nTrkPV2_spl_[ivtx] = v->tracksSize();
1380	>	nTrkPV2_spl_[ivtx] = v->tracksSize();
1381	>	sumptsq2_spl_[ivtx] = sumPtSquared(*v);
1382	>	isValid2_spl_[ivtx] = v->isValid();
1383	>	isFake2_spl_[ivtx] = v->isFake();
1384		recx2_spl_[ivtx] = v->x();
1385		recy2_spl_[ivtx] = v->y();
1386		recz2_spl_[ivtx] = v->z();
#	Line 1161 | Line 1388 | void PVStudy::fillTrackHistoInPV(const r
1388		recy2_err_spl_[ivtx] = v->yError();
1389		recz2_err_spl_[ivtx] = v->zError();
1390		}
1391	+	if(datatype == 4) { // for pixelVertices
1392	+	nTrkPV_pxlpvtx_[ivtx] = v->tracksSize();
1393	+	sumptsq_pxlpvtx_[ivtx] = sumPtSquared(*v);
1394	+	isValid_pxlpvtx_[ivtx] = v->isValid();
1395	+	isFake_pxlpvtx_[ivtx] = v->isFake();
1396	+	recx_pxlpvtx_[ivtx] = v->x();
1397	+	recy_pxlpvtx_[ivtx] = v->y();
1398	+	recz_pxlpvtx_[ivtx] = v->z();
1399	+	recx_err_pxlpvtx_[ivtx] = v->xError();
1400	+	recy_err_pxlpvtx_[ivtx] = v->yError();
1401	+	recz_err_pxlpvtx_[ivtx] = v->zError();
1402	+	}
1403		}
1404	<
1405	<	if(fillHisto) {
1406	<	h["nTrkPV"+suffix]->Fill(v->tracksSize());
1407	<	try {
1408	<	for(reco::Vertex::trackRef_iterator t = v->tracks_begin();
1409	<	t!=v->tracks_end(); t++) {
1410	<	// illegal charge
1411	<	if ( (**t).charge() < -1 || (**t).charge() > 1 ) {
1412	<	//    h["trkPtPV"]->Fill(0.);
1413	<	}
1414	<	else {
1415	<	h["trkPtPV"+suffix]->Fill((**t).pt());
1416	<	h["trkDxyPV"+suffix]->Fill((**t).dxy());
1417	<	h["trkDzPV"+suffix]->Fill((**t).dz());
1418	<	h["trkEtaPV"+suffix]->Fill((**t).eta());
1419	<	h["trkPhiPV"+suffix]->Fill((**t).phi());
1420	<	}
1404	>	}
1405	>	// For histogram only fill the leading pvtx parameters
1406	>	if(fillHisto) {
1407	>	h["nTrkPV"+suffix]->Fill(vertexColl->begin()->tracksSize());
1408	>	try {
1409	>	for(reco::Vertex::trackRef_iterator t = vertexColl->begin()->tracks_begin();
1410	>	t!=vertexColl->begin()->tracks_end(); t++) {
1411	>	// illegal charge
1412	>	if ( (**t).charge() < -1 || (**t).charge() > 1 ) {
1413	>	//      h["trkPtPV"]->Fill(0.);
1414	>	}
1415	>	else {
1416	>	h["trkPtPV"+suffix]->Fill((**t).pt());
1417	>	h["trkDxyPV"+suffix]->Fill((**t).dxy());
1418	>	h["trkDzPV"+suffix]->Fill((**t).dz());
1419	>	h["trkEtaPV"+suffix]->Fill((**t).eta());
1420	>	h["trkPhiPV"+suffix]->Fill((**t).phi());
1421		}
1183	–	}
1184	–	catch (...) {
1185	–	// exception thrown when trying to use linked track
1186	–	//        h["trkPtPV"]->Fill(0.);
1422		}
1423		}
1424	+	catch (...) {
1425	+	// exception thrown when trying to use linked track
1426	+	//          h["trkPtPV"]->Fill(0.);
1427	+	}
1428		}
1429	+
1430		}
1431
1432		void PVStudy::fillMCHisto(std::vector<simPrimaryVertex>::iterator vsim, int isimpv, const reco::VertexCollection *vtxColl, int datatype)
#	Line 1310 | Line 1550 | void PVStudy::SetVarToZero() {
1550		}
1551
1552		//Yanyan's variables
1553	<	// Number of reconstructed vertices  // for GeneralTracks
1553	>	// Number of reconstructed vertices
1554		nrecPV_ = 0;
1555		nrecPV1_spl_ = 0;
1556		nrecPV2_spl_ = 0;
1557		nrecPV_twovtx_ = 0;
1558		nsimPV_ = 0;
1559	+
1560		nrecPV_mct_ = 0;
1561		nrecPV_spl1_mct_ = 0;
1562		nrecPV_spl2_mct_ = 0;
1563
1564	+	// Mininum separation between the secondary pvtxes and leading pvtx
1565	+	min_zsep_ = 9999.0;
1566	+	min_ntrksep_ = 9999.0;
1567	+	min_sumpt2sep_ = -9999.0;
1568	+
1569	+
1570		for (int i = 0; i < nMaxPVs_; i++) {
1571	<	nTrkPV_[i] = 0; // Number of tracks in the pvtx
1571	>	// recoVertices with all tracks
1572	>	nTrkPV_[i] = 0; // Number of tracks in the pvtx
1573	>	sumptsq_[i] = 0;
1574	>	isValid_[i] = -1;
1575	>	isFake_[i] = -1;
1576		recx_[i] = 0;
1577		recy_[i] = 0;
1578		recz_[i] = 0;
1579		recx_err_[i] = 0;
1580		recy_err_[i] = 0;
1581		recz_err_[i] = 0;
1582	<
1583	<	nTrkPV1_spl_[i] = 0; // Number of tracks in the pvtx
1582	>
1583	>	// recoVertices with splitTrack1
1584	>	nTrkPV1_spl_[i] = 0; // Number of tracks in the pvtx
1585	>	sumptsq1_spl_[i] = 0;
1586	>	isValid1_spl_[i] = -1;
1587	>	isFake1_spl_[i] = -1;
1588		recx1_spl_[i] = 0;
1589		recy1_spl_[i] = 0;
1590		recz1_spl_[i] = 0;
1591		recx1_err_spl_[i] = 0;
1592		recy1_err_spl_[i] = 0;
1593		recz1_err_spl_[i] = 0;
1339	–
1340	–	nTrkPV2_spl_[i] = 0; // Number of tracks in the pvtx
1594
1595	+	// recoVertices with splitTrack2
1596	+	nTrkPV2_spl_[i] = 0; // Number of tracks in the pvtx
1597	+	sumptsq2_spl_[i] = 0;
1598	+	isValid2_spl_[i] = -1;
1599	+	isFake2_spl_[i] = -1;
1600		recx2_spl_[i] = 0;
1601		recy2_spl_[i] = 0;
1602		recz2_spl_[i] = 0;
#	Line 1346 | Line 1604 | void PVStudy::SetVarToZero() {
1604		recy2_err_spl_[i] = 0;
1605		recz2_err_spl_[i] = 0;
1606
1607	+	//pixelVertices
1608	+	nTrkPV_pxlpvtx_[i] = 0; // Number of tracks in the pvtx
1609	+	sumptsq_pxlpvtx_[i] = 0;
1610	+	isValid_pxlpvtx_[i] = -1;
1611	+	isFake_pxlpvtx_[i] = -1;
1612	+	recx_pxlpvtx_[i] = 0;
1613	+	recy_pxlpvtx_[i] = 0;
1614	+	recz_pxlpvtx_[i] = 0;
1615	+	recx_err_pxlpvtx_[i] = 0;
1616	+	recy_err_pxlpvtx_[i] = 0;
1617	+	recz_err_pxlpvtx_[i] = 0;
1618	+
1619	+	// matched two-vertices
1620		nTrkPV1_spl_twovtx_[i] = 0;
1621		nTrkPV2_spl_twovtx_[i] = 0;
1622		nTrkPV_twovtx_[i] = 0;
#	Line 1361 | Line 1632 | void PVStudy::SetVarToZero() {
1632		pully_twovtx_[i] = 0;
1633		pullz_twovtx_[i] = 0;
1634
1635	+	//simpv
1636		nsimTrkPV_[i] = 0;
1637		simx_[i] = 0;
1638		simy_[i] = 0;
1639		simz_[i] = 0;
1640		simptsq_[i] = 0;
1641
1642	+	// residual and pulls with mc truth required
1643		deltax_mct_[i] = 0;
1644		deltay_mct_[i] = 0;
1645		deltaz_mct_[i] = 0;
#	Line 1400 | Line 1673 | double PVStudy::sumPtSquared(const reco:
1673		}
1674		return sum;
1675		}
1676	+

Diff Legend

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