From 757734410c972032a93396a61dc698936783fe2d Mon Sep 17 00:00:00 2001 From: Antonio De Maria <antonio.de.maria@cern.ch> Date: Sat, 20 Jul 2024 07:38:44 +0200 Subject: [PATCH 1/2] Remove unused function in MMC --- .../DiTauMassTools/MissingMassCalculatorV2.h | 11 -- .../Root/MissingMassCalculatorV2.cxx | 110 +----------------- 2 files changed, 1 insertion(+), 120 deletions(-) diff --git a/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/MissingMassCalculatorV2.h b/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/MissingMassCalculatorV2.h index a3dd64f6856d..e5981213c6b8 100644 --- a/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/MissingMassCalculatorV2.h +++ b/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/MissingMassCalculatorV2.h @@ -243,8 +243,6 @@ class MissingMassCalculatorV2 { DitauStuff m_fDitauStuffFit; // results based on fit method DitauStuff m_fDitauStuffHisto; // results based on histo method - int m_fApplyMassScale; // switch to apply mass scale correction - int m_niter_fit1; // number of iterations for dR-dPhi scan int m_niter_fit2; // number of iterations for MET-scan int m_niter_fit3; // number of iterations for Mnu-scan @@ -253,8 +251,6 @@ class MissingMassCalculatorV2 { int m_RMSStop; int m_RndmSeedAltering; // reset seed (not necessary by default) - - int m_fJERsyst; // switch for JER systematics double m_dTheta3d_binMin; // minimal step size for dTheta3D double m_dTheta3d_binMax; // maximum step size for dTheta3D double m_dRmax_tau; // maximum dR(nu-visTau) @@ -267,9 +263,6 @@ class MissingMassCalculatorV2 { void DoOutputInfo(); void PrintOtherInput(); void PrintResults(); - int NuPsolution(TVector2 met_vec, double theta1, double phi1, - double theta2, double phi2, double &P1, double &P2); // keep this version for simple tests - inline int NuPsolutionV3(const double & mNu1, const double & mNu2, const double & phi1, const double & phi2, int & nsol1, int & nsol2); @@ -358,9 +351,6 @@ public: void SetdTheta3d_binMin(const double val) { m_dTheta3d_binMin=val; } // minimal step size for dTheta3D void SetEventNumber(const int eventNumber) { m_eventNumber = eventNumber; } - void SetJERsyst(const int val) { m_fJERsyst=val; } - void SetApplyMassScale(const int val) { m_fApplyMassScale=val; } - void SetMnuScanRange(const double val) { m_MnuScanRange=val; } void SetProposalTryMEt(const double val) {m_proposalTryMEt=val; } @@ -410,7 +400,6 @@ public: //-------- Get results; - int StandardCollApprox(const TLorentzVector & tau_vec1, const TLorentzVector & tau_vec2, const TVector2 & met_vec, double &Mrec); // standard collinear approximation Double_t maxFitting(Double_t *x, Double_t *par); // compute maximum from histo diff --git a/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassCalculatorV2.cxx b/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassCalculatorV2.cxx index b8b90eb62939..29cca5bc756b 100644 --- a/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassCalculatorV2.cxx +++ b/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassCalculatorV2.cxx @@ -74,9 +74,7 @@ MissingMassCalculatorV2::MissingMassCalculatorV2( Prob->SetUseDphiLL(false); // added by Tomas Davidek for lep-lep m_dTheta3d_binMin = 0.0025; m_dTheta3d_binMax = 0.02; - m_fJERsyst = 0; // no JER systematics by default (+/-1: up/down 1 sigma) preparedInput.m_METresSyst = 0; // no MET resolution systematics by default (+/-1: up/down 1 sigma) - m_fApplyMassScale = 0; // don't apply mass scale correction by default preparedInput.m_dataType = 1; // set to "data" by default preparedInput.m_fUseTailCleanup = 1; // cleanup by default for lep-had Moriond 2012 analysis preparedInput.m_fUseDefaults = 0; // use pre-set defaults for various configurations; if set it to 0 @@ -401,13 +399,7 @@ void MissingMassCalculatorV2::DoOutputInfo() { OutputInfo.m_nuvec2[MMCFitMethodV2::MAXW].Py()); OutputInfo.m_FittedMetVec[MMCFitMethodV2::MAXW] = metmaxw; - // MLM method : can only get MMC, rest is dummy - double scale = MassScale(MMCFitMethodV2::MAXW, m_fDitauStuffHisto.Mditau_best, - preparedInput.m_type_visTau1, - preparedInput.m_type_visTau2); // only for histo method for now. In - // practice disabled by default - - OutputInfo.m_FittedMass[MMCFitMethodV2::MLM] = scale * m_fDitauStuffHisto.Mditau_best; + OutputInfo.m_FittedMass[MMCFitMethodV2::MLM] = m_fDitauStuffHisto.Mditau_best; OutputInfo.m_FittedMassLowerError[MMCFitMethodV2::MLM] = yq[0]; OutputInfo.m_FittedMassUpperError[MMCFitMethodV2::MLM] = yq[1]; @@ -2106,63 +2098,6 @@ int MissingMassCalculatorV2::TailCleanUp(const TLorentzVector &vis1, return pass_code; } -//-------- This function applies correction to compensate for the off-set -double MissingMassCalculatorV2::MassScale(int method, double mass, const int &tau_type1, - const int &tau_type2) { - double Fscale = 1.0; - // calibration for rel16 lep-had analysis only - if (m_fApplyMassScale == 1) { - if (preparedInput.m_tauTypes == TauTypes::lh) { - if (method != 1) - return 1.0; - // float p0, p1, p2, p3; - // if(tau_type1==1 || tau_type2==1) // 1-prong tau's - // { - // p0=3.014; p1=-71.86; p2=1.018; p3=0.8912; - // if(mass>91.2) Fscale=p0/(p1+p2*mass)+p3; - // else Fscale=p0/(p1+p2*91.2)+p3; - // } - // if(tau_type1==3 || tau_type2==3) // 3-prong tau's - // { - // p0=0.4576; p1=-84.22; p2=0.9783; p3=0.9136; - // if(mass>91.2) Fscale=p0/(p1+p2*mass)+p3; - // else Fscale=p0/(p1+p2*91.2)+p3; - // } - // if(Fscale>1.0) Fscale=1.0; - // if(Fscale<0.89) Fscale=0.89; - - float p0, p1, p2, p3, p4, p5, p6, p7; - if ((tau_type1 >= 0 && tau_type1 <= 2) || (tau_type2 >= 0 && tau_type2 <= 2)) - return 1.0; // 1-prong tau's - if ((tau_type1 >= 3 && tau_type1 <= 5) || (tau_type2 >= 3 && tau_type2 <= 5)) // 3-prong tau's - { - p0 = 3.014; - p1 = -71.86; - p2 = 1.018; - p3 = 0.8912; - p4 = 0.4576; - p5 = -84.22; - p6 = 0.9783; - p7 = 0.9136; - double scale1 = p0 / (p1 + p2 * mass) + p3; - double scale3 = p4 / (p5 + p6 * mass) + p7; - if (mass > 91.2) - Fscale = scale3 / scale1; - else { - scale1 = p0 / (p1 + p2 * 91.2) + p3; - scale3 = p4 / (p5 + p6 * 91.2) + p7; - Fscale = scale3 / scale1; - } - } - if (Fscale > 1.0) - Fscale = 1.0; - if (Fscale < 0.95) - Fscale = 0.95; - } - } - return 1.0 / Fscale; -} - // note that if MarkovChain the input solutions can be modified void MissingMassCalculatorV2::handleSolutions() @@ -2935,49 +2870,6 @@ double MissingMassCalculatorV2::dTheta3DLimit(const int &tau_type, const int &li return limit; } -// returns P(nu1) & P(nu2) -int MissingMassCalculatorV2::NuPsolution(TVector2 met_vec, double theta1, double phi1, - double theta2, double phi2, double &P1, - double &P2) { - int solution_code = 0; // 0== no solution, 1==with solution - P1 = 0.0; - P2 = 0.0; - double D = sin(theta1) * sin(theta2) * sin(phi2 - phi1); - if (std::abs(D) > 0.0) // matrix deteriminant is non-zero - { - P1 = (met_vec.Px() * sin(phi2) - met_vec.Py() * cos(phi2)) * sin(theta2) / D; - P2 = (met_vec.Py() * cos(phi1) - met_vec.Px() * sin(phi1)) * sin(theta1) / D; - if (P1 > 0.0 && P2 > 0.0) - solution_code = 1; - } - return solution_code; -} - -// standard collinear approximation -// it returns code=0 if collinear approximation can't be applied -// and code=1 and Mrec if collinear approximation was applied -int MissingMassCalculatorV2::StandardCollApprox(const TLorentzVector &tau_vec1, - const TLorentzVector &tau_vec2, - const TVector2 &met_vec, double &Mrec) { - int code = 0; - Mrec = 0.0; - double P_nu1 = 0.0; - double P_nu2 = 0.0; - int coll_code = NuPsolution(met_vec, tau_vec1.Theta(), tau_vec1.Phi(), tau_vec2.Theta(), - tau_vec2.Phi(), P_nu1, P_nu2); - if (coll_code == 1) { - code = 1; - TLorentzVector nu1(P_nu1 * sin(tau_vec1.Theta()) * cos(tau_vec1.Phi()), - P_nu1 * sin(tau_vec1.Theta()) * sin(tau_vec1.Phi()), - P_nu1 * cos(tau_vec1.Theta()), P_nu1); - TLorentzVector nu2(P_nu2 * sin(tau_vec2.Theta()) * cos(tau_vec2.Phi()), - P_nu2 * sin(tau_vec2.Theta()) * sin(tau_vec2.Phi()), - P_nu2 * cos(tau_vec2.Theta()), P_nu2); - Mrec = (nu1 + nu2 + tau_vec1 + tau_vec2).M(); - } - return code; -} - // checks units of input variables, converts into [GeV] if needed, make all // possible corrections DR new : now a second structure preparedInput is derived // from the input one which only has direct user input -- GitLab From 76a8f5c23cc0d48abad6082f9a66f89457abc8b9 Mon Sep 17 00:00:00 2001 From: Antonio De Maria <antonio.de.maria@cern.ch> Date: Thu, 25 Jul 2024 16:04:15 +0200 Subject: [PATCH 2/2] further cleaning --- .../DiTauMassTools/HelperFunctions.h | 2 - .../DiTauMassTools/Root/HelperFunctions.cxx | 26 ------ .../DiTauMassTools/Root/MissingMassProb.cxx | 90 ------------------- 3 files changed, 118 deletions(-) diff --git a/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h b/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h index 3b69f2bc36fc..b8ffd41197b6 100644 --- a/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h +++ b/PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h @@ -53,8 +53,6 @@ namespace TauTypes // see source file of MissingMassProb for further reasoning template <typename T> void ignore(T &&){} -int getFirstBinBelowMax(const std::shared_ptr<TH1F>& hist, double max, double targetVal); -int getFirstBinAboveMax(const std::shared_ptr<TH1F>& hist, double max, double targetVal); double Angle(const TLorentzVector & vec1, const TLorentzVector & vec2); bool updateDouble (const double in, double & out) ; void fastSinCos (const double & phi, double & sinPhi, double & cosPhi); diff --git a/PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx b/PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx index 2a0f187adc78..a4cdf27c87c4 100644 --- a/PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx +++ b/PhysicsAnalysis/TauID/DiTauMassTools/Root/HelperFunctions.cxx @@ -12,32 +12,6 @@ using namespace DiTauMassTools; -int DiTauMassTools::getFirstBinAboveMax(const std::shared_ptr<TH1F>& hist, double max, double targetVal) -{ - int maxBin = hist->FindBin(max); - int targetBin = maxBin; - for(int i=maxBin; i<=hist->GetNbinsX(); i++){ - if(hist->GetBinContent(i)<targetVal){ - targetBin = i-1; - break; - } - } - return targetBin; -} - -int DiTauMassTools::getFirstBinBelowMax(const std::shared_ptr<TH1F>& hist, double max, double targetVal) -{ - int maxBin = hist->FindBin(max); - int targetBin = maxBin; - for(int i=1; i<=maxBin; i++){ - if(hist->GetBinContent(i)>=targetVal){ - targetBin = i; - break; - } - } - return targetBin; -} - double DiTauMassTools::MaxDelPhi(int tau_type, double Pvis, double dRmax_tau) { return dRmax_tau+0*Pvis*tau_type; // hack to avoid warning diff --git a/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassProb.cxx b/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassProb.cxx index 142a220d7ca8..9978cebb868f 100644 --- a/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassProb.cxx +++ b/PhysicsAnalysis/TauID/DiTauMassTools/Root/MissingMassProb.cxx @@ -1153,96 +1153,6 @@ double MissingMassProb::dTheta3Dparam(const int & parInd, const int & tau_type, return 0.; } -/* -// Sasha: keep this for now, may need in the future -// returns analytical P(theta)-probability for given tau-topology -// decayType==1 for leptonic decays and 0 for hadronic decays -// uses product of probabilities -double MissingMassProb::AngularProbability(TLorentzVector nu_vec, TLorentzVector vis_vec, int decayType) { - double prob=0.0; - double M=1.777; - double angl=0.0; - double P=(vis_vec+nu_vec).P(); - double V=P/sqrt(P*P+M*M); // tau speed - double dA=dRmax_tau/(2.0*Niter_fit1); - - if(decayType==1) // leptonic tau for now - { - // exact formular for energy probability is sqrt(1-V^2)/(2*V*p_0)*dE - double m_1=nu_vec.M(); - double m_2=vis_vec.M(); - double E_nu=(M*M+m_1*m_1-m_2*m_2)/(2.0*M); - if(E_nu<=m_1) return 0.0; - double P_nu=sqrt(pow(E_nu,2)-pow(m_1,2)); - double prob1=0.5*sqrt(1-V*V)/(P_nu*V); // energy probability - - angl=Angle(vis_vec,vis_vec+nu_vec); // using lepton direction - double det1=1.0-V*cos(angl+dA); - double det2= (angl-dA)>0.0 ? 1.0-V*cos(angl-dA) : 1.0-V; - double prob2=std::abs(1.0/det1-1.0/det2)*(1.0-V*V)/(2.0*V); // using massless approximation for leptons - prob=prob1*prob2; - } - if(decayType==0) // hadronic tau - { - // exact formula for energy probability is sqrt(1-V^2)/(2*V*p_0)*dE - // drop p_0 because it's a contstant for a given hadronic tau - double prob1=0.5*sqrt(1-V*V)/V; // "energy" probability - - angl=Angle(nu_vec,vis_vec+nu_vec); // using neutrino direction - double det1=1.0-V*cos(angl+dA); - double det2= (angl-dA)>0.0 ? 1.0-V*cos(angl-dA) : 1.0-V; - double prob2=std::abs(1.0/det1-1.0/det2)*(1.0-V*V)/(2.0*V); - prob=prob1*prob2; - } - return prob; -} - -// returns probability of angle between two tau's -// assuming massless tau's for now, should be small effect for M>M_Z -double MissingMassProb::ResonanceProbability(TLorentzVector vec1, TLorentzVector vec2) { - - double prob=1.0; - double boson_P=(vec1+vec2).P(); - if(boson_P==0.0) return 1.0; - double boson_E=(vec1+vec2).E(); - double boson_V=0.0; - if(boson_E>0.0) boson_V=boson_P/boson_E; - else return 1.0E-10; - - double testMass=(vec1+vec2).M(); - double m=1.777; // tau mass - double E_tau=testMass/2.0; - double P_tau=sqrt(pow(E_tau,2)-m*m); - double Evis_lim[2]; - Evis_lim[0]=(E_tau-boson_V*P_tau)/sqrt(1.0-boson_V*boson_V); - Evis_lim[1]=(E_tau+boson_V*P_tau)/sqrt(1.0-boson_V*boson_V); - if(vec1.E()<Evis_lim[0] || vec1.E()>Evis_lim[1]) return 1.0E-20; - if(vec2.E()<Evis_lim[0] || vec2.E()>Evis_lim[1]) return 1.0E-20; - - double prob1=0.5*sqrt(1-boson_V*boson_V)/(P_tau*boson_V); - - if(vec1.P()*vec2.P()>0) - { - double theta=acos((vec1.Px()*vec2.Px()+vec1.Py()*vec2.Py()+vec1.Pz()*vec2.Pz())/(vec1.P()*vec2.P())); - if(boson_V>0.0 && boson_V<1.0) - { - if(boson_V<cos(theta/2)) return 1.0E-10; - double num=(1.0-boson_V*boson_V)*cos(theta/2); - double denom=4*boson_V*sin(theta/2)*sin(theta/2)*sqrt(boson_V*boson_V-cos(theta/2)*cos(theta/2)); - prob=num/denom; - } - else - { - if(std::abs(theta-TMath::Pi())>0.0001) return 1.0E-10; - } - } - else return 1.0E-10; - prob=prob*prob1; - if(prob<1.0E-20) prob=1.0E-20; - return prob; -} -*/ - void MissingMassProb::MET(MissingMassInput& preparedInput){ // compute MET resolution (eventually use HT) if(preparedInput.m_METsigmaP<0.1 || preparedInput.m_METsigmaL<0.1) -- GitLab