Skip to content
Snippets Groups Projects
Commit 6c157fef authored by Marc Pelizaeus's avatar Marc Pelizaeus
Browse files

Add f2 amplitude

parent c6567a9f
No related branches found
No related tags found
No related merge requests found
......@@ -11,25 +11,14 @@ JpsiGamKsKlKKProdLh::JpsiGamKsKlKKProdLh(boost::shared_ptr<const EvtDataBaseList
AbsLh(theEvtList)
,_eta2225Hyp(true)
,_f02020Hyp(true)
,_f22300Hyp(true)
{
initializeHypothesisMap( hypMap);
std::map<const std::string, bool>::const_iterator iter= hypMap.find("eta2225Hyp");
if (iter !=hypMap.end()){
_eta2225Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _eta2225Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis eta2225Hyp not set!!!" <<endmsg;
iter= hypMap.find("f02020Hyp");
if (iter !=hypMap.end()){
_f02020Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _f02020Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis f02020Hyp not set!!!" <<endmsg;
}
......@@ -37,24 +26,28 @@ JpsiGamKsKlKKProdLh::JpsiGamKsKlKKProdLh( boost::shared_ptr<AbsLh> theLhPtr, con
AbsLh(theLhPtr->getEventList())
,_eta2225Hyp(true)
,_f02020Hyp(true)
,_f22300Hyp(true)
{
std::map<const std::string, bool>::const_iterator iter= hypMap.find("eta2225Hyp");
initializeHypothesisMap( hypMap);
if (iter !=hypMap.end()){
_eta2225Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _eta2225Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis eta2225Hyp not set!!!" <<endmsg;
// std::map<const std::string, bool>::const_iterator iter= hypMap.find("eta2225Hyp");
iter= hypMap.find("f02020Hyp");
if (iter !=hypMap.end()){
_f02020Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _f02020Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis f02020Hyp not set!!!" <<endmsg;
// if (iter !=hypMap.end()){
// _eta2225Hyp= iter->second;
// Info<< "hypothesis " << iter->first << "\t" << _eta2225Hyp <<endmsg;
// _hypMap[iter->first]= iter->second;
// }
// else Alert << "hypothesis eta2225Hyp not set!!!" <<endmsg;
// iter= hypMap.find("f02020Hyp");
// if (iter !=hypMap.end()){
// _f02020Hyp= iter->second;
// Info<< "hypothesis " << iter->first << "\t" << _f02020Hyp <<endmsg;
// _hypMap[iter->first]= iter->second;
// }
// else Alert << "hypothesis f02020Hyp not set!!!" <<endmsg;
}
......@@ -109,7 +102,24 @@ double JpsiGamKsKlKKProdLh::calcEvtIntensity(EvtData* theData, fitParams& thePar
JmmGmp+=f0GammaAmp(-1, 1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,mass,width );
JmmGmm+=f0GammaAmp(-1, -1, theData, PsiTof02020GamMag, PsiTof02020GamPhi,F02020ToPhiPhiMag,F02020ToPhiPhiPhi,mass,width );
}
if (_f22300Hyp){
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiTof2GamMag=theParamVal.Mags[paramEnumJpsiGamKsKlKK::PsiToF22300Gamma];
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiTof2GamPhi=theParamVal.Phis[paramEnumJpsiGamKsKlKK::PsiToF22300Gamma];
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];
JmpGmp+=f2GammaAmp(1, 1, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width );
JmpGmp+=f2GammaAmp(1, -1, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width );
JmpGmp+=f2GammaAmp(-1, 1, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width );
JmpGmp+=f2GammaAmp(-1, -1, theData, PsiTof2GamMag, PsiTof2GamPhi,F2ToPhiPhiMag,F2ToPhiPhiPhi,mass,width );
}
result=norm(JmpGmp)+norm(JmpGmm)+norm(JmmGmp)+norm(JmmGmm);
return result;
}
......@@ -195,10 +205,6 @@ return result;
complex<double> JpsiGamKsKlKKProdLh::etaToPhiPhiTo4KAmp(EvtData* theData){
complex<double> result(0.,0.);
......@@ -238,6 +244,93 @@ complex<double> JpsiGamKsKlKKProdLh::f0ToPhiPhiTo4KAmp(EvtData* theData,
complex<double> JpsiGamKsKlKKProdLh::f2GammaAmp(Spin Minit, Spin Mgamma, 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 ){
complex<double> result(0.,0.);
Vector4<double> fv2Phi= theData->FourVecs[enumJpsiGamKsKlKKData::V4_KsKlKpKm_HeliPsi];
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itPsi;
for ( itPsi=ampf2ProdMag.begin(); itPsi!=ampf2ProdMag.end(); ++itPsi){
boost::shared_ptr<const JPCLS> PsiState=itPsi->first;
double thePsiMag=itPsi->second;
double thePsiPhi=ampf2ProdMag[PsiState];
complex<double> expiphiPsi(cos(thePsiPhi), sin(thePsiPhi));
for(Spin f2Lambda=-2;f2Lambda<=2; f2Lambda++){
Spin lambda = f2Lambda-Mgamma;
if( fabs(lambda)>PsiState->J || fabs(lambda)>PsiState->S) continue;
complex<double> amp = thePsiMag*expiphiPsi*sqrt(2*PsiState->L+1)
*Clebsch(PsiState->L, 0, PsiState->S, lambda, PsiState->J, lambda)
*Clebsch(2, f2Lambda, 1, -Mgamma, PsiState->S, lambda )
*conj( theData->WignerDs[enumJpsiGamKsKlKKData::Df_Psi][PsiState->J][Minit][lambda] );
amp=amp*f2ToPhiPhiTo4KAmp(theData, f2Lambda, ampf2DecMag,ampf2DecPhi );
result+= amp;
}
}
result*=BreitWigner( fv2Phi, mass, width);
return result;
}
complex<double> JpsiGamKsKlKKProdLh::f2ToPhiPhiTo4KAmp( EvtData* theData, Spin f2Lambda,
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &ampf2DecMag ,
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > &ampf2DecPhi ){
complex<double> result(0.,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;
double theMag=itf2->second;
double thePhi=ampf2DecMag[f2State];
complex<double> expiphi(cos(thePhi), sin(thePhi));
for(Spin lambdaPhi1=-1; lambdaPhi1<=1; lambdaPhi1++){
for(Spin lambdaPhi2=-1; lambdaPhi2<=1; lambdaPhi2++){
Spin lambda = lambdaPhi1-lambdaPhi2;
if( fabs(lambda)>f2State->J || fabs(lambda)>f2State->S) continue;
complex<double> amp = theMag*expiphi*sqrt(2*f2State->L+1)
*Clebsch(f2State->L, 0, f2State->S, lambda, f2State->J, lambda)
*Clebsch(1, lambdaPhi1, 1, -lambdaPhi2, f2State->S, lambda )
*conj( theData->WignerDs[enumJpsiGamKsKlKKData::Df_Spin2][f2State->J][f2Lambda][lambda] );
amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 );
result +=amp;
}
}
}
return result;
}
complex<double> JpsiGamKsKlKKProdLh::phiphiTo4KAmp( EvtData* theData, Spin lambdaPhi1, Spin lambdaPhi2 ){
complex<double> result(0.,0.);
result = 3. * conj(theData->WignerDs[enumJpsiGamKsKlKKData::Df_KsKl][1][lambdaPhi1][0])
* conj(theData->WignerDs[enumJpsiGamKsKlKKData::Df_KpKm][1][lambdaPhi2][0]);
return result;
}
void JpsiGamKsKlKKProdLh::print(std::ostream& os) const{
os << "JpsiGamKsKlKKProdLh::print\n";
}
......@@ -257,6 +350,13 @@ void JpsiGamKsKlKKProdLh::getDefaultParams(fitParams& fitVal, fitParams& fitErr)
theAmpMap[paramEnumJpsiGamKsKlKK::PsiToF02020Gamma] = theFitParams.jpclsVec(paramEnumJpsiGamKsKlKK::PsiToF02020Gamma);
theAmpMap[paramEnumJpsiGamKsKlKK::F02020ToPhiPhi] = theFitParams.jpclsVec(paramEnumJpsiGamKsKlKK::F02020ToPhiPhi);
}
if(_f22300Hyp){
theAmpMap[paramEnumJpsiGamKsKlKK::PsiToF22300Gamma] = theFitParams.jpclsVec(paramEnumJpsiGamKsKlKK::PsiToF22300Gamma);
theAmpMap[paramEnumJpsiGamKsKlKK::F22300ToPhiPhi] = theFitParams.jpclsVec(paramEnumJpsiGamKsKlKK::F22300ToPhiPhi);
}
std::map<int, std::vector< boost::shared_ptr<const JPCLS> > >::iterator itAmpMap;
......@@ -302,6 +402,46 @@ void JpsiGamKsKlKKProdLh::getDefaultParams(fitParams& fitVal, fitParams& fitErr)
fitErr.Widths[paramEnumJpsiGamKsKlKK::f02020]=0.1;
}
if(_f22300Hyp){
fitVal.Masses[paramEnumJpsiGamKsKlKK::f22300]=2.02;
fitErr.Masses[paramEnumJpsiGamKsKlKK::f22300]=0.3;
fitVal.Widths[paramEnumJpsiGamKsKlKK::f22300]=.30;
fitErr.Widths[paramEnumJpsiGamKsKlKK::f22300]=0.1;
}
}
bool
JpsiGamKsKlKKProdLh::initializeHypothesisMap( const std::map<const std::string, bool>& hypMap ){
std::map<const std::string, bool>::const_iterator iter= hypMap.find("eta2225Hyp");
if (iter !=hypMap.end()){
_eta2225Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _eta2225Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis eta2225Hyp not set!!!" <<endmsg;
iter= hypMap.find("f02020Hyp");
if (iter !=hypMap.end()){
_f02020Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _f02020Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis f02020Hyp not set!!!" <<endmsg;
iter= hypMap.find("f22300Hyp");
if (iter !=hypMap.end()){
_f22300Hyp= iter->second;
Info<< "hypothesis " << iter->first << "\t" << _f22300Hyp <<endmsg;
_hypMap[iter->first]= iter->second;
}
else Alert << "hypothesis f22300Hyp not set!!!" <<endmsg;
return true;
}
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