diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKHist.cc b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKHist.cc index 018d98522550755fbcb610d44fa1f5c0f004e0ad..0d66542e46b57c5cff1abce12c5bbb90ef51736e 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKHist.cc +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKHist.cc @@ -354,7 +354,18 @@ void JpsiGamKsKlKKHist::fillTuple( TNtuple* theTuple, EvtData* theData, double Vector4<double>& V4_Kp_HeliPsi= theData->FourVecs[enumJpsiGamKsKlKKData::V4_Kp_HeliPsi] ; Vector4<double>& V4_Km_HeliPsi= theData->FourVecs[enumJpsiGamKsKlKKData::V4_Km_HeliPsi] ; - double thePhiPhiDecayPlaneAngle = decayAngleChi( V4_KsKlKpKm_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi, V4_Ks_HeliPsi, V4_Kl_HeliPsi ); + + Vector4<double>& V4_normKpKmDecHeliKsKlKpKm = theData->FourVecs[enumJpsiGamKsKlKKData::V4_normKpKmDecHeliKsKlKpKm]; + Vector4<double>& V4_normKsKlDecHeliKsKlKpKm = theData->FourVecs[enumJpsiGamKsKlKKData::V4_normKsKlDecHeliKsKlKpKm]; + double cosChi=(V4_normKpKmDecHeliKsKlKpKm.Px()*V4_normKsKlDecHeliKsKlKpKm.Px() + +V4_normKpKmDecHeliKsKlKpKm.Py()*V4_normKsKlDecHeliKsKlKpKm.Py() + +V4_normKpKmDecHeliKsKlKpKm.Pz()*V4_normKsKlDecHeliKsKlKpKm.Pz()) + / (V4_normKpKmDecHeliKsKlKpKm.P()*V4_normKsKlDecHeliKsKlKpKm.P()); + + double chi=acos(fabs(cosChi)); + double thePhiPhiDecayPlaneAngle = chi/TMath::Pi()*180.; + + //double thePhiPhiDecayPlaneAngle = decayAngleChi( V4_KsKlKpKm_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi, V4_Ks_HeliPsi, V4_Kl_HeliPsi ); double testHeli = costDecHeli( V4_KsKlKpKm_HeliPsi+V4_gamma_HeliPsi, V4_Ks_HeliPsi+V4_Kl_HeliPsi+V4_Km_HeliPsi+V4_Kp_HeliPsi, V4_Ks_HeliPsi+V4_Kl_HeliPsi ); theTuple->Fill( V4_KsKlKpKm_HeliPsi.M(), diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.cc b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.cc index c1b965dab4d12902f80a79ce157ad786e02e7d5e..b5a48c1849bde4dc854333e0c766dfa5eda28396 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.cc +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.cc @@ -39,15 +39,19 @@ bool JpsiGamKsKlKKParser::parseCommandLine(int argc, char **argv) ("paramFile",po::value<string>(&_paramFile), "file with start parameters for fit or QA (full path)") ("startHypo",po::value<string>(&_startHypo), "choose the hyopthesis to start") ("enableHyp",po::value< vector<string> >(&_enabledHyps), "enable hypotheses") - ("mode",po::value<string>(&_mode), "enable/diable QA mode") + ("mode",po::value<string>(&_mode), "modes are: pwa, dumpDefaulParams, qaMode") + ("massIndependentFit", po::value<bool>(&_massIndependentFit), "enable/disable mass independence in fit") + ("commonProdPhases",po::value<bool>(&_useCommonProductionPhases), "enable/disable common production phases") ; - + po::options_description config("Configuration file options"); config.add_options() ("verbose",po::value<bool>(&verbose)->default_value(true), "Determines whether additional information should be emitted") ("mnParFix",po::value< vector<string> >(&_mnParFixs), "minuit parameters can be fixed here") + ("massRangeMin",po::value<double>(&_massMin), "min of phi phi mass range for mass indep. fit") + ("massRangeMax",po::value<double>(&_massMax), "max of phi phi mass range for mass indep. fit") ; - + po::options_description cmdline_options; cmdline_options.add(desc).add(common); diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.hh b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.hh index e102f8349031acd30c74866f11bebe1e54726348..12b4db7ab20601fa9ba1df135ce7dc8cf3416f0d 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.hh +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKParser.hh @@ -6,6 +6,7 @@ #include <sstream> #include <string> #include <vector> +#include <utility> // Boost headers go here @@ -45,7 +46,11 @@ class JpsiGamKsKlKKParser , _paramFile("/data/sleipnir1/bertram/JpsiGamKsKlKKData/startParamSpin02.dat") , _startHypo("base") , _mode("qaMode") - { + , _massIndependentFit(false) + , _useCommonProductionPhases(false) + , _massMin(0.0) + , _massMax(10.) + { if (!parseCommandLine(argc, argv)) throw false; } @@ -57,7 +62,12 @@ class JpsiGamKsKlKKParser const std::vector<std::string>& enabledHyps() const { return _enabledHyps; } const std::string startHypo() const {return _startHypo;} const std::string mode() const {return _mode;} - const std::vector<std::string>& fixedParams() const { return _mnParFixs; } + const std::vector<std::string>& fixedParams() const { return _mnParFixs; } + const bool massIndependentFit() const {return _massIndependentFit; } + const bool useCommonProductionPhases() const {return _useCommonProductionPhases; } + const std::pair<double, double> massRange() const { return std::make_pair( _massMin, _massMax ) ; } + + protected: bool parseCommandLine(int argc,char **argv); @@ -71,6 +81,12 @@ protected: std::string _mode; std::vector<std::string> _enabledHyps; std::vector<std::string> _mnParFixs; + bool _massIndependentFit; + bool _useCommonProductionPhases; + + double _massMin; + double _massMax; + }; diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.cc b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.cc index 6afb35c7782ed7aabb98330ece9bf2b6a99e032f..c5897d6ba4713072870c0eac5f6a37d7880bd8ff 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.cc +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.cc @@ -17,6 +17,10 @@ JpsiGamKsKlKKProdLh::JpsiGamKsKlKKProdLh(boost::shared_ptr<const EvtDataBaseList ,_eta21870Hyp(false) ,_f1Hyp(false) ,_usePhasespace(false) + ,_useCommonProductionPhase(true) + ,_massIndependentFit(false) + ,_phiMass( 1.019455) + ,_kaonMass(0.493677) { @@ -37,6 +41,10 @@ JpsiGamKsKlKKProdLh::JpsiGamKsKlKKProdLh( boost::shared_ptr<AbsLh> theLhPtr, con ,_eta21870Hyp(false) ,_f1Hyp(false) ,_usePhasespace(false) + ,_useCommonProductionPhase(true) + ,_massIndependentFit(false) + ,_phiMass( 1.019455) + ,_kaonMass(0.493677) { initializeHypothesisMap( hypMap); @@ -52,31 +60,60 @@ JpsiGamKsKlKKProdLh::~JpsiGamKsKlKKProdLh() double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& theParamVal){ double result=0.; + + Vector4<double> fv2Phi = theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; + complex<double> JmpGmp(0.,0.); + complex<double> JmpGmm(0.,0.); + complex<double> JmmGmp(0.,0.); + complex<double> JmmGmm(0.,0.); + - Vector4<double> fvKsKlKpKm = theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; - double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::etac]; - double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::etac]; + if(_massIndependentFit){ + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > etaMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > etaPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; + JmpGmp+=etaGammaAmp(1, 0, 1, theData, etaMag, etaPhi ); + JmpGmm+=etaGammaAmp(1, 0, -1, theData, etaMag, etaPhi ); + JmmGmp+=etaGammaAmp(-1, 0, 1, theData, etaMag, etaPhi ); + JmmGmm+=etaGammaAmp(-1, 0, -1, theData, etaMag, etaPhi ); + + result=norm(JmpGmp)+norm(JmpGmm)+norm(JmmGmp)+norm(JmmGmm); + return result; + } - complex<double> JmpGmp=etaGammaAmp(1, 0, 1, theData, PsiToEtacGamMag, PsiToEtacGamPhi, mass, width ); - complex<double> JmpGmm=etaGammaAmp(1, 0, -1, theData, PsiToEtacGamMag, PsiToEtacGamPhi, mass, width ); - complex<double> JmmGmp=etaGammaAmp(-1, 0, 1, theData, PsiToEtacGamMag, PsiToEtacGamPhi, mass, width ); - complex<double> JmmGmm=etaGammaAmp(-1, 0, -1, theData, PsiToEtacGamMag, PsiToEtacGamPhi, mass, width ); + + if(_etacHyp){ + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::PsiToEtacGamma]; + + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::etac]; + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::etac]; + complex<double> theDynamicPart = BreitWigner(fv2Phi, mass, width); + //complex<double> theDynamicPart = BreitWignerBlattW(fv2Phi, _phiMass, _phiMass, mass, width, 1); + + JmpGmp+=etaGammaAmp(1, 0, 1, theData, PsiToEtacGamMag, PsiToEtacGamPhi )*theDynamicPart; + JmpGmm+=etaGammaAmp(1, 0, -1, theData, PsiToEtacGamMag, PsiToEtacGamPhi )*theDynamicPart; + JmmGmp+=etaGammaAmp(-1, 0, 1, theData, PsiToEtacGamMag, PsiToEtacGamPhi )*theDynamicPart; + JmmGmm+=etaGammaAmp(-1, 0, -1, theData, PsiToEtacGamMag, PsiToEtacGamPhi )*theDynamicPart; + + } + if (_eta2225Hyp){ std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEta2225GamMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToEta2225Gamma]; std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEta2225GamPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::PsiToEta2225Gamma]; - mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::eta2225]; - width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::eta2225]; + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::eta2225]; + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::eta2225]; + complex<double> theDynamicPart = BreitWigner(fv2Phi, mass, width); + //complex<double> theDynamicPart = BreitWignerBlattW(fv2Phi, _phiMass, _phiMass, mass, width, 1); + + JmpGmp+=etaGammaAmp(1, 0, 1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi )*theDynamicPart; + JmpGmm+=etaGammaAmp(1, 0, -1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi )*theDynamicPart; + JmmGmp+=etaGammaAmp(-1, 0, 1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi )*theDynamicPart; + JmmGmm+=etaGammaAmp(-1, 0, -1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi )*theDynamicPart; - JmpGmp+=etaGammaAmp(1, 0, 1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi, mass, width ); - JmpGmm+=etaGammaAmp(1, 0, -1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi, mass, width ); - JmmGmp+=etaGammaAmp(-1, 0, 1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi, mass, width ); - JmmGmm+=etaGammaAmp(-1, 0, -1, theData, PsiToEta2225GamMag, PsiToEta2225GamPhi, mass, width ); } if (_f02020Hyp || _f02020FlatteHyp ){ @@ -86,34 +123,24 @@ double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& thePar std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F02020ToPhiPhiMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::F02020ToPhiPhi]; std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F02020ToPhiPhiPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::F02020ToPhiPhi]; - - complex<double> dynamicPart(0.,0.); - Vector4<double> fv2Phi = theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; - mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f02020]; + dynamicModelParams theDynModel; + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f02020]; + if( _f02020FlatteHyp ){ double gKK = theParamVal.gFactors[paramEnumJpsiGamKsKlKK::f02020gKK]; double gPhiPhi = theParamVal.gFactors[paramEnumJpsiGamKsKlKK::f02020gPhiPhi]; - - const double phiMass = 1.019455; - const double kaonMass = 0.493677; - std::pair <const double, const double> kkPair=make_pair(kaonMass, kaonMass); - std::pair <const double, const double> phiphiPair=make_pair(phiMass, phiMass); - dynamicPart = Flatte( fv2Phi , phiphiPair, kkPair, mass, gPhiPhi, gKK ); + initializeFlatteModel( theDynModel, fv2Phi, mass, gPhiPhi, gKK ); - } else{ - width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f02020]; - dynamicPart = BreitWigner( fv2Phi, mass, width); + }else if( _f02020Hyp ){ + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f02020]; + initializeBreitWignerModel( theDynModel, fv2Phi, mass, width, dynamicModelParams::BreitWigner ); } - bool useCommonProductionPhase=true; - - JmpGmp+=f0GammaAmp(1, 1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,useCommonProductionPhase ) * dynamicPart; - JmpGmm+=f0GammaAmp(1, -1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,useCommonProductionPhase ) * dynamicPart; - JmmGmp+=f0GammaAmp(-1, 1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,useCommonProductionPhase ) * dynamicPart; - JmmGmm+=f0GammaAmp(-1, -1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,useCommonProductionPhase ) * dynamicPart; - - - } + JmpGmp+=f0GammaAmp(1, 1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,_useCommonProductionPhase, theDynModel ); + JmpGmm+=f0GammaAmp(1, -1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,_useCommonProductionPhase, theDynModel ); + JmmGmp+=f0GammaAmp(-1, 1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,_useCommonProductionPhase, theDynModel ); + JmmGmm+=f0GammaAmp(-1, -1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,_useCommonProductionPhase, theDynModel ); + } if (_f22300Hyp){ std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiTof2GamMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToF22300Gamma]; @@ -122,15 +149,18 @@ double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& thePar std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F2ToPhiPhiMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::F22300ToPhiPhi]; std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F2ToPhiPhiPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::F22300ToPhiPhi]; - mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f22300]; - width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f22300]; + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f22300]; + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f22300]; Spin f2Spin=2; - bool useCommonProductionPhase=true; - JmpGmp+=f2GammaAmp(1, 1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmpGmm+=f2GammaAmp(1, -1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmp+=f2GammaAmp(-1, 1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmm+=f2GammaAmp(-1, -1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); + dynamicModelParams theDynModel; + initializeBreitWignerModel( theDynModel, fv2Phi, mass, width, dynamicModelParams::BreitWigner ); + + + JmpGmp+=f2GammaAmp(1, 1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmpGmm+=f2GammaAmp(1, -1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmmGmp+=f2GammaAmp(-1, 1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmmGmm+=f2GammaAmp(-1, -1, f2Spin, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); } @@ -141,16 +171,17 @@ double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& thePar std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > Eta2ToPhiPhiMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::Eta21870ToPhiPhi]; std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > Eta2ToPhiPhiPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::Eta21870ToPhiPhi]; - mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::eta21870]; - width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::eta21870]; + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::eta21870]; + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::eta21870]; Spin eta2Spin=2; - bool useCommonProductionPhase=true; + dynamicModelParams theDynModel; + initializeBreitWignerModel( theDynModel, fv2Phi, mass, width, dynamicModelParams::BreitWigner ); - JmpGmp+=f2GammaAmp(1, 1, eta2Spin, theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmpGmm+=f2GammaAmp(1, -1, eta2Spin,theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmp+=f2GammaAmp(-1, 1, eta2Spin,theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmm+=f2GammaAmp(-1, -1, eta2Spin, theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width , useCommonProductionPhase); + JmpGmp+=f2GammaAmp(1, 1, eta2Spin, theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, _useCommonProductionPhase,theDynModel ); + JmpGmm+=f2GammaAmp(1, -1, eta2Spin,theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, _useCommonProductionPhase,theDynModel ); + JmmGmp+=f2GammaAmp(-1, 1, eta2Spin,theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width, _useCommonProductionPhase,theDynModel ); + JmmGmm+=f2GammaAmp(-1, -1, eta2Spin, theData, PsiToEta2GamMag, PsiToEta2GamPhi,Eta2ToPhiPhiMag,Eta2ToPhiPhiPhi,mass,width , _useCommonProductionPhase,theDynModel ); } @@ -162,15 +193,17 @@ double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& thePar std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F1ToPhiPhiMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::F1ToPhiPhi]; std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > F1ToPhiPhiPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::F1ToPhiPhi]; - mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f1]; - width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f1]; + double mass = theParamVal.Masses[paramEnumJpsiGamKsKlKK::f1]; + double width = theParamVal.Widths[paramEnumJpsiGamKsKlKK::f1]; Spin f1Spin=1; - bool useCommonProductionPhase=true; - JmpGmp+=f2GammaAmp(1, 1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmpGmm+=f2GammaAmp(1, -1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmp+=f2GammaAmp(-1, 1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, useCommonProductionPhase ); - JmmGmm+=f2GammaAmp(-1, -1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, useCommonProductionPhase ); + dynamicModelParams theDynModel; + initializeBreitWignerModel( theDynModel, fv2Phi, mass, width, dynamicModelParams::BreitWigner ); + + JmpGmp+=f2GammaAmp(1, 1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmpGmm+=f2GammaAmp(1, -1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmmGmp+=f2GammaAmp(-1, 1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); + JmmGmm+=f2GammaAmp(-1, -1, f1Spin, theData, PsiTof1GamMag, PsiTof1GamPhi,F1ToPhiPhiMag,F1ToPhiPhiPhi,mass,width, _useCommonProductionPhase, theDynModel ); } @@ -193,8 +226,7 @@ complex<double> JpsiGamKsKlKKProdLh::calcCoherentAmp(Spin Minit, Spin lamGam, fi complex<double> JpsiGamKsKlKKProdLh::etaGammaAmp(Spin Minit, Spin Metac, Spin Mgamma, EvtData* theData, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaMag, - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaPhi, - double mass, double width){ + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaPhi ){ complex<double> result(0.,0.); @@ -218,7 +250,7 @@ complex<double> JpsiGamKsKlKKProdLh::etaGammaAmp(Spin Minit, Spin Metac, Spin Mg result+= amp; } - result*=BreitWigner( fv2Phi, mass, width)*etaToPhiPhiTo4KAmp(theData); + result*=etaToPhiPhiTo4KAmp(theData); return result; } @@ -234,7 +266,7 @@ complex<double> JpsiGamKsKlKKProdLh::f0GammaAmp(Spin Minit, Spin Mgamma, EvtData std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0ProdPhi, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecPhi, - bool useCommonProductionPhase){ + bool useCommonProductionPhase, dynamicModelParams& dynModPars){ complex<double> result(0.,0.); Spin f0Spin=0; @@ -271,7 +303,7 @@ complex<double> JpsiGamKsKlKKProdLh::f0GammaAmp(Spin Minit, Spin Mgamma, EvtData } //result*=BreitWigner( fv2Phi, mass, width)*f0ToPhiPhiTo4KAmp(theData, ampf0DecMag,ampf0DecPhi ); - result*=f0ToPhiPhiTo4KAmp(theData, ampf0DecMag,ampf0DecPhi ); + result*=f0ToPhiPhiTo4KAmp(theData, ampf0DecMag,ampf0DecPhi, dynModPars ); return result; @@ -294,25 +326,37 @@ complex<double> JpsiGamKsKlKKProdLh::etaToPhiPhiTo4KAmp(EvtData* theData){ } complex<double> JpsiGamKsKlKKProdLh::f0ToPhiPhiTo4KAmp(EvtData* theData, - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f0DecMag , - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f0DecPhi ){ + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f0DecMag , + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f0DecPhi, + dynamicModelParams& dynModPars){ complex<double> result(0.,0.); - + complex<double> dynamicPart(0.,0.); + Vector4<double> fv2Phi = theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; + if(dynModPars.dynamicModel == dynamicModelParams::Flatte || dynModPars.dynamicModel==dynamicModelParams::BreitWigner){ + dynamicPart = dynModPars.value; + }else if( dynModPars.dynamicModel==dynamicModelParams::MassIndependent ){ + dynamicPart = complex<double>(1,0); + } + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itf0; for ( itf0=ampf0DecMag.begin(); itf0!=ampf0DecMag.end(); ++itf0){ boost::shared_ptr<const JPCLS> f0ToPhiPhiState=itf0->first; double theMag=itf0->second; double thePhi=ampf0DecPhi[f0ToPhiPhiState]; complex<double> expiphi(cos(thePhi), sin(thePhi)); - + complex<double> tmp2PhiDecAmp(0.,0.); for (Spin lambdaPhi1=-1; lambdaPhi1<=1; lambdaPhi1++){ Spin lambdaPhi2=lambdaPhi1; tmp2PhiDecAmp+=Clebsch(1, lambdaPhi1, 1, -lambdaPhi2, f0ToPhiPhiState->S, lambdaPhi1-lambdaPhi2) *3.*conj(theData->WignerDs[enumJpsiGamKsKlKKData::Df_KsKl][1][lambdaPhi1][0])*conj(theData->WignerDs[enumJpsiGamKsKlKKData::Df_KpKm][1][lambdaPhi2][0]); //3=sqrt(2L+1)*sqrt(2L+1) Cls=1 } - - result+= theMag*expiphi*sqrt(2*f0ToPhiPhiState->L+1)*tmp2PhiDecAmp; + + if(dynModPars.dynamicModel==dynamicModelParams::BreitWignerBlattW){ + dynamicPart=BreitWignerBlattW(fv2Phi, _phiMass, _phiMass, dynModPars.mass, dynModPars.width, f0ToPhiPhiState->L ); + } + + result+= theMag*expiphi*sqrt(2*f0ToPhiPhiState->L+1)*tmp2PhiDecAmp*dynamicPart; } return result; } @@ -320,20 +364,19 @@ complex<double> JpsiGamKsKlKKProdLh::f0ToPhiPhiTo4KAmp(EvtData* theData, - complex<double> JpsiGamKsKlKKProdLh::f2GammaAmp(Spin Minit, Spin Mgamma, Spin fJSpin, EvtData* theData, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2ProdMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2ProdPhi, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2DecMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2DecPhi, - double mass, double width, bool useCommonProductionPhase ){ + double mass, double width, bool useCommonProductionPhase, dynamicModelParams& dynModPars){ complex<double> result(0.,0.); Vector4<double> fv2Phi= theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; std::map<Spin, complex<double> > decAmp; for (Spin lambdaf2=-fJSpin; lambdaf2<=fJSpin; ++lambdaf2){ - decAmp[lambdaf2]=f2ToPhiPhiTo4KAmp(theData, lambdaf2, ampf2DecMag,ampf2DecPhi ); + decAmp[lambdaf2]=f2ToPhiPhiTo4KAmp(theData, lambdaf2, ampf2DecMag,ampf2DecPhi, dynModPars ); } double theCommonPhasePhi=999.; @@ -371,7 +414,6 @@ complex<double> JpsiGamKsKlKKProdLh::f2GammaAmp(Spin Minit, Spin Mgamma, Spin fJ } } - result*=BreitWigner( fv2Phi, mass, width); return result; } @@ -379,11 +421,19 @@ complex<double> JpsiGamKsKlKKProdLh::f2GammaAmp(Spin Minit, Spin Mgamma, Spin fJ complex<double> JpsiGamKsKlKKProdLh::f2ToPhiPhiTo4KAmp( EvtData* theData, Spin f2Lambda, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecMag , - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecPhi ){ + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecPhi, + dynamicModelParams& dynModPars ){ complex<double> result(0.,0.); - - + + complex<double> dynamicPart(0.,0.); + Vector4<double> fv2Phi = theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi]; + if(dynModPars.dynamicModel == dynamicModelParams::Flatte || dynModPars.dynamicModel==dynamicModelParams::BreitWigner){ + dynamicPart = dynModPars.value; + }else if( dynModPars.dynamicModel==dynamicModelParams::MassIndependent ){ + dynamicPart = complex<double>(1,0); + } + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itf2; for ( itf2=ampf2DecMag.begin(); itf2!=ampf2DecMag.end(); ++itf2){ boost::shared_ptr<const JPCLS> f2State=itf2->first; @@ -401,7 +451,11 @@ complex<double> JpsiGamKsKlKKProdLh::f2ToPhiPhiTo4KAmp( EvtData* theData, Spin f *Clebsch(1, lambdaPhi1, 1, -lambdaPhi2, f2State->S, lambda ) *conj( theData->WignerDs[enumJpsiGamKsKlKKData::Df_Spin2][f2State->J][f2Lambda][lambda] ); - amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 ); + if(dynModPars.dynamicModel==dynamicModelParams::BreitWignerBlattW){ + dynamicPart=BreitWignerBlattW(fv2Phi, _phiMass, _phiMass, dynModPars.mass, dynModPars.width, f2State->L ); + } + + amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 ) * dynamicPart; result +=amp; } @@ -424,6 +478,45 @@ complex<double> JpsiGamKsKlKKProdLh::phiphiTo4KAmp( EvtData* theData, Spin lambd + +bool JpsiGamKsKlKKProdLh::initializeFlatteModel( dynamicModelParams &theDynModel, const Vector4<double> &fv2Phi, double mass, double gPhiPhi, double gKK ){ + + theDynModel.dynamicModel = dynamicModelParams::Flatte; + theDynModel.gFactor1 = gPhiPhi; + theDynModel.gFactor2 = gKK; + pair<const double, const double> mp1 = make_pair(_phiMass, _phiMass); + pair<const double, const double> mp2 = make_pair(_kaonMass, _kaonMass); + theDynModel.massPair1 = mp1; + theDynModel.massPair2 = mp2; + theDynModel.value = Flatte( fv2Phi, mp1, mp2, theDynModel.mass, theDynModel.gFactor1, theDynModel.gFactor2 ); + + return true; +} + + +bool JpsiGamKsKlKKProdLh::initializeBreitWignerModel(dynamicModelParams &theDynModel,const Vector4<double> &fv2Phi, double mass, double width, dynamicModelParams::enumDynamicModel theModel ){ + theDynModel.mass = mass; + theDynModel.width=width; + theDynModel.dynamicModel = theModel; + if(theModel == dynamicModelParams::BreitWigner){ + theDynModel.value = BreitWigner(fv2Phi, mass, width ); + } + return true; +} + + + + + + + + + + + + + + void JpsiGamKsKlKKProdLh::print(std::ostream& os) const{ os << "JpsiGamKsKlKKProdLh::print\n"; } @@ -623,8 +716,6 @@ JpsiGamKsKlKKProdLh::initializeHypothesisMap( const std::map<const std::string, } else Alert << "using phasespace not set!!!" <<endmsg; - - return true; } diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.hh b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.hh index 8d8c8a6830c753b7d3ff855e8afd0f738c3202c8..f08d82d6c73d28eb76456d55163e88344204879b 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.hh +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKProdLh.hh @@ -16,12 +16,14 @@ #include "PwaUtils/AbsLh.hh" #include "Examples/JpsiGamKsKlKK/JpsiGamKsKlKKData.hh" +#include "Examples/JpsiGamKsKlKK/JpsiGamKsKlKKDynamicModelParams.hh" #include "PwaUtils/DataUtils.hh" #include "Minuit2/MnUserParameters.h" + class JpsiGamKsKlKKProdLh : public AbsLh{ public: @@ -47,7 +49,9 @@ public: virtual void getDefaultParams(fitParams& fitVal, fitParams& fitErr); virtual void print(std::ostream& os) const; - + void useCommonProductionPhase( bool commonPhase ){_useCommonProductionPhase=commonPhase;} + void massIndependentFit( bool massIndep){ _massIndependentFit=massIndep; } + protected: virtual complex<double> calcCoherentAmp(Spin Minit, Spin lamGam, fitParams& theParamVal, EvtData* theData); @@ -55,20 +59,19 @@ protected: virtual complex<double> etaGammaAmp(Spin Minit, Spin Metac, Spin Mgamma, EvtData* theData, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaMag, - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaPhi, - double mass, double width); + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& PsiToEtaPhi ); virtual complex<double> f0GammaAmp(Spin Minit, Spin Mgamma, EvtData* theData, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0ProdMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0ProdPhi, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecPhi, - bool useCommonProductionPhase ); + bool useCommonProductionPhase, dynamicModelParams& dynModPars ); virtual complex<double> etaToPhiPhiTo4KAmp(EvtData* theData); virtual complex<double> f0ToPhiPhiTo4KAmp( EvtData* theData, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecMag, - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecPhi ); + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf0DecPhi, dynamicModelParams& dynModPars ); virtual complex<double> f2GammaAmp(Spin Minit, Spin Mgamma, Spin fJSpin, EvtData* theData, @@ -76,15 +79,17 @@ protected: std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2ProdPhi, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2DecMag, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >& ampf2DecPhi, - double mass, double width, bool useCommonProductionPhase ); + double mass, double width, bool useCommonProductionPhase, dynamicModelParams& dynModPars ); complex<double> f2ToPhiPhiTo4KAmp( EvtData* theData, Spin f2Lambda, std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecMag , - std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecPhi ); + std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &f2DecPhi, + dynamicModelParams& dynModPars); complex<double> phiphiTo4KAmp( EvtData* theData, Spin lambdaPhi1, Spin lambdaPhi2 ); - + bool initializeFlatteModel( dynamicModelParams &theDynModel, const Vector4<double> &fv2Phi, double mass, double gPhiPhi, double gKK ); + bool initializeBreitWignerModel(dynamicModelParams &theDynModel,const Vector4<double> &fv2Phi, double mass, double width, dynamicModelParams::enumDynamicModel theModel ); bool _etacHyp; bool _eta2225Hyp; @@ -96,11 +101,22 @@ protected: bool _f1Hyp; bool _usePhasespace; + bool _useCommonProductionPhase; + bool _massIndependentFit; std::map<const std::string, bool> _hypMap; + + + + private: bool initializeHypothesisMap( const std::map<const std::string, bool>& hypMap ); + + const double _phiMass;// = 1.019455; + const double _kaonMass;// = 0.493677; + + }; #endif diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.cc b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.cc index 7ae0ff29273c7f5ed7c8dbc5704a02988c0f3e9f..c6d767e343d73715c191baefc6d2863b5a445275 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.cc +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.cc @@ -9,7 +9,9 @@ JpsiGamKsKlKKReader::JpsiGamKsKlKKReader() JpsiGamKsKlKKReader::JpsiGamKsKlKKReader(const std::vector<std::string>& files, int particles, int skip): numParticles(particles), - linesToSkip(skip) + linesToSkip(skip), + _useMassRange(false), + _massRange(std::make_pair(0.,10.) ) { if (0 == files.size()) { Alert << "empty list of event files" ; // << endmsg; @@ -38,10 +40,19 @@ bool JpsiGamKsKlKKReader::fillAll(EventList& evtList) double e,px,py,pz; Event* newEvent = new Event(); int parts; + + Vector4<double> fv2Phi(0,0,0,0); //phi phi four-vector for (parts = 0; parts < numParticles; parts++) { currentStream >> px >> py >> pz >> e; newEvent->addParticle(e,px,py,pz); + Vector4<double> tmp = newEvent->p4(parts); + if(parts>0) fv2Phi= fv2Phi+tmp; } + + if(_useMassRange){ + if(fv2Phi.Mass()<_massRange.first || fv2Phi.Mass()>_massRange.second ) continue; + } + if (!currentStream.fail()) { evtList.add(newEvent); for (parts = 0; parts < linesToSkip; parts++) diff --git a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.hh b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.hh index c4a027a3e7d8a7d27cf9a026a5d84491ba421ae8..dec93e4a87c020ebe898e443f46d03afeba84fea 100644 --- a/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.hh +++ b/Examples/JpsiGamKsKlKK/JpsiGamKsKlKKReader.hh @@ -8,6 +8,7 @@ #include <iostream> #include <fstream> #include <cstdlib> +#include <utility> class EventList; @@ -21,13 +22,21 @@ public: virtual ~JpsiGamKsKlKKReader(); virtual bool fillAll(EventList& evtList); - + bool fillMassRange(EventList& evtList, std::pair<double,double> massRange){ + _useMassRange=true; + _massRange = massRange; + return fillAll(evtList); + }; + private: std::vector<std::string> fileNames; std::vector<std::string>::const_iterator currentFile; std::ifstream currentStream; int numParticles; int linesToSkip; + bool _useMassRange; + std::pair<double,double> _massRange; + }; #endif diff --git a/Examples/JpsiGamKsKlKK/MJpsiGamKlKsKKApp.cc b/Examples/JpsiGamKsKlKK/MJpsiGamKlKsKKApp.cc index 57763ded71351c9f855e3e72a3fd25a28ca8c191..548c6101e22c095f5e185277c4e27c4459894392 100644 --- a/Examples/JpsiGamKsKlKK/MJpsiGamKlKsKKApp.cc +++ b/Examples/JpsiGamKsKlKK/MJpsiGamKlKsKKApp.cc @@ -82,7 +82,7 @@ int main(int __argc,char *__argv[]){ std::string theCfgFile = theAppParams.getConfigFile(); - + const std::string datFile=theAppParams.dataFile(); const std::string mcFile=theAppParams.mcFile(); Info << "data file: " << datFile ; // << endmsg; @@ -96,13 +96,15 @@ int main(int __argc,char *__argv[]){ Alert << "Error: could not parse " << pdtFile ; // << endmsg; exit(1); } - + + std::pair<double, double> massRange = theAppParams.massRange(); + Info << "Mass range: " << massRange.first << " " << massRange.second ; + std::vector<std::string> fileNames; - fileNames.push_back(datFile); JpsiGamKsKlKKReader eventReader(fileNames, 5, 0); EventList eventsData; - eventReader.fillAll(eventsData); + eventReader.fillMassRange(eventsData, massRange ); if (!eventsData.findParticleTypes(pTable)) Warning << "could not find all particles" ; // << endmsg; @@ -131,7 +133,7 @@ int main(int __argc,char *__argv[]){ fileNamesMc.push_back(mcFile); JpsiGamKsKlKKReader eventReaderMc(fileNamesMc, 5, 0); EventList eventsMc; - eventReaderMc.fillAll(eventsMc); + eventReaderMc.fillMassRange(eventsMc, massRange); eventsMc.rewind(); // @@ -141,7 +143,7 @@ int main(int __argc,char *__argv[]){ // - //disable hypotheses, currently not in use + //disable hypotheses // std::map<const std::string, bool> hypMap; hypMap["etacHyp"] =false; @@ -172,16 +174,18 @@ int main(int __argc,char *__argv[]){ std::string startWithHyp=theAppParams.startHypo(); if (startWithHyp=="production"){ - theLhPtr = boost::shared_ptr<AbsLh> (new JpsiGamKsKlKKProdLh(theJpsiGamKsKlKKEventListPtr, hypMap )); + JpsiGamKsKlKKProdLh* tmplh = new JpsiGamKsKlKKProdLh(theJpsiGamKsKlKKEventListPtr, hypMap); + tmplh->massIndependentFit( theAppParams.massIndependentFit() ); + tmplh->useCommonProductionPhase( theAppParams.useCommonProductionPhases() ); + theLhPtr = boost::shared_ptr<AbsLh> (tmplh); + //theLhPtr = boost::shared_ptr<AbsLh> (new JpsiGamKsKlKKProdLh(theJpsiGamKsKlKKEventListPtr, hypMap )); } else { Alert << "start with hypothesis " << startWithHyp << " not supported!!!!" ; // << endmsg; exit(1); } - - - + std::string mode=theAppParams.mode(); @@ -259,8 +263,14 @@ int main(int __argc,char *__argv[]){ fitParams finalFitParams=theFitParamBase->getFitParamVal(finalParamVec); + + //calculate intensity contributions + //EvtData* evtdata = eventsData.getDataVecs(); + //theLhPtr->calcEvtIntensity( evtdata , finalFitParams ); + + JpsiGamKsKlKKHist theHist(theLhPtr, finalFitParams); theFitParamBase->printParams(finalFitParams); diff --git a/SetEnv_rub b/SetEnv_rub index 5fd464050a7dfa717fd78e140690d2faab639f69..2d39cf792aead8b70f664f732a313a1b40b432ee 100755 --- a/SetEnv_rub +++ b/SetEnv_rub @@ -1,7 +1,7 @@ setenv ROOTSYS /opt/root/5.26-00b/Linux26SL5_x86_64_gcc412 #setenv CMAKE_SOURCE_DIR /data/pegasus2/bertram/PawianPsi2S_V2/Pawian -setenv CMAKE_SOURCE_DIR `pwd` +setenv CMAKE_SOURCE_DIR `pwd | sed -e 's/\/nfs//'` setenv LD_LIBRARY_PATH ${ROOTSYS}/lib setenv PATH ${PATH}:${ROOTSYS}/bin:/data/pegasus1/PANDA/bin/ setenv BOOST_BUILD_PATH /data/pegasus1/PANDA/PWA/boost-build/