#include <MuonAcceptor.h>
Inheritance diagram for Atlfast::MuonAcceptor:
Public Member Functions | |
MuonAcceptor (int muSmearKey, MsgStream &log) | |
virtual | ~MuonAcceptor () |
virtual bool | accept (const ReconstructedParticle &particle, MsgStream &log) |
int | muReadEff (std::ifstream &InputFile, vector< muEffdata > &VectorOfmuEffdata) |
vector< pair< int, int > > | muEffGetPaire (double pt, double etamin, double etamax) |
pair< double, double > | muGetPtLowHigh (double pt) |
int | muGetPaireNumeDenoPosi (double pt, double etamin, double etamax) |
pair< double, double > | muEfficiency (double ptIn, double etaIn, double phi) |
pair< double, double > | muEff (double ptIn, double etaIn, double phi) |
int | muGetEtaBin (double etaIn) |
int | muGetPhiBin (double phi) |
int | muGetPhiBin1 (double phi) |
pair< double, double > | muPhiEff (double pt, double etamin, double etamax, double phimin, double phimax) |
pair< double, double > | getPhiEff (double pt, double effin, double phig, double etamin, double etamax) |
void | muPhiEfficiency (double pt, double etaIn, double phi, double &eff, double &eeff) |
Private Attributes | |
vector< muEffdata > | m_muEffdataObject |
vector< double > | m_ptValues |
vector< double > | m_etaminValues |
vector< double > | m_etamaxValues |
MsgStream & | m_log |
HepRandomEngine * | m_randomEngine |
Definition at line 102 of file MuonAcceptor.h.
Atlfast::MuonAcceptor::MuonAcceptor | ( | int | muSmearKey, | |
MsgStream & | log | |||
) |
Definition at line 92 of file MuonAcceptor.cxx.
00092 : m_log(log) { 00093 00094 m_log << MSG::DEBUG << "MuonAcceptor::MuonAcceptor" << endreq; 00095 // Initialize random number generator 00096 m_randomEngine = new HepJamesRandom(12345); 00097 00098 //Open input file 00099 00100 std::string filename = muSmearKey == 1 ? 00101 "atlfastDatafiles/MuonSpectrometerEff.txt" : 00102 "atlfastDatafiles/MuonCombinedEff.txt"; 00103 filename = PathResolver::find_file (filename, "DATAPATH"); 00104 m_log << MSG::DEBUG << "Opening file: " << filename << endreq; 00105 std::ifstream InputFile; 00106 InputFile.open(filename.c_str()); 00107 if(!InputFile) { 00108 m_log << MSG::ERROR << "Error: file could not be opened" << endreq; 00109 exit(1); 00110 } 00111 // 00112 std::vector<muEffdata>VectorOfmuEffdata; 00113 std::vector<std::pair<int, int> > ptEtaNumeDeno; 00114 VectorOfmuEffdata.clear(); 00115 ptEtaNumeDeno.clear(); 00116 // 00117 m_muEffdataObject.clear(); 00118 m_ptValues.clear(); 00119 m_etaminValues.clear(); 00120 m_etamaxValues.clear(); 00121 00122 00123 //Loop on evts 00124 int Istatus = 1; 00125 while ( Istatus == 1 ) { 00126 00127 // Get evt data 00128 // 00129 Istatus = muReadEff(InputFile,VectorOfmuEffdata); 00130 // Do something 00131 if (Istatus == 1){ 00132 00133 std::vector<int> ptEtaNum; 00134 int SizeOfVectorOfmuEffdata = VectorOfmuEffdata.size() ; 00135 for (int Item=0; Item <SizeOfVectorOfmuEffdata ; Item++){ 00136 double pt;double etamin;double etamax; 00137 VectorOfmuEffdata[Item].GetPtEtaNumeDeno(pt,etamin,etamax,ptEtaNumeDeno); 00138 m_muEffdataObject.push_back(VectorOfmuEffdata[Item]); 00139 m_ptValues.push_back(pt); 00140 m_etaminValues.push_back(etamin); 00141 m_etamaxValues.push_back(etamax); 00142 } 00143 } 00144 } 00145 //Close input file 00146 InputFile.close(); 00147 m_log << MSG::DEBUG << "Done with MuonAcceptor constructor" << endreq; 00148 }
virtual Atlfast::MuonAcceptor::~MuonAcceptor | ( | ) | [inline, virtual] |
bool Atlfast::MuonAcceptor::accept | ( | const ReconstructedParticle & | particle, | |
MsgStream & | log | |||
) | [virtual] |
Accept function, must be overridden in derived class
Implements Atlfast::IAcceptor.
Definition at line 150 of file MuonAcceptor.cxx.
00150 { 00151 00152 log << MSG::DEBUG << "MuonAcceptor::accept" << endreq; 00153 00154 // The efficiencies are a function of pT in GeV 00155 double pTinGeV = particle.pT() / GeV; 00156 00157 double muEff = muEfficiency(pTinGeV, 00158 particle.eta(), 00159 particle.phi()).first; 00160 log << MSG::DEBUG << "pT,eta,phi,efficiency = " << pTinGeV 00161 << "," << particle.eta() << "," 00162 << particle.phi() << "," << muEff << endreq; 00163 double random = m_randomEngine->flat(); 00164 log << MSG::DEBUG << "random = " << random << endreq; 00165 return ( random < muEff) ? true : false; 00166 }
int Atlfast::MuonAcceptor::muReadEff | ( | std::ifstream & | InputFile, | |
vector< muEffdata > & | VectorOfmuEffdata | |||
) |
std::vector< std::pair< int, int > > Atlfast::MuonAcceptor::muEffGetPaire | ( | double | pt, | |
double | etamin, | |||
double | etamax | |||
) |
Definition at line 776 of file MuonAcceptor.cxx.
00776 { 00777 00778 std::vector<std::pair<int, int> > vecPaireNumeDeno; 00779 vecPaireNumeDeno.clear(); 00780 // 00781 std::vector<std::pair<int, int> > ptEtaNumeDeno; 00782 ptEtaNumeDeno.clear(); 00783 // 00784 int iposition = muGetPaireNumeDenoPosi(pt,etamin,etamax); 00785 00786 if(iposition != -1){ 00787 int i=iposition; 00788 m_muEffdataObject[i].GetPtEtaNumeDeno(pt,etamin,etamax,ptEtaNumeDeno); 00789 vecPaireNumeDeno=ptEtaNumeDeno; 00790 } 00791 // 00792 return vecPaireNumeDeno; 00793 }
std::pair< double, double > Atlfast::MuonAcceptor::muGetPtLowHigh | ( | double | pt | ) |
Definition at line 810 of file MuonAcceptor.cxx.
00810 { 00811 00812 std::pair<double,double> ptLowHigh; 00813 00814 int sizeOfptvalues = m_ptValues.size(); 00815 double ptmin= m_ptValues.at(0); 00816 double ptmax= m_ptValues.at(sizeOfptvalues-1); 00817 00818 for (int i=0; i <sizeOfptvalues-1 ; i++){ 00819 ptmin= m_ptValues[i]; 00820 ptmax= m_ptValues[i+1]; 00821 if(pt >=ptmin && ptmin<ptmax){ 00822 double ptlow =m_ptValues[i] ; 00823 double pthigh =m_ptValues[i+1] ; 00824 ptLowHigh=make_pair(ptlow,pthigh); 00825 } 00826 } 00827 if(pt>=m_ptValues.at(sizeOfptvalues-1)){ 00828 double ptlow = m_ptValues.at(sizeOfptvalues-2); 00829 double pthigh = m_ptValues.at(sizeOfptvalues-1); 00830 ptLowHigh=make_pair(ptlow,pthigh); 00831 } 00832 if(pt <=m_ptValues.at(0)){ 00833 double ptlow = m_ptValues.at(0); 00834 double pthigh = m_ptValues.at(0) ; 00835 ptLowHigh=make_pair(ptlow,pthigh); 00836 } 00837 return ptLowHigh ; 00838 }
int Atlfast::MuonAcceptor::muGetPaireNumeDenoPosi | ( | double | pt, | |
double | etamin, | |||
double | etamax | |||
) |
Definition at line 795 of file MuonAcceptor.cxx.
00795 { 00796 00797 int iposition =-1; 00798 const double epsilon=1.e-6; 00799 for(unsigned int i= 0;i<m_ptValues.size();i++){ 00800 if( fabs(pt-m_ptValues[i]) < epsilon && 00801 fabs(etamin - m_etaminValues[i]) < epsilon && 00802 fabs( etamax - m_etamaxValues[i]) < epsilon 00803 ) 00804 iposition=i; 00805 } 00806 return iposition; 00807 }
std::pair< double, double > Atlfast::MuonAcceptor::muEfficiency | ( | double | ptIn, | |
double | etaIn, | |||
double | phi | |||
) |
Definition at line 427 of file MuonAcceptor.cxx.
00428 { 00429 00430 double efficiency =0; 00431 double errEfficiency =0; 00432 00433 double eta = fabs(etaIn); 00434 double pt = fabs(ptIn); 00435 00436 // Interpolate between pt points 00437 double pt1 = 0.; 00438 double pt2 = 0.; 00439 00440 // Get the two pt points 00441 std::pair<double,double>ptLowHigh; 00442 ptLowHigh= muGetPtLowHigh(pt); 00443 pt1 = ptLowHigh.first; 00444 pt2 = ptLowHigh.second; 00445 00446 if(pt >= pt1 && pt < pt2){ 00447 std::pair<double,double> eff1_errff1; 00448 std::pair<double,double> eff2_errff2; 00449 00450 eff1_errff1 = muEff(pt1,eta,phi); 00451 eff2_errff2 = muEff(pt2,eta,phi); 00452 00453 double eff1 = eff1_errff1.first ; 00454 double eeff1 = eff1_errff1.second; 00455 00456 double eff2 = eff2_errff2.first; 00457 double eeff2 = eff2_errff2.second; 00458 00459 efficiency = eff2 *log(pt/pt1)/log(pt2/pt1); 00460 efficiency = efficiency + eff1 *log(pt2/pt)/log(pt2/pt1); 00461 errEfficiency = pow((eeff2*log(pt/pt1)/log(pt2/pt1)),2); 00462 errEfficiency = errEfficiency +pow((eeff1*log(pt2/pt)/log(pt2/pt1)),2); 00463 errEfficiency = sqrt(errEfficiency); 00464 } 00465 00466 if(efficiency < 0){ 00467 efficiency =0.; 00468 errEfficiency=0.; 00469 } 00470 if(efficiency > 1){ 00471 efficiency =1.; 00472 errEfficiency=0.; 00473 } 00474 00475 return make_pair(efficiency,errEfficiency); 00476 }
std::pair< double, double > Atlfast::MuonAcceptor::muEff | ( | double | ptIn, | |
double | etaIn, | |||
double | phi | |||
) |
Definition at line 479 of file MuonAcceptor.cxx.
00480 { 00481 double eta = fabs(etaIn); 00482 double pt = fabs(ptIn); 00483 00484 double xEtaMinOut = 0.; 00485 double xEtaMaxOut = 2.5; 00486 00487 int nNumOut = 0; 00488 int nDenoOut = 0; 00489 00490 int nNumTotOut = 0; 00491 int nDenoTotOut= 0; 00492 // 00493 int iEtaBin = muGetEtaBin(eta); 00494 // 00495 std::vector<std::pair<int, int> > ptEtaNumeDeno; 00496 ptEtaNumeDeno.clear(); 00497 ptEtaNumeDeno= muEffGetPaire(pt,xEtaMinOut,xEtaMaxOut); 00498 00499 nNumOut=ptEtaNumeDeno[iEtaBin].first; 00500 nDenoOut=ptEtaNumeDeno[iEtaBin].second; 00501 00502 int nBinOut = ptEtaNumeDeno.size(); 00503 00504 for(unsigned int i=0 ; i<ptEtaNumeDeno.size() ; i++){ 00505 nNumTotOut = nNumTotOut+ptEtaNumeDeno[i].first; 00506 nDenoTotOut = nDenoTotOut+ptEtaNumeDeno[i].second; 00507 } 00508 00509 double eff1 = double(nNumOut)/double(nDenoOut); 00510 double eeff1 = sqrt(eff1*(1.-eff1)/double(nDenoOut)); 00511 double delta = (xEtaMaxOut-xEtaMinOut)/double(nBinOut) ; 00512 double xcbin1 = xEtaMinOut+delta/2.+delta*double(iEtaBin) ; 00513 // 00514 muPhiEfficiency(pt,xcbin1,phi,eff1,eeff1); 00515 // 00516 return make_pair(eff1,eeff1); 00517 }
int Atlfast::MuonAcceptor::muGetEtaBin | ( | double | etaIn | ) |
Definition at line 519 of file MuonAcceptor.cxx.
00519 { 00520 00521 int iEtaBin = -1; 00522 double eta = fabs(etaIn); 00523 00524 double xMaxOut =2.5; 00525 double xMinOut = 0.; 00526 double pt=50.; 00527 // 00528 std::vector<std::pair<int, int> > ptEtaNumeDeno; 00529 ptEtaNumeDeno.clear(); 00530 ptEtaNumeDeno= muEffGetPaire(pt,xMinOut,xMaxOut); 00531 // 00532 int NbinOut = ptEtaNumeDeno.size(); 00533 00534 double delta=double(xMaxOut-xMinOut)/double(NbinOut); 00535 iEtaBin = int(floor(eta/delta)); 00536 if( iEtaBin >= NbinOut) iEtaBin= NbinOut-1; 00537 00538 return iEtaBin; 00539 }
int Atlfast::MuonAcceptor::muGetPhiBin | ( | double | phi | ) |
Definition at line 617 of file MuonAcceptor.cxx.
00617 { 00618 00619 double pi = 3.1415927; 00620 int iPhiBin = -1; 00621 00622 double xMaxOut =360; 00623 double xMinOut = 0.; 00624 int NbinOut = 20; 00625 00626 double phig = phi; 00627 if (phi < 0.) { 00628 phig = pi*2 + phi; 00629 } 00630 double phideg = phig*180./ pi; 00631 00632 double delta=double(xMaxOut-xMinOut)/double(NbinOut); 00633 iPhiBin = int(floor(phideg/delta)); 00634 if( iPhiBin >= NbinOut) iPhiBin= NbinOut-1; 00635 00636 return iPhiBin; 00637 }
int Atlfast::MuonAcceptor::muGetPhiBin1 | ( | double | phi | ) |
Definition at line 639 of file MuonAcceptor.cxx.
00639 { 00640 00641 double pi = 3.1415927; 00642 int iPhiBin = -1; 00643 00644 double xMaxOut = 22.5; 00645 double xMinOut = 0.; 00646 int NbinOut = 20; 00647 00648 double phig = phi; 00649 if (phi < 0.) { 00650 phig = pi*2 + phi; 00651 } 00652 double phideg = phig*180./ pi; 00653 00654 double delta=double(xMaxOut-xMinOut)/double(NbinOut); 00655 iPhiBin = int(floor(phideg/delta)); 00656 if( iPhiBin >= NbinOut) iPhiBin= NbinOut-1; 00657 00658 return iPhiBin; 00659 }
std::pair< double, double > Atlfast::MuonAcceptor::muPhiEff | ( | double | pt, | |
double | etamin, | |||
double | etamax, | |||
double | phimin, | |||
double | phimax | |||
) |
Definition at line 740 of file MuonAcceptor.cxx.
00740 { 00741 00742 int NnumTot=0 ; 00743 int NdenoTot=0; 00744 00745 double eff =0; 00746 double eeff=0; 00747 00748 int ibinmin; 00749 int ibinmax; 00750 00751 if(etamin==0. && etamax==0.15){ 00752 ibinmin = muGetPhiBin1(phimin) ; 00753 ibinmax = muGetPhiBin1(phimax) ; 00754 } 00755 else{ 00756 ibinmin = muGetPhiBin(phimin) ; 00757 ibinmax = muGetPhiBin(phimax) ; 00758 } 00759 std::vector<std::pair<int, int> > ptEtaNumeDeno; 00760 ptEtaNumeDeno.clear(); 00761 ptEtaNumeDeno= muEffGetPaire(pt,etamin,etamax); 00762 00763 for (int ibin=ibinmin;ibin<=ibinmax;ibin++){ 00764 NnumTot =NnumTot+ptEtaNumeDeno[ibin].first; 00765 NdenoTot=NdenoTot+ptEtaNumeDeno[ibin].second; 00766 } 00767 00768 eff =double(NnumTot)/double(NdenoTot); 00769 eeff=sqrt(eff*(1.-eff)/double(NdenoTot)); 00770 // 00771 return make_pair(eff,eeff); 00772 }
std::pair< double, double > Atlfast::MuonAcceptor::getPhiEff | ( | double | pt, | |
double | effin, | |||
double | phig, | |||
double | etamin, | |||
double | etamax | |||
) |
Definition at line 663 of file MuonAcceptor.cxx.
00664 { 00665 double eff = 1.; 00666 double eeff =0.; 00667 00668 double phivalue[21]; 00669 double efficiency[21]; 00670 double errefficiency[21]; 00671 00672 for (int j = 0; j < 21; j++) { 00673 efficiency[j] = 0.; 00674 errefficiency[j] = 0.; 00675 phivalue[j] = 0.; 00676 } 00677 // 00678 double pi = 3.1415927; 00679 double degrad = pi/180.; 00680 // 00681 for (int j = 0; j < 21; j++) { 00682 phivalue[j] = j*18.*degrad; 00683 if(etamin ==0. && etamax ==0.15)phivalue[j] = j*1.125*degrad; 00684 } 00685 // 00686 // Mean Efficiency in phi 00687 00688 std::pair<double,double> eff_errff; 00689 double phiLow = phivalue[0]; 00690 double phiHigh = phivalue[20]; 00691 eff_errff = muPhiEff(pt,etamin,etamax,phiLow,phiHigh); 00692 efficiency[20] = eff_errff.first ; 00693 errefficiency[20] = eff_errff.second; 00694 00695 // The other bins 00696 double delta2 = 8.9*degrad; 00697 if(etamin ==0. && etamax ==0.15)delta2 = 0.55*degrad; 00698 for (int j = 0; j < 20; j++) { 00699 double phicenter = (phivalue[j+1]+phivalue[j])/2.; 00700 double phimin = phicenter - delta2; 00701 double phimax = phicenter + delta2; 00702 eff_errff = muPhiEff(pt,etamin,etamax,phimin,phimax); 00703 efficiency[j] = eff_errff.first ; 00704 errefficiency[j] = eff_errff.second; 00705 } 00706 // 00707 double averaefficiency = effin/efficiency[20]; 00708 double factornorm = averaefficiency; 00709 if(factornorm > 1. || factornorm < 0.9)factornorm =1; 00710 // 00711 if (phig >= phivalue[0] && phig <= phivalue[1]) { 00712 eff = efficiency[0] * factornorm; 00713 eeff = errefficiency[0]; 00714 } 00715 for (int j = 1; j < 20; j++) { 00716 if (phig > phivalue[j] && phig <= phivalue[j+1]) { 00717 eff = efficiency[j] * factornorm; 00718 eeff= errefficiency[j]; 00719 } 00720 } 00721 if (phig > phivalue[20] && phig <= phivalue[21]) { 00722 eff = efficiency[19] * factornorm; 00723 eeff= errefficiency[19]; 00724 } 00725 // 00726 if (eff < 0.) { 00727 eff = 0.; 00728 eeff = 0.; 00729 } 00730 if (eff > 1.) { 00731 eff = 1.; 00732 eeff =0.; 00733 } 00734 // 00735 return make_pair(eff,eeff); 00736 }
void Atlfast::MuonAcceptor::muPhiEfficiency | ( | double | pt, | |
double | etaIn, | |||
double | phi, | |||
double & | eff, | |||
double & | eeff | |||
) |
Definition at line 541 of file MuonAcceptor.cxx.
00541 { 00542 00543 double pi = 3.1415927; 00544 double phig = phi; 00545 if (phi < 0.) { 00546 phig = pi*2 + phi; 00547 } 00548 00549 double effin = eff; 00550 //double eeffin = eeff; 00551 00552 double eta = fabs(etaIn); 00553 00554 //--------------------------- 00555 //Eta region 0. <|eta|<0.15 00556 //--------------------------- 00557 00558 double etamin1 = 0.; 00559 double etamax1 = 0.15; 00560 if (eta >= etamin1 && eta < etamax1) { 00561 double pion4 = pi/4.; 00562 double pion8 = pi/8.; 00563 double phig2 = fmod(phig,pion4); 00564 if (phig2 >= pion8) phig2 = pion4-phig2; 00565 // 00566 std::pair<double,double> eff_erreff; 00567 eff_erreff = getPhiEff(pt,effin,phig2,etamin1,etamax1); 00568 eff = eff_erreff.first ; 00569 eeff = eff_erreff.second; 00570 } 00571 00572 //------------------------------------------------------ 00573 //Feets: Eta region 0.3 <|eta|<0.4 && 0.8 <|eta|<0.9 00574 //------------------------------------------------------ 00575 00576 double etamin2 = 0.3; 00577 double etamax2 = 0.4; 00578 double etamin22 = 0.8; 00579 double etamax22 = 0.9; 00580 if ((eta >= etamin2 && eta < etamax2) || 00581 (eta >= etamin22 && eta < etamax22) ) { 00582 std::pair<double,double> eff_erreff; 00583 eff_erreff = getPhiEff(pt,effin,phig,etamin2,etamax2); 00584 eff = eff_erreff.first ; 00585 eeff = eff_erreff.second; 00586 } 00587 //--------------------------- 00588 // Eta region 1.05 <|eta|<1.15 00589 //--------------------------- 00590 double etamin3 =1.05; 00591 double etamax3 =1.15; 00592 if (eta >= etamin3 && eta < etamax3) { 00593 00594 std::pair<double,double> eff_erreff; 00595 eff_erreff = getPhiEff(pt,effin,phig,etamin3,etamax3); 00596 eff = eff_erreff.first ; 00597 eeff = eff_erreff.second; 00598 } 00599 // 00600 //--------------------------- 00601 // Eta region 1.20 <|eta|<1.25 00602 //--------------------------- 00603 double etamin4 =1.20; 00604 double etamax4 =1.25; 00605 00606 if (eta >= etamin4 && eta < etamax4) { 00607 00608 std::pair<double,double> eff_erreff; 00609 eff_erreff = getPhiEff(pt,effin,phig,etamin4,etamax4); 00610 eff = eff_erreff.first ; 00611 eeff = eff_erreff.second; 00612 } 00613 00614 }
vector<muEffdata> Atlfast::MuonAcceptor::m_muEffdataObject [private] |
Definition at line 136 of file MuonAcceptor.h.
vector<double> Atlfast::MuonAcceptor::m_ptValues [private] |
Definition at line 137 of file MuonAcceptor.h.
vector<double> Atlfast::MuonAcceptor::m_etaminValues [private] |
Definition at line 138 of file MuonAcceptor.h.
vector<double> Atlfast::MuonAcceptor::m_etamaxValues [private] |
Definition at line 139 of file MuonAcceptor.h.
MsgStream& Atlfast::MuonAcceptor::m_log [private] |
Definition at line 140 of file MuonAcceptor.h.
HepRandomEngine* Atlfast::MuonAcceptor::m_randomEngine [private] |
Definition at line 141 of file MuonAcceptor.h.