CMS 3D CMS Logo

BeamSpotAnalyzer Class Reference

_________________________________________________________________ class: BeamSpotAnalyzer.h package: RecoVertex/BeamSpotProducer More...

#include <RecoVertex/BeamSpotProducer/interface/BeamSpotAnalyzer.h>

Inheritance diagram for BeamSpotAnalyzer:

edm::EDAnalyzer

List of all members.

Public Member Functions

 BeamSpotAnalyzer (const edm::ParameterSet &)
 _________________________________________________________________ class: BeamSpotAnalyzer.cc package: RecoVertex/BeamSpotProducer
 ~BeamSpotAnalyzer ()

Private Member Functions

virtual void analyze (const edm::Event &, const edm::EventSetup &)
virtual void beginJob (const edm::EventSetup &)
virtual void endJob ()

Private Attributes

std::string ckfSeedProducerLabel_
std::string ckfTrackCandidateProducerLabel_
std::string ckfTrackProducerLabel_
std::ofstream fasciiFile
std::string fasciiFileName
std::string fasciifilename
std::vector< BSTrkParametersfBSvector
int fcharge
double fchi2
double fcov [7][7]
double fd0
double fd0phi_chi2
double fd0phi_d0
double feta
TFile * file_
int fmaxNtracks
double fndof
unsigned int fnHit
unsigned int fnPixelHit
unsigned int fnPXBHit
unsigned int fnPXFHit
unsigned int fnStripHit
unsigned int fnTECHit
unsigned int fnTIBHit
unsigned int fnTIDHit
unsigned int fnTOBHit
double fphi0
double fpt
float fptmin
double fsigmad0
double fsigmaz0
double ftheta
int ftotal_tracks
int ftotalevents
TTree * ftree_
double fz0
double inputBeamWidth_
std::string outputfilename_
bool runallfitters_
unsigned int sameNumberOfTracks
bool write2DB_


Detailed Description

_________________________________________________________________ class: BeamSpotAnalyzer.h package: RecoVertex/BeamSpotProducer

author: Francisco Yumiceva, Fermilab (yumiceva@fnal.gov)

version

Id
BeamSpotAnalyzer.h,v 1.5 2008/05/14 15:31:57 yumiceva Exp

________________________________________________________________

Definition at line 32 of file BeamSpotAnalyzer.h.


Constructor & Destructor Documentation

BeamSpotAnalyzer::BeamSpotAnalyzer ( const edm::ParameterSet iConfig  )  [explicit]

_________________________________________________________________ class: BeamSpotAnalyzer.cc package: RecoVertex/BeamSpotProducer

author: Francisco Yumiceva, Fermilab (yumiceva@fnal.gov)

version

Id
BeamSpotAnalyzer.cc,v 1.4 2008/05/14 15:31:57 yumiceva Exp

________________________________________________________________

Definition at line 41 of file BeamSpotAnalyzer.cc.

References ckfTrackProducerLabel_, fasciiFile, fasciiFileName, fBSvector, fcharge, fchi2, fcov, fd0, feta, file_, fmaxNtracks, fndof, fnHit, fnPixelHit, fnPXBHit, fnPXFHit, fnStripHit, fnTECHit, fnTIBHit, fnTIDHit, fnTOBHit, fphi0, fpt, fptmin, fsigmad0, fsigmaz0, ftheta, ftotal_tracks, ftotalevents, ftree_, fz0, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), inputBeamWidth_, outputfilename_, runallfitters_, and write2DB_.

00042 {
00043 
00044         outputfilename_ = iConfig.getUntrackedParameter<std::string>("OutputFileName");
00045         
00046   file_ = TFile::Open(outputfilename_.c_str(),"RECREATE");
00047 
00048   ftree_ = new TTree("mytree","mytree");
00049   ftree_->AutoSave();
00050   
00051   ftree_->Branch("pt",&fpt,"fpt/D");
00052   ftree_->Branch("d0",&fd0,"fd0/D");
00053   ftree_->Branch("sigmad0",&fsigmad0,"fsigmad0/D");
00054   ftree_->Branch("phi0",&fphi0,"fphi0/D");
00055   ftree_->Branch("z0",&fz0,"fz0/D");
00056   ftree_->Branch("sigmaz0",&fsigmaz0,"fsigmaz0/D");
00057   ftree_->Branch("theta",&ftheta,"ftheta/D");
00058   ftree_->Branch("eta",&feta,"feta/D");
00059   ftree_->Branch("charge",&fcharge,"fcharge/I");
00060   ftree_->Branch("chi2",&fchi2,"fchi2/D");
00061   ftree_->Branch("ndof",&fndof,"fndof/D");
00062   ftree_->Branch("nHit",&fnHit,"fnHit/i");
00063   ftree_->Branch("nStripHit",&fnStripHit,"fnStripHit/i");
00064   ftree_->Branch("nPixelHit",&fnPixelHit,"fnPixelHit/i");
00065   ftree_->Branch("nTIBHit",&fnTIBHit,"fnTIBHit/i");
00066   ftree_->Branch("nTOBHit",&fnTOBHit,"fnTOBHit/i");
00067   ftree_->Branch("nTIDHit",&fnTIDHit,"fnTIDHit/i");
00068   ftree_->Branch("nTECHit",&fnTECHit,"fnTECHit/i");
00069   ftree_->Branch("nPXBHit",&fnPXBHit,"fnPXBHit/i");
00070   ftree_->Branch("nPXFHit",&fnPXFHit,"fnPXFHit/i");
00071   ftree_->Branch("cov",&fcov,"fcov[7][7]/D");
00072   
00073    
00074   fBSvector.clear();
00075 
00076   //dump to file
00077   fasciiFileName = outputfilename_.replace(outputfilename_.size()-4,outputfilename_.size(),"txt");
00078   fasciiFile.open(fasciiFileName.c_str());
00079 
00080   
00081   // get parameter
00082  
00083   //ckfSeedProducerLabel_ = iConfig.getUntrackedParameter<std::string>("ckfSeedProducerLabel");
00084   //ckfTrackCandidateProducerLabel_ = iConfig.getUntrackedParameter<std::string>("ckfTrackCandidateProducerLabel");
00085   //ckfTrackProducerLabel_ = iConfig.getUntrackedParameter<std::string>("ckfTrackProducerLabel");
00086 
00087   //sameNumberOfTracks = iConfig.getUntrackedParameter<unsigned int>("sameNumberOfTracks");
00088 
00089   //
00090   fptmin = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<double>("MinimumPt");
00091   fmaxNtracks = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<int>("MaximumNtracks");
00092   ckfTrackProducerLabel_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getUntrackedParameter<std::string>("TrackCollection");
00093    
00094   write2DB_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<bool>("WriteToDB");
00095   runallfitters_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getParameter<bool>("RunAllFitters");
00096   inputBeamWidth_ = iConfig.getParameter<edm::ParameterSet>("BSAnalyzerParameters").getUntrackedParameter<double>("InputBeamWidth",-1.);
00097   
00098   ftotal_tracks = 0;
00099   ftotalevents = 0;
00100   
00101 }

BeamSpotAnalyzer::~BeamSpotAnalyzer (  ) 

Definition at line 104 of file BeamSpotAnalyzer.cc.

00105 {
00106         
00107 }


Member Function Documentation

void BeamSpotAnalyzer::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
) [private, virtual]

Implements edm::EDAnalyzer.

Definition at line 111 of file BeamSpotAnalyzer.cc.

References funct::abs(), ckfTrackProducerLabel_, fBSvector, fcharge, fchi2, fcov, fd0, feta, fndof, fnHit, fnPixelHit, fnPXBHit, fnPXFHit, fnStripHit, fnTECHit, fnTIBHit, fnTIDHit, fnTOBHit, fphi0, fpt, fsigmad0, fsigmaz0, ftheta, ftotal_tracks, ftotalevents, ftree_, fz0, edm::Event::getByLabel(), i, id, j, PixelSubdetector::PixelBarrel, PixelSubdetector::PixelEndcap, edm::Handle< T >::product(), StripSubdetector::TEC, StripSubdetector::TIB, StripSubdetector::TID, StripSubdetector::TOB, and track.

00112 {
00113         
00114         ftree_->SetBranchAddress("theta",&ftheta);
00115         ftree_->SetBranchAddress("pt",&fpt);
00116         ftree_->SetBranchAddress("eta",&feta);
00117         ftree_->SetBranchAddress("charge",&fcharge);
00118         ftree_->SetBranchAddress("chi2",&fchi2);
00119         ftree_->SetBranchAddress("ndof",&fndof);
00120         ftree_->SetBranchAddress("d0",&fd0);
00121         ftree_->SetBranchAddress("sigmad0",&fsigmad0);
00122         ftree_->SetBranchAddress("phi0",&fphi0);
00123         ftree_->SetBranchAddress("z0",&fz0);
00124         ftree_->SetBranchAddress("sigmaz0",&fsigmaz0);
00125         ftree_->SetBranchAddress("nHit",&fnHit);
00126         ftree_->SetBranchAddress("nStripHit",&fnStripHit);
00127         ftree_->SetBranchAddress("nPixelHit",&fnPixelHit);
00128         ftree_->SetBranchAddress("nTIBHit",&fnTIBHit);
00129         ftree_->SetBranchAddress("nTOBHit",&fnTOBHit);
00130         ftree_->SetBranchAddress("nTIDHit",&fnTIDHit);
00131         ftree_->SetBranchAddress("nTECHit",&fnTECHit);
00132         ftree_->SetBranchAddress("nPXBHit",&fnPXBHit);
00133         ftree_->SetBranchAddress("nPXFHit",&fnPXFHit);
00134         ftree_->SetBranchAddress("cov",&fcov);
00135         
00136   
00137         // get collections
00138         
00139         //edm::Handle<TrackCandidateCollection> ckfTrackCandidateCollectionHandle;
00140         //iEvent.getByLabel(ckfTrackCandidateProducerLabel_,ckfTrackCandidateCollectionHandle);
00141         //const TrackCandidateCollection *ckfTrackCandidateCollection = ckfTrackCandidateCollectionHandle.product();
00142         
00143         edm::Handle<reco::TrackCollection> ckfTrackCollectionHandle;
00144         //iEvent.getByLabel(TkTag,tkCollection);
00145         iEvent.getByLabel(ckfTrackProducerLabel_,ckfTrackCollectionHandle);
00146         const reco::TrackCollection *ckfTrackCollection = ckfTrackCollectionHandle.product();
00147 
00148 
00149   // Ckf tracks
00150   //for (unsigned int itrack=0;itrack<tkCollection->size();++itrack){
00151           //build the ref to the track
00152           //reco::TrackRef track=reco::TrackRef(ckfTrackCollectionHandly,itrack);
00153 
00154         //std::cout << "track collection size: "<< ckfTrackCollection->size() << std::endl;
00155         
00156         for ( reco::TrackCollection::const_iterator track = ckfTrackCollection->begin();
00157                   track != ckfTrackCollection->end();
00158                   ++track ) {
00159 
00160                 
00161           fpt = track->pt();
00162           //std::cout << "pt= "<< track->pt() << std::endl;
00163           //std::cout << "eta= "<< track->eta() << std::endl;
00164           
00165           feta = track->eta();
00166           fphi0 = track->momentum().phi();
00167           fcharge = track->charge();
00168           fchi2 = track->chi2();
00169           fndof = track->ndof();
00170           
00171           //fsigmaphi0 = track->phi0Error();
00172           fd0 = track->d0();
00173           fsigmad0 = track->d0Error();
00174           fz0 = track->dz();
00175           fsigmaz0 = track->dzError();
00176           ftheta = track->theta();
00177 
00178           for (int i=0; i<5; ++i) {
00179                   for (int j=0; j<5; ++j) {
00180                           fcov[i][j] = track->covariance(i,j);
00181                   }
00182           }
00183           
00184           // loop over hits in tracks, count
00185           fnHit      = 0;
00186           fnStripHit = 0;
00187           fnPixelHit = 0;
00188           fnTIBHit   = 0;
00189           fnTOBHit   = 0;
00190           fnTIDHit   = 0;
00191           fnTECHit   = 0;
00192           fnPXBHit   = 0;
00193           fnPXFHit   = 0;
00194 
00195           for ( trackingRecHit_iterator recHit = track->recHitsBegin();
00196                         recHit != track->recHitsEnd();
00197                         ++ recHit ) {
00198 
00199                   
00200                   ++fnHit;
00201                   //std::cout << "fnHit="<< fnHit << std::endl;
00202                   
00203                   DetId id((*recHit)->geographicalId());
00204 
00205                   if ( (unsigned int)id.subdetId() == StripSubdetector::TIB ) {
00206                           ++fnStripHit;
00207                           ++fnTIBHit;
00208                   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TOB ) {
00209                           ++fnStripHit;
00210                           ++fnTOBHit;
00211                   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TID ) {
00212                           ++fnStripHit;
00213                           ++fnTIDHit;
00214                   } else if ( (unsigned int)id.subdetId() == StripSubdetector::TEC ) {
00215                           ++fnStripHit;
00216                           ++fnTECHit;
00217                   } else if ( (unsigned int)id.subdetId() == PixelSubdetector::PixelBarrel ) {
00218                           ++fnPixelHit;
00219                           ++fnPXBHit;
00220                   } else if ( (unsigned int)id.subdetId() == PixelSubdetector::PixelEndcap ) {
00221                           ++fnPixelHit;
00222                           ++fnPXFHit;
00223                   }
00224           }
00225 
00226           ftree_->Fill();
00227 
00228           ftotal_tracks++;
00229           // track quality
00230           if (fnStripHit>=8 && fnPixelHit >= 2 && fchi2/fndof<5 && fpt>2 && std::abs(fd0)<0.9) {
00231                   //fchi2/fndof<5 && fpt>4 && std::abs(fd0)<0.1 && TMath::Prob(fchi2,((int)fndof))>0.02 ) {
00232                   //if 
00233                   fBSvector.push_back(BSTrkParameters(fz0,fsigmaz0,fd0,fsigmad0,fphi0,fpt,0.,0.));
00234           }
00235           
00236     
00237         }
00238   
00239         ftotalevents++;
00240 }

void BeamSpotAnalyzer::beginJob ( const edm::EventSetup  )  [private, virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 245 of file BeamSpotAnalyzer.cc.

00246 {
00247 }

void BeamSpotAnalyzer::endJob ( void   )  [private, virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 250 of file BeamSpotAnalyzer.cc.

References reco::BeamSpot::BeamWidth(), GenMuonPlsPt100GeV_cfg::cout, reco::BeamSpot::covariance(), reco::BeamSpot::dxdz(), reco::BeamSpot::dydz(), e, lat::endl(), fasciiFile, fBSvector, fd0phi_chi2, fd0phi_d0, file_, BSFitter::Fit(), BSFitter::Fit_d0phi(), BSFitter::Fit_ited0phi(), ftotal_tracks, ftotalevents, ftree_, BSFitter::GetData(), BSFitter::GetResPar0(), BSFitter::GetResPar0Err(), BSFitter::GetResPar1(), BSFitter::GetResPar1Err(), i, inputBeamWidth_, edm::Service< T >::isAvailable(), j, python::TagTree::newtree, reco::BeamSpot::position(), funct::pow(), runallfitters_, BeamSpotObjects::SetBeamWidth(), BeamSpotObjects::SetCovariance(), BSFitter::Setd0Cut_d0phi(), BeamSpotObjects::Setdxdz(), BeamSpotObjects::Setdydz(), BSFitter::SetFitType(), BSFitter::SetFitVariable(), BeamSpotObjects::SetPosition(), BeamSpotObjects::SetSigmaZ(), reco::BeamSpot::sigmaZ(), write2DB_, reco::BeamSpot::x0(), reco::BeamSpot::y0(), and reco::BeamSpot::z0().

00250                          {
00251 
00252         std::cout << "\n-------------------------------------\n" << std::endl;
00253         std::cout << "\n Total number of events processed: "<< ftotalevents << std::endl;
00254         std::cout << "\n-------------------------------------\n\n" << std::endl;
00255         std::cout << " calculating beam spot..." << std::endl;
00256         std::cout << " we will use " << fBSvector.size() << " good tracks out of " << ftotal_tracks << std::endl;
00257 
00258         // default fit to extract beam spot info
00259         BSFitter *myalgo = new BSFitter( fBSvector );
00260         reco::BeamSpot beam_default = myalgo->Fit();
00261         std::cout << "\n RESULTS OF DEFAULT FIT:" << std::endl;
00262         std::cout << beam_default << std::endl;
00263 
00264         // dump to file
00265         fasciiFile << "X " << beam_default.x0() << std::endl;
00266         fasciiFile << "Y " << beam_default.y0() << std::endl;
00267         fasciiFile << "Z " << beam_default.z0() << std::endl;
00268         fasciiFile << "sigmaZ " << beam_default.sigmaZ() << std::endl;
00269         fasciiFile << "dxdz " << beam_default.dxdz() << std::endl;
00270         fasciiFile << "dydz " << beam_default.dydz() << std::endl;
00271         if (inputBeamWidth_ > 0 ) {
00272                 fasciiFile << "BeamWidth " << inputBeamWidth_ << std::endl;
00273         } else {
00274                 fasciiFile << "BeamWidth " << beam_default.BeamWidth() << std::endl;
00275         }
00276         
00277         for (int i = 0; i<6; ++i) {
00278                 fasciiFile << "Cov("<<i<<",j) ";
00279                 for (int j=0; j<7; ++j) {
00280                         fasciiFile << beam_default.covariance(i,j) << " ";
00281                 }
00282                 fasciiFile << std::endl;
00283         }
00284         // beam width error
00285         if (inputBeamWidth_ > 0 ) {
00286                 fasciiFile << "Cov(6,j) 0 0 0 0 0 0 " << pow(2.e-4,2) << std::endl;
00287         } else {
00288                 fasciiFile << "Cov(6,j) 0 0 0 0 0 0 " << beam_default.covariance(6,6) << std::endl;
00289         }
00290         
00291         
00292         if (write2DB_) {
00293                 std::cout << "\n-------------------------------------\n\n" << std::endl;
00294                 std::cout << " write results to DB..." << std::endl;
00295 
00296                 BeamSpotObjects *pBSObjects = new BeamSpotObjects();
00297 
00298                 //pBSObjects->Put(beam_default);
00299                 pBSObjects->SetPosition(beam_default.position().X(),beam_default.position().Y(),beam_default.position().Z());
00300                 //std::cout << " wrote: x= " << beam_default.position().X() << " y= "<< beam_default.position().Y() << " z= " << beam_default.position().Z() << std::endl;
00301                 pBSObjects->SetSigmaZ(beam_default.sigmaZ());
00302                 pBSObjects->Setdxdz(beam_default.dxdz());
00303                 pBSObjects->Setdydz(beam_default.dydz());
00304                 if (inputBeamWidth_ > 0 ) {
00305                         std::cout << " beam width value forced to be " << inputBeamWidth_ << std::endl;
00306                         pBSObjects->SetBeamWidth(inputBeamWidth_);
00307                 } else {
00308                         std::cout << " using default value, 15e-4, for beam width!"<<std::endl;
00309                         pBSObjects->SetBeamWidth(15.0e-4);
00310                 }
00311                 
00312                 for (int i = 0; i<7; ++i) {
00313                   for (int j=0; j<7; ++j) {
00314                     pBSObjects->SetCovariance(i,j,beam_default.covariance(i,j));
00315                   }
00316                 }
00317                 edm::Service<cond::service::PoolDBOutputService> poolDbService;
00318                 if( poolDbService.isAvailable() ) {
00319                   std::cout << "poolDBService available"<<std::endl;
00320                   if ( poolDbService->isNewTagRequest( "BeamSpotObjectsRcd" ) ) {
00321                     std::cout << "new tag requested" << std::endl;
00322                     poolDbService->createNewIOV<BeamSpotObjects>( pBSObjects, poolDbService->beginOfTime(),poolDbService->endOfTime(),
00323                                                                   "BeamSpotObjectsRcd"  );
00324                   }
00325                   else {
00326                     std::cout << "no new tag requested" << std::endl;
00327                     poolDbService->appendSinceTime<BeamSpotObjects>( pBSObjects, poolDbService->currentTime(),
00328                                                                      "BeamSpotObjectsRcd" );
00329                   }
00330 
00331                 
00332                 }
00333         }
00334 
00335         if (runallfitters_) {
00336         
00337           
00338         // add new branches
00339         std::cout << " add new branches to output file " << std::endl;
00340         beam_default = myalgo->Fit_d0phi();
00341         file_->cd();
00342         TTree *newtree = new TTree("mytreecorr","mytreecorr");
00343         newtree->Branch("d0phi_chi2",&fd0phi_chi2,"fd0phi_chi2/D");
00344         newtree->Branch("d0phi_d0",&fd0phi_d0,"fd0phi_d0/D");
00345         newtree->SetBranchAddress("d0phi_chi2",&fd0phi_chi2);
00346         newtree->SetBranchAddress("d0phi_d0",&fd0phi_d0);
00347         std::vector<BSTrkParameters>  tmpvector = myalgo->GetData();
00348         
00349         std::vector<BSTrkParameters>::iterator iparam = tmpvector.begin();
00350         for( iparam = tmpvector.begin() ;
00351                  iparam != tmpvector.end() ; ++iparam) {
00352                 fd0phi_chi2 = iparam->d0phi_chi2();
00353                 fd0phi_d0   = iparam->d0phi_d0();
00354                 newtree->Fill();
00355         }
00356         newtree->Write();
00357 
00358         // iterative
00359         std::cout << " d0-phi Iterative:" << std::endl;
00360         BSFitter *myitealgo = new BSFitter( fBSvector );
00361         myitealgo->Setd0Cut_d0phi(4.0);
00362         reco::BeamSpot beam_ite = myitealgo->Fit_ited0phi();
00363         std::cout << beam_ite << std::endl;
00364 
00365         
00366         std::cout << "\n Now run tests of the different fits\n";
00367         // from here are just tests
00368         std::string fit_type = "chi2";
00369         myalgo->SetFitVariable(std::string("z"));
00370         myalgo->SetFitType(std::string("chi2"));
00371         reco::BeamSpot beam_fit_z_chi2 = myalgo->Fit();
00372         std::cout << " z Chi2 Fit ONLY:" << std::endl;
00373         std::cout << beam_fit_z_chi2 << std::endl;
00374         
00375         
00376         fit_type = "combined";
00377         myalgo->SetFitVariable("z");
00378         myalgo->SetFitType("combined");
00379         reco::BeamSpot beam_fit_z_lh = myalgo->Fit();
00380         std::cout << " z Log-Likelihood Fit ONLY:" << std::endl;
00381         std::cout << beam_fit_z_lh << std::endl;
00382 
00383         
00384         myalgo->SetFitVariable("d");
00385         myalgo->SetFitType("d0phi");
00386         reco::BeamSpot beam_fit_dphi = myalgo->Fit();
00387         std::cout << " d0-phi0 Fit ONLY:" << std::endl;
00388         std::cout << beam_fit_dphi << std::endl;
00389 
00390                 
00391         myalgo->SetFitVariable(std::string("d*z"));
00392         myalgo->SetFitType(std::string("likelihood"));
00393         reco::BeamSpot beam_fit_dz_lh = myalgo->Fit();
00394         std::cout << " Log-Likelihood Fit:" << std::endl;
00395         std::cout << beam_fit_dz_lh << std::endl;
00396 
00397         
00398         myalgo->SetFitVariable(std::string("d*z"));
00399         myalgo->SetFitType(std::string("resolution"));
00400         reco::BeamSpot beam_fit_dresz_lh = myalgo->Fit();
00401         std::cout << " IP Resolution Fit" << std::endl;
00402         std::cout << beam_fit_dresz_lh << std::endl;
00403 
00404         std::cout << "c0 = " << myalgo->GetResPar0() << " +- " << myalgo->GetResPar0Err() << std::endl;
00405         std::cout << "c1 = " << myalgo->GetResPar1() << " +- " << myalgo->GetResPar1Err() << std::endl;
00406 
00407         }
00408 
00409         // let's close everything
00410     file_->cd();
00411         ftree_->Write();
00412         file_->Close();
00413                 
00414         std::cout << "[BeamSpotAnalyzer] endJob done \n" << std::endl;
00415 }


Member Data Documentation

std::string BeamSpotAnalyzer::ckfSeedProducerLabel_ [private]

Definition at line 77 of file BeamSpotAnalyzer.h.

std::string BeamSpotAnalyzer::ckfTrackCandidateProducerLabel_ [private]

Definition at line 78 of file BeamSpotAnalyzer.h.

std::string BeamSpotAnalyzer::ckfTrackProducerLabel_ [private]

Definition at line 79 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

std::ofstream BeamSpotAnalyzer::fasciiFile [private]

Definition at line 44 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer(), and endJob().

std::string BeamSpotAnalyzer::fasciiFileName [private]

Definition at line 45 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer().

std::string BeamSpotAnalyzer::fasciifilename [private]

Definition at line 43 of file BeamSpotAnalyzer.h.

std::vector< BSTrkParameters > BeamSpotAnalyzer::fBSvector [private]

Definition at line 75 of file BeamSpotAnalyzer.h.

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

int BeamSpotAnalyzer::fcharge [private]

Definition at line 54 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fchi2 [private]

Definition at line 55 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fcov[7][7] [private]

Definition at line 73 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fd0 [private]

Definition at line 58 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fd0phi_chi2 [private]

Definition at line 71 of file BeamSpotAnalyzer.h.

Referenced by endJob().

double BeamSpotAnalyzer::fd0phi_d0 [private]

Definition at line 72 of file BeamSpotAnalyzer.h.

Referenced by endJob().

double BeamSpotAnalyzer::feta [private]

Definition at line 53 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

TFile* BeamSpotAnalyzer::file_ [private]

Definition at line 47 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer(), and endJob().

int BeamSpotAnalyzer::fmaxNtracks [private]

Definition at line 84 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer().

double BeamSpotAnalyzer::fndof [private]

Definition at line 56 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnHit [private]

Definition at line 62 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnPixelHit [private]

Definition at line 64 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnPXBHit [private]

Definition at line 69 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnPXFHit [private]

Definition at line 70 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnStripHit [private]

Definition at line 63 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnTECHit [private]

Definition at line 68 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnTIBHit [private]

Definition at line 65 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnTIDHit [private]

Definition at line 67 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

unsigned int BeamSpotAnalyzer::fnTOBHit [private]

Definition at line 66 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fphi0 [private]

Definition at line 57 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fpt [private]

Definition at line 52 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

float BeamSpotAnalyzer::fptmin [private]

Definition at line 83 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer().

double BeamSpotAnalyzer::fsigmad0 [private]

Definition at line 59 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::fsigmaz0 [private]

Definition at line 61 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::ftheta [private]

Definition at line 51 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

int BeamSpotAnalyzer::ftotal_tracks [private]

Definition at line 88 of file BeamSpotAnalyzer.h.

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

int BeamSpotAnalyzer::ftotalevents [private]

Definition at line 50 of file BeamSpotAnalyzer.h.

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

TTree* BeamSpotAnalyzer::ftree_ [private]

Definition at line 48 of file BeamSpotAnalyzer.h.

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

double BeamSpotAnalyzer::fz0 [private]

Definition at line 60 of file BeamSpotAnalyzer.h.

Referenced by analyze(), and BeamSpotAnalyzer().

double BeamSpotAnalyzer::inputBeamWidth_ [private]

Definition at line 89 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer(), and endJob().

std::string BeamSpotAnalyzer::outputfilename_ [private]

Definition at line 42 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer().

bool BeamSpotAnalyzer::runallfitters_ [private]

Definition at line 87 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer(), and endJob().

unsigned int BeamSpotAnalyzer::sameNumberOfTracks [private]

Definition at line 81 of file BeamSpotAnalyzer.h.

bool BeamSpotAnalyzer::write2DB_ [private]

Definition at line 86 of file BeamSpotAnalyzer.h.

Referenced by BeamSpotAnalyzer(), and endJob().


The documentation for this class was generated from the following files:
Generated on Tue Jun 9 18:15:00 2009 for CMSSW by  doxygen 1.5.4