MINT2
Public Member Functions | Protected Attributes | List of all members
AmpsPdfFlexiFastCPV_time_integrated Class Reference
Inheritance diagram for AmpsPdfFlexiFastCPV_time_integrated:
MINT::PdfBase< IDalitzEvent > IDalitzPdf MINT::IPdf< IDalitzEvent > MINT::IReturnRealForEvent< IDalitzEvent > MINT::IReturnRealForEvent< IDalitzEvent > MINT::IPdf< IDalitzEvent > MINT::IReturnRealForEvent< IDalitzEvent > MINT::IReturnRealForEvent< IDalitzEvent >

Public Member Functions

void parametersChanged ()
 
void beginFit ()
 
void endFit ()
 
void printIntegralVals ()
 
double un_normalised_noPs (IDalitzEvent &evt)
 
virtual double getVal (IDalitzEvent &evt)
 
virtual double getVal_withPs (IDalitzEvent &evt)
 
virtual double getVal_noPs (IDalitzEvent &evt)
 
virtual double getVal (IDalitzEvent *evt)
 
virtual double getVal_withPs (IDalitzEvent *evt)
 
virtual double getVal_noPs (IDalitzEvent *evt)
 
virtual DalitzHistoSet histoSet ()
 
void doFinalStatsAndSaveForAmp12 (MINT::Minimiser *min=0, const std::string &fname="FitAmpResults", const std::string &fnameROOT="fitFractions")
 
 AmpsPdfFlexiFastCPV_time_integrated (AmpsPdfFlexiFast *amps1, AmpsPdfFlexiFast *amps2, AmpsPdfFlexiFast *ampsSum, MINT::FitParameter &r, MINT::FitParameter &delta, MINT::FitParameter &gamma, MINT::FitParameter &tau, MINT::FitParameter &dGamma, MINT::FitParameter &dm, MINT::FitParameter &eff_tag, MINT::FitParameter &w)
 
- Public Member Functions inherited from MINT::PdfBase< IDalitzEvent >
 PdfBase ()
 
 PdfBase (const PdfBase< IDalitzEvent > &)
 
virtual double RealVal (IDalitzEvent &evt)
 
virtual double getNewVal (IDalitzEvent &evt)
 
virtual ~PdfBase ()
 
- Public Member Functions inherited from MINT::IPdf< IDalitzEvent >
virtual void Gradient (IDalitzEvent &evt, std::vector< double > &grad, MINT::MinuitParameterSet *mps)
 
virtual bool useAnalyticGradient ()
 
- Public Member Functions inherited from MINT::IReturnRealForEvent< IDalitzEvent >
virtual ~IReturnRealForEvent ()
 

Protected Attributes

AmpsPdfFlexiFast_amps1
 
AmpsPdfFlexiFast_amps2
 
AmpsPdfFlexiFast_ampsSum
 
FitParameter_r
 
FitParameter_delta
 
FitParameter_gamma
 
FitParameter_tau
 
FitParameter_dGamma
 
FitParameter_dm
 
FitParameter_eff_tag
 
FitParameter_w
 
double _intA
 
double _intAbar
 
complex< double > _intAAbar
 

Additional Inherited Members

- Protected Member Functions inherited from MINT::IReturnRealForEvent< IDalitzEvent >
 IReturnRealForEvent ()
 

Detailed Description

Definition at line 993 of file BsDsKpipi_TD_ampFit.cpp.

Constructor & Destructor Documentation

◆ AmpsPdfFlexiFastCPV_time_integrated()

AmpsPdfFlexiFastCPV_time_integrated::AmpsPdfFlexiFastCPV_time_integrated ( AmpsPdfFlexiFast amps1,
AmpsPdfFlexiFast amps2,
AmpsPdfFlexiFast ampsSum,
MINT::FitParameter r,
MINT::FitParameter delta,
MINT::FitParameter gamma,
MINT::FitParameter tau,
MINT::FitParameter dGamma,
MINT::FitParameter dm,
MINT::FitParameter eff_tag,
MINT::FitParameter w 
)
inline

Definition at line 1110 of file BsDsKpipi_TD_ampFit.cpp.

1112  :
1113  _amps1(amps1),_amps2(amps2),_ampsSum(ampsSum),_r(r),_delta(delta),_gamma(gamma),_tau(tau),_dGamma(dGamma),_dm(dm),_eff_tag(eff_tag),_w(w),
1114  _intA(-1),_intAbar(-1),_intAAbar(-1)
1115  {;}

Member Function Documentation

◆ beginFit()

void AmpsPdfFlexiFastCPV_time_integrated::beginFit ( )
inlinevirtual

◆ doFinalStatsAndSaveForAmp12()

void AmpsPdfFlexiFastCPV_time_integrated::doFinalStatsAndSaveForAmp12 ( MINT::Minimiser min = 0,
const std::string &  fname = "FitAmpResults",
const std::string &  fnameROOT = "fitFractions" 
)
inline

Definition at line 1103 of file BsDsKpipi_TD_ampFit.cpp.

1103  {
1106  _amps1->doFinalStatsAndSave(min,((string)fname+".txt").c_str(),((string)fnameROOT+"_TI.root").c_str());
1107  _amps2->doFinalStatsAndSave(min,((string)fname+"_CC.txt").c_str(),((string)fnameROOT+"_TI_CC.root").c_str());
1108  }
void doFinalStatsAndSave(MINT::Minimiser *min=0, const std::string &fname="FitAmpResults.txt", const std::string &fnameROOT="fitFractions.root")

◆ endFit()

void AmpsPdfFlexiFastCPV_time_integrated::endFit ( )
inlinevirtual

◆ getVal() [1/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1066 of file BsDsKpipi_TD_ampFit.cpp.

1066  {
1067 
1068  const double f = static_cast<double>((int)evt.getValueFromVector(4));
1069  const double val = un_normalised_noPs(evt);
1070 
1071  double r = (double)_r; // * sqrt(_intA/_intAbar);
1072  const double Gamma = 1./((double) _tau);
1073  const complex<double> phase_diff = polar((double)r,((double) _delta -(double)_gamma*f)/360.*2*pi);
1074  const double int_interference = (phase_diff*_intAAbar).real();
1075 
1076  if(_intA == -1 ){
1077  cout << "AmpsPdfFlexiFastCPV:: _norm = -1, should not have happened." << endl;
1078  throw "can't deal with that";
1079  }
1080 
1081  return val/(2.* ((_intA + r* r * _intAbar) * 4.*Gamma - int_interference * 2. * _dGamma )/ (4.*Gamma*Gamma-_dGamma*_dGamma));
1082  }
static const double pi
virtual double getValueFromVector(unsigned int i) const =0
double un_normalised_noPs(IDalitzEvent &evt)

◆ getVal() [2/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1087 of file BsDsKpipi_TD_ampFit.cpp.

1087  {
1088  if(0 == evt) return 0;
1089  return getVal(*evt);
1090  }
virtual double getVal(IDalitzEvent &evt)

◆ getVal_noPs() [1/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal_noPs ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1085 of file BsDsKpipi_TD_ampFit.cpp.

1085 {return getVal(evt);}
virtual double getVal(IDalitzEvent &evt)

◆ getVal_noPs() [2/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal_noPs ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1095 of file BsDsKpipi_TD_ampFit.cpp.

1095  {
1096  if(0 == evt) return 0;
1097  return getVal_noPs(*evt);
1098  }
virtual double getVal_noPs(IDalitzEvent &evt)

◆ getVal_withPs() [1/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal_withPs ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1084 of file BsDsKpipi_TD_ampFit.cpp.

1084 {return getVal(evt);}
virtual double getVal(IDalitzEvent &evt)

◆ getVal_withPs() [2/2]

virtual double AmpsPdfFlexiFastCPV_time_integrated::getVal_withPs ( IDalitzEvent evt)
inlinevirtual

Implements IDalitzPdf.

Definition at line 1091 of file BsDsKpipi_TD_ampFit.cpp.

1091  {
1092  if(0 == evt) return 0;
1093  return getVal_withPs(*evt);
1094  }
virtual double getVal_withPs(IDalitzEvent &evt)

◆ histoSet()

virtual DalitzHistoSet AmpsPdfFlexiFastCPV_time_integrated::histoSet ( )
inlinevirtual

Implements IDalitzPdf.

Definition at line 1100 of file BsDsKpipi_TD_ampFit.cpp.

1100 {return _ampsSum->histoSet();}
virtual DalitzHistoSet histoSet() const

◆ parametersChanged()

void AmpsPdfFlexiFastCPV_time_integrated::parametersChanged ( )
inlinevirtual

◆ printIntegralVals()

void AmpsPdfFlexiFastCPV_time_integrated::printIntegralVals ( )
inline

Definition at line 1031 of file BsDsKpipi_TD_ampFit.cpp.

1031  {
1032  cout << "intSum = " << _ampsSum->getIntegralValue() << endl;
1033  cout << "intA = " << _intA << endl;
1034  cout << "intAbar = " << _intAbar << endl;
1035  cout << "intAAbar = " << _intAAbar.real() << endl;
1036  }

◆ un_normalised_noPs()

double AmpsPdfFlexiFastCPV_time_integrated::un_normalised_noPs ( IDalitzEvent evt)
inline

Definition at line 1038 of file BsDsKpipi_TD_ampFit.cpp.

1038  {
1039  const double t = (double) evt.getValueFromVector(0);
1040  const double dt = (double) evt.getValueFromVector(1);
1041  const double q = static_cast<double>((int)evt.getValueFromVector(2));
1042  const double w = (double) evt.getValueFromVector(3);
1043  const double f = static_cast<double>((int)evt.getValueFromVector(4));
1044 
1045  const double e_eff = fabs(q)*_eff_tag + (1.-fabs(q))*(1.-_eff_tag);
1046 
1047  double r = (double)_r; // * sqrt(_intA/_intAbar);
1048  const complex<double> phase_diff = polar((double)r,((double) _delta -(double)_gamma*f)/360.*2*pi);
1049 
1050  const std::complex<double> amp = _amps1->ComplexVal_un_normalised_noPs(evt) ;
1051  const std::complex<double> amp_bar = _amps2->ComplexVal_un_normalised_noPs(evt) * phase_diff;
1052 
1053  double Gamma = 1./((double) _tau);
1054 
1055  const double val =
1056  (
1057  (2.-fabs(q))*(norm(amp) + norm(amp_bar))*(4.*Gamma/(4.*Gamma*Gamma-_dGamma*_dGamma))
1058  +f*q*(1.-2.*w)*(norm(amp) - norm(amp_bar)) *(Gamma/(Gamma*Gamma+_dm*_dm))
1059  -(2.-fabs(q))*2.0*real(amp_bar*conj(amp))*(2.*_dGamma/(4.*Gamma*Gamma-_dGamma*_dGamma))
1060  -f*2.0*q*(1.-2.*w)*imag(amp_bar*conj(amp))*(_dm/(Gamma*Gamma+_dm*_dm))
1061  )*e_eff;
1062 
1063  return val;
1064  }
static const double pi
virtual double getValueFromVector(unsigned int i) const =0
std::complex< double > ComplexVal_un_normalised_noPs(IDalitzEvent &evt)

Member Data Documentation

◆ _amps1

AmpsPdfFlexiFast* AmpsPdfFlexiFastCPV_time_integrated::_amps1
protected

Definition at line 997 of file BsDsKpipi_TD_ampFit.cpp.

◆ _amps2

AmpsPdfFlexiFast* AmpsPdfFlexiFastCPV_time_integrated::_amps2
protected

Definition at line 998 of file BsDsKpipi_TD_ampFit.cpp.

◆ _ampsSum

AmpsPdfFlexiFast* AmpsPdfFlexiFastCPV_time_integrated::_ampsSum
protected

Definition at line 999 of file BsDsKpipi_TD_ampFit.cpp.

◆ _delta

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_delta
protected

Definition at line 1002 of file BsDsKpipi_TD_ampFit.cpp.

◆ _dGamma

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_dGamma
protected

Definition at line 1006 of file BsDsKpipi_TD_ampFit.cpp.

◆ _dm

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_dm
protected

Definition at line 1007 of file BsDsKpipi_TD_ampFit.cpp.

◆ _eff_tag

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_eff_tag
protected

Definition at line 1008 of file BsDsKpipi_TD_ampFit.cpp.

◆ _gamma

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_gamma
protected

Definition at line 1003 of file BsDsKpipi_TD_ampFit.cpp.

◆ _intA

double AmpsPdfFlexiFastCPV_time_integrated::_intA
protected

Definition at line 1011 of file BsDsKpipi_TD_ampFit.cpp.

◆ _intAAbar

complex<double> AmpsPdfFlexiFastCPV_time_integrated::_intAAbar
protected

Definition at line 1013 of file BsDsKpipi_TD_ampFit.cpp.

◆ _intAbar

double AmpsPdfFlexiFastCPV_time_integrated::_intAbar
protected

Definition at line 1012 of file BsDsKpipi_TD_ampFit.cpp.

◆ _r

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_r
protected

Definition at line 1001 of file BsDsKpipi_TD_ampFit.cpp.

◆ _tau

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_tau
protected

Definition at line 1005 of file BsDsKpipi_TD_ampFit.cpp.

◆ _w

FitParameter& AmpsPdfFlexiFastCPV_time_integrated::_w
protected

Definition at line 1009 of file BsDsKpipi_TD_ampFit.cpp.


The documentation for this class was generated from the following file: