Skip to content
Snippets Groups Projects
Commit 8a078c86 authored by Bertram Kopf's avatar Bertram Kopf
Browse files

added missing classes in JpsiToPhiPhiGam

parent 2b3ee9f9
No related branches found
No related tags found
No related merge requests found
#include <getopt.h>
#include <fstream>
#include <string>
#include "Examples/JpsiToPhiPhiGam/JpsiToPhiPhiGamLh.hh"
#include "Examples/JpsiToPhiPhiGam/XToPhiPhiDecAmps.hh"
#include "PwaUtils/FitParamsBaseNew.hh"
#include "PwaUtils/AbsXdecAmp.hh"
#include "ErrLogger/ErrLogger.hh"
#include <boost/bind.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>
JpsiToPhiPhiGamLh::JpsiToPhiPhiGamLh(boost::shared_ptr<const EvtDataBaseListNew> theEvtList, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates) :
PsiProdBaseLhNew(theEvtList, hypVec, theStates)
,_jpsiToPhiPhiGamStatesPtr(theStates)
{
initializeHypothesis();
}
JpsiToPhiPhiGamLh::JpsiToPhiPhiGamLh( boost::shared_ptr<AbsLhNew> theLhPtr, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates ) :
PsiProdBaseLhNew(theLhPtr->getEventList(), hypVec, theStates)
,_jpsiToPhiPhiGamStatesPtr(theStates)
{
initializeHypothesis();
}
JpsiToPhiPhiGamLh::~JpsiToPhiPhiGamLh()
{;
}
void JpsiToPhiPhiGamLh::print(std::ostream& os) const{
}
void JpsiToPhiPhiGamLh::initializeHypothesis(){
std::vector<std::string>::const_iterator it;
for (it=_GammaEtaHyps.begin(); it!=_GammaEtaHyps.end(); ++it){
size_t pos=it->find(_EtaKey);
std::string etaDecAmpName=it->substr(pos);
_etaDecAmpMap[*it]=boost::shared_ptr<AbsXdecAmp>(new XToPhiPhiDecAmps( etaDecAmpName, _hypVec, _jpsiToPhiPhiGamStatesPtr, Spin(0)) );
_allAmpMap[*it]=_etaDecAmpMap[*it];
}
for (it=_GammaEta2Hyps.begin(); it!=_GammaEta2Hyps.end(); ++it){
size_t pos=it->find(_Eta2Key);
std::string eta2DecAmpName=it->substr(pos);
_eta2DecAmpMap[*it]=boost::shared_ptr<AbsXdecAmp>(new XToPhiPhiDecAmps( eta2DecAmpName, _hypVec, _jpsiToPhiPhiGamStatesPtr, Spin(2)) );
_allAmpMap[*it]=_eta2DecAmpMap[*it];
}
for (it=_GammaF0Hyps.begin(); it!=_GammaF0Hyps.end(); ++it){
size_t pos=it->find(_F0Key);
std::string f0DecAmpName=it->substr(pos);
_f0DecAmpMap[*it]=boost::shared_ptr<AbsXdecAmp>(new XToPhiPhiDecAmps( f0DecAmpName, _hypVec, _jpsiToPhiPhiGamStatesPtr, Spin(0)) );
_allAmpMap[*it]=_f0DecAmpMap[*it];
}
for (it=_GammaF1Hyps.begin(); it!=_GammaF1Hyps.end(); ++it){
size_t pos=it->find(_F1Key);
std::string f1DecAmpName=it->substr(pos);
_f1DecAmpMap[*it]=boost::shared_ptr<AbsXdecAmp>(new XToPhiPhiDecAmps( f1DecAmpName, _hypVec, _jpsiToPhiPhiGamStatesPtr, Spin(1)) );
_allAmpMap[*it]=_f1DecAmpMap[*it];
}
for (it=_GammaF2Hyps.begin(); it!=_GammaF2Hyps.end(); ++it){
size_t pos=it->find(_F2Key);
std::string f2DecAmpName=it->substr(pos);
_f2DecAmpMap[*it]=boost::shared_ptr<AbsXdecAmp>(new XToPhiPhiDecAmps( f2DecAmpName, _hypVec, _jpsiToPhiPhiGamStatesPtr, Spin(2)) );
_allAmpMap[*it]=_f2DecAmpMap[*it];
}
}
#pragma once
#include <iostream>
#include <fstream>
#include <string>
#include <vector>
#include <complex>
#include <cassert>
#include <boost/shared_ptr.hpp>
#include <boost/function.hpp>
#include "TROOT.h"
// #include <TSystem.h>
#include "qft++/topincludes/relativistic-quantum-mechanics.hh"
#include "PwaUtils/PsiProdBaseLhNew.hh"
#include "Examples/JpsiToPhiPhiGam/JpsiToPhiPhiGamStates.hh"
#include "PwaUtils/DataUtils.hh"
#include "Minuit2/MnUserParameters.h"
class AbsXdecAmp;
class JpsiToPhiPhiGamLh : public PsiProdBaseLhNew{
public:
// create/copy/destroy:
///Constructor
JpsiToPhiPhiGamLh(boost::shared_ptr<const EvtDataBaseListNew>, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates);
JpsiToPhiPhiGamLh(boost::shared_ptr<AbsLhNew>, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates);
/** Destructor */
virtual ~JpsiToPhiPhiGamLh();
virtual AbsLhNew* clone_() const {
return new JpsiToPhiPhiGamLh(_evtListPtr, _hypVec, _jpsiToPhiPhiGamStatesPtr);
}
// virtual double calcEvtIntensity( EvtDataNew* theData, fitParamsNew& theParamVal);
//Getters:
// virtual void getDefaultParams(fitParamsNew& fitVal, fitParamsNew& fitErr);
virtual void print(std::ostream& os) const;
protected:
boost::shared_ptr<JpsiToPhiPhiGamStates> _jpsiToPhiPhiGamStatesPtr;
private:
void initializeHypothesis();
};
#include <getopt.h>
#include <fstream>
#include <string>
#include "Examples/JpsiToPhiPhiGam/XToPhiPhiDecAmps.hh"
#include "qft++/relativistic-quantum-mechanics/Utils.hh"
#include "ErrLogger/ErrLogger.hh"
// #include "PwaUtils/EvtDataBaseListNew.hh"
#include "Examples/JpsiToPhiPhiGam/JpsiToPhiPhiGamEventList.hh"
XToPhiPhiDecAmps::XToPhiPhiDecAmps(const std::string& name, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates, Spin spinX) :
AbsXdecAmp(name, hypVec, spinX)
,_phiPhiKey(name+"ToPhiPhi")
,_xBWKey(name+"BreitWigner")
,_massIndependent(true)
,_theStatesPtr(theStates)
{
initialize();
}
XToPhiPhiDecAmps::~XToPhiPhiDecAmps()
{
}
complex<double> XToPhiPhiDecAmps::XdecAmp(Spin lamX, EvtDataNew* theData, fitParamsNew& theParamVal){
complex<double> result(0.,0.);
result+=XToPhiPhiAmp(lamX, theData, theParamVal);
complex<double> dynModel(1.,0.);
if (!_massIndependent){
double xMass=theParamVal.Masses[_name];
double xWidth=theParamVal.Widths[_name];
Vector4<double> p4PhiPhi = theData->FourVecsDec[enumJpsiGamX4V::V4_KsKlKpKm_HeliPsi];
dynModel=BreitWigner(p4PhiPhi, xMass, xWidth);
}
result *=dynModel;
return result;
}
complex<double> XToPhiPhiDecAmps::XToPhiPhiAmp(Spin lamX, EvtDataNew* theData, fitParamsNew& theParamVal){
complex<double> result(0.,0.);
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > XToPhiPhiMag=theParamVal.Mags[_phiPhiKey];
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > XToPhiPhiPhi=theParamVal.Phis[_phiPhiKey];
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itXMag;
for ( itXMag=XToPhiPhiMag.begin(); itXMag!=XToPhiPhiMag.end(); ++itXMag){
boost::shared_ptr<const JPCLS> XState=itXMag->first;
double theXMag=itXMag->second;
double theXPhi=XToPhiPhiPhi[XState];
complex<double> expiphiX(cos(theXPhi), sin(theXPhi));
for(Spin lambdaPhi1=-1; lambdaPhi1<=1; lambdaPhi1++){
for(Spin lambdaPhi2=-1; lambdaPhi2<=1; lambdaPhi2++){
Spin lambda = lambdaPhi1-lambdaPhi2;
if( fabs(lambda)>XState->J || fabs(lambda)>XState->S) continue;
complex<double> amp = theXMag*expiphiX*sqrt(2*XState->L+1)
*Clebsch(XState->L, 0, XState->S, lambda, XState->J, lambda)
*Clebsch(1, lambdaPhi1, 1, -lambdaPhi2, XState->S, lambda )
*conj( theData->WignerDsDec[enumJpsiGamXDfunc::Df_XToPhiPhi][XState->J][_J_X][lambda] );
amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 );
result +=amp;
}
}
}
result*=conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_XToPhiPhi][_J_X][lamX][0]);
return result;
}
complex<double> XToPhiPhiDecAmps::phiphiTo4KAmp( EvtDataNew* theData, Spin lambdaPhi1, Spin lambdaPhi2 ){
complex<double> result(0.,0.);
result = 3. * conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_PhiToKsKl][1][lambdaPhi1][0])
* 3.* conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_PhiToKpKm][1][lambdaPhi2][0]);
return result;
}
void XToPhiPhiDecAmps::getDefaultParams(fitParamsNew& fitVal, fitParamsNew& fitErr){
std::vector< boost::shared_ptr<const JPCLS> > PhiPhiStates;
if(_J_X==0) PhiPhiStates=_theStatesPtr->EtaToPhiPhiStates();
else if(_J_X==1) PhiPhiStates=_theStatesPtr->F1ToPhiPhiStates();
else if(_J_X==2) PhiPhiStates=_theStatesPtr->Eta2ToPhiPhiStates();
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > currentMagValMap;
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > currentPhiValMap;
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > currentMagErrMap;
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > currentPhiErrMap;
std::vector< boost::shared_ptr<const JPCLS> >::const_iterator itLS;
for(itLS=PhiPhiStates.begin(); itLS!=PhiPhiStates.end(); ++itLS){
currentMagValMap[*itLS]=0.2;
currentPhiValMap[*itLS]=0.;
currentMagErrMap[*itLS]=0.8;
currentPhiErrMap[*itLS]=0.3;
}
fitVal.Mags[_phiPhiKey]=currentMagValMap;
fitVal.Phis[_phiPhiKey]=currentPhiValMap;
fitErr.Mags[_phiPhiKey]=currentMagErrMap;
fitErr.Phis[_phiPhiKey]=currentPhiErrMap;
if (!_massIndependent){
size_t pos=_name.find("_");
std::string massMeVString=_name.substr(pos+1);
stringstream massMeVStrStream(massMeVString);
int MassMeV;
massMeVStrStream >> MassMeV;
double MassGeV= ( (double) MassMeV)/1000.;
fitVal.Masses[_name]=MassGeV;
fitErr.Masses[_name]=0.01;
fitVal.Widths[_name]=0.2;
fitErr.Widths[_name]=0.02;
}
}
void XToPhiPhiDecAmps::print(std::ostream& os) const{
return; //dummy
}
void XToPhiPhiDecAmps::initialize(){
std::vector<std::string>::const_iterator it;
for (it=_hypVec.begin(); it!=_hypVec.end(); ++it){
if (it->compare(0, _xBWKey.size(), _xBWKey) ==0){
_massIndependent=false;
}
}
}
#pragma once
#include <iostream>
#include <vector>
#include <complex>
#include <map>
#include <string>
#include <cassert>
#include <boost/shared_ptr.hpp>
#include "PwaUtils/AbsXdecAmp.hh"
#include "Examples/JpsiToPhiPhiGam/JpsiToPhiPhiGamStates.hh"
class XToPhiPhiDecAmps : public AbsXdecAmp{
public:
// create/copy/destroy:
///Constructor
XToPhiPhiDecAmps(const std::string& name, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToPhiPhiGamStates> theStates, Spin spinX);
/** Destructor */
virtual ~XToPhiPhiDecAmps();
// Getters:
virtual complex<double> XdecAmp(Spin lamX, EvtDataNew* theData, fitParamsNew& theParamVal);
virtual void getDefaultParams(fitParamsNew& fitVal, fitParamsNew& fitErr);
virtual void print(std::ostream& os) const;
protected:
const std::string _phiPhiKey;
const std::string _xBWKey;
bool _massIndependent;
boost::shared_ptr<JpsiToPhiPhiGamStates> _theStatesPtr;
complex<double> XToPhiPhiAmp(Spin lamX, EvtDataNew* theData, fitParamsNew& theParamVal);
complex<double> phiphiTo4KAmp( EvtDataNew* theData, Spin lambdaPhi1, Spin lambdaPhi2 );
virtual void initialize();
private:
};
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment