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

improved again performance for calculation LS amplitudes

parent 7ba674c2
No related branches found
No related tags found
No related merge requests found
......@@ -154,6 +154,7 @@ void AbsLSDecAmps::updateFitParams(fitParams& theParamVal){
_currentParamPhis[*it]=thePhi;
complex<double> expi(cos(thePhi), sin(thePhi));
_currentParamMagExpi[*it]=theMag*expi;
_currentParamPreFacMagExpi[*it]=_preFactor*_isospinCG*theMag*expi;
}
_absDyn->updateFitParams(theParamVal);
......
......@@ -77,6 +77,7 @@ protected:
std::map< std::shared_ptr<const LScomb>, double, pawian::Collection::SharedPtrLess > _currentParamMags;
std::map< std::shared_ptr<const LScomb>, double, pawian::Collection::SharedPtrLess > _currentParamPhis;
std::map< std::shared_ptr<const LScomb>, complex<double>, pawian::Collection::SharedPtrLess > _currentParamMagExpi;
std::map< std::shared_ptr<const LScomb>, complex<double>, pawian::Collection::SharedPtrLess > _currentParamPreFacMagExpi;
Spin _lam1Min;
Spin _lam1Max;
......@@ -85,6 +86,7 @@ protected:
void fillCgPreFactor();
// virtual complex<double> lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtData* theData, Spin& lam1Min, Spin& lam1Max, Spin& lam2Min, Spin& lam2Max, bool withDecs, Spin lamFs=0 )=0;
private:
......
......@@ -111,3 +111,20 @@ void AbsXdecAmp::cacheAmplitudes(){
if(!_daughter2IsStable) _decAmpDaughter2->cacheAmplitudes();
}
void AbsXdecAmp::calcDynamics(EvtData* theData, AbsXdecAmp* grandmaAmp){
if(!_recalculate) return;
// Info << "threadID: " << std::this_thread::get_id() << endmsg;
// if(!_recalculate) return;
if(!_absDyn->isLdependent()){
theMutex.lock();
_cachedDynMap[std::this_thread::get_id()][_absDyn->grandMaKey(grandmaAmp)] = _absDyn->eval( theData, grandmaAmp);
theMutex.unlock();
}
if(!_daughter1IsStable) _decAmpDaughter1->calcDynamics(theData, this);
if(!_daughter2IsStable) _decAmpDaughter2->calcDynamics(theData, this);
return;
}
......@@ -34,6 +34,7 @@
#include <mutex>
#include <memory>
#include <boost/unordered_map.hpp>
#include <thread>
#include "PwaUtils/AbsChannelEnv.hh"
#include "PwaUtils/EvtDataBaseList.hh"
......@@ -68,6 +69,7 @@ public:
const std::string wignerDKey() const {return _wignerDKey;}
const std::string refKey() const {return _refKey;}
virtual void cacheAmplitudes();
virtual void calcDynamics(EvtData* theData, AbsXdecAmp* grandmaAmp=0);
protected:
......@@ -98,6 +100,7 @@ protected:
Spin _J;
intStringShortComplFloatMap _cachedAmpMap;
intStringStringShortComplFloatMap _cachedAmpMapNew;
std::map<std::thread::id, std::map<std::string, complex<float> > > _cachedDynMap;
virtual void initialize();
};
......@@ -123,8 +123,8 @@ complex<double> LSDecNonRefAmps::XdecAmp(Spin& lamX, EvtData* theData, Spin& lam
result=lsLoop(grandmaAmp, lamX, theData, _lam1Min, _lam1Max, _lam2Min, _lam2Max, true, lamFs);
result*=_preFactor*_isospinCG;
if (!_absDyn->isLdependent()) result *=_absDyn->eval(theData, grandmaAmp);
// result*=_preFactor*_isospinCG;
// if (!_absDyn->isLdependent()) result *=_absDyn->eval(theData, grandmaAmp);
if ( _cacheAmps){
theMutex.lock();
......@@ -148,12 +148,12 @@ complex<double> LSDecNonRefAmps::lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtD
complex<double> result(0.,0.);
std::vector< std::shared_ptr<const LScomb> >::iterator it;
std::map< Spin, complex<double> > dynLSs;
if (_absDyn->isLdependent()){
for (it=_LSs.begin(); it!=_LSs.end(); ++it){
dynLSs[(*it)->L]=_absDyn->eval(theData, grandmaAmp, (*it)->L);
}
}
// std::map< Spin, complex<double> > dynLSs;
// if (_absDyn->isLdependent()){
// for (it=_LSs.begin(); it!=_LSs.end(); ++it){
// dynLSs[(*it)->L]=_absDyn->eval(theData, grandmaAmp, (*it)->L);
// }
// }
std::map<Id3StringType, complex<double> >& currentWignerDMap=theData->WignerDStringId.at(_wignerDKey);
......@@ -174,7 +174,8 @@ complex<double> LSDecNonRefAmps::lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtD
// else amp+=theMag*expi*_cgPreFactor.at(*it).at(lambda1).at(lambda2);
// if (_absDyn->isLdependent()) amp+=theMag*expi*cgPre_LSMap.at(*it)*dynLSs.at((*it)->L);
// else amp+=theMag*expi*cgPre_LSMap.at(*it);
if (_absDyn->isLdependent()) amp+=_currentParamMagExpi.at(*it)*cgPre_LSMap.at(*it)*dynLSs.at((*it)->L);
// if (_absDyn->isLdependent()) amp+=_currentParamPreFacMagExpi.at(*it)*cgPre_LSMap.at(*it)*dynLSs.at((*it)->L);
if (_absDyn->isLdependent()) amp+=_currentParamPreFacMagExpi.at(*it)*cgPre_LSMap.at(*it)*_cachedDynLSMap[std::this_thread::get_id()][(*it)->L];
else amp+=_currentParamMagExpi.at(*it)*cgPre_LSMap.at(*it);
}
Id3StringType IdJLamXLam12=FunctionUtils::spin3Index(_J, lamX, lambda);
......@@ -183,8 +184,39 @@ complex<double> LSDecNonRefAmps::lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtD
result+=amp;
}
}
// if (!_absDyn->isLdependent()) result *=_absDyn->eval(theData, grandmaAmp);
if (!_absDyn->isLdependent()) result *=_cachedDynMap.at(std::this_thread::get_id()).at(_absDyn->grandMaKey(grandmaAmp));
return result;
}
void LSDecNonRefAmps::calcDynamics(EvtData* theData, AbsXdecAmp* grandmaAmp){
if(!_recalculate) return;
if(!_absDyn->isLdependent()){
AbsXdecAmp::calcDynamics(theData, grandmaAmp);
return;
}
std::vector< std::shared_ptr<const LScomb> >::iterator it;
for (it=_LSs.begin(); it!=_LSs.end(); ++it){
theMutex.lock();
_cachedDynLSMap[std::this_thread::get_id()][(*it)->L]=_absDyn->eval(theData, grandmaAmp, (*it)->L);
theMutex.unlock();
}
if(!_daughter1IsStable) _decAmpDaughter1->calcDynamics(theData, this);
if(!_daughter2IsStable) _decAmpDaughter2->calcDynamics(theData, this);
return;
}
// void LSDecNonRefAmps::resetRecalcDynamics(EvtData* theData, AbsXdecAmp* grandmaAmp){
// if(!_absDyn->isLdependent()){
// AbsXdecAmp::resetRecalcDynamics(theData, grandmaAmp);
// return;
// }
// _cachedDynRecalcLSMap[std::this_thread::get_id()][(*it)->L]=true;
// if(!_daughter1IsStable) _decAmpDaughter1->resetRecalcDynamics(theData, this);
// if(!_daughter2IsStable) _decAmpDaughter2->resetRecalcDynamics(theData, this);
// return;
// }
......@@ -59,11 +59,12 @@ public:
virtual complex<double> XdecPartAmp(Spin& lamX, Spin& lamDec, short fixDaughterNr,
EvtData* theData, Spin& lamFs, AbsXdecAmp* grandmaAmp);
virtual void calcDynamics(EvtData* theData, AbsXdecAmp* grandmaAmp=0);
protected:
virtual complex<double> lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtData* theData, Spin& lam1Min, Spin& lam1Max, Spin& lam2Min, Spin& lam2Max, bool withDecs, Spin lamFs=0 );
Spin _Smax;
std::map<std::thread::id, std::map<Spin, complex<double> > > _cachedDynLSMap;
private:
......
......@@ -152,7 +152,7 @@ complex<double> LSDecRefAmps::lsLoopRef(AbsXdecAmp* grandmaAmp, std::string& ref
// double theMag=_currentParamMags.at(*it);
// double thePhi=_currentParamPhis.at(*it);
// complex<double> expi(cos(thePhi), sin(thePhi));
complex<double> currentMagExpi=_currentParamMagExpi.at(*it);
complex<double> currentPreMagExpi=_currentParamPreFacMagExpi.at(*it);
complex<double> tmpResult(0.,0.);
for(Spin lambda1=lam1Min; lambda1<=lam1Max; ++lambda1){
......@@ -161,7 +161,7 @@ complex<double> LSDecRefAmps::lsLoopRef(AbsXdecAmp* grandmaAmp, std::string& ref
if( fabs(lambda)>_JPCPtr->J || fabs(lambda)>(*it)->S) continue;
Id3StringType IdJLamXLam12=FunctionUtils::spin3Index(_J, lamX, lambda);
// complex<double> amp = theMag*expi*currentCgFactor.at(lambda1).at(lambda2)*conj(currentWignerDMap.at(IdJLamXLam12));
complex<double> amp = currentMagExpi*currentCgFactor.at(lambda1).at(lambda2)*conj(currentWignerDMap.at(IdJLamXLam12));
complex<double> amp = currentPreMagExpi*currentCgFactor.at(lambda1).at(lambda2)*conj(currentWignerDMap.at(IdJLamXLam12));
if(withDecs) amp *=daughterAmp(lambda1, lambda2, theData, lamFs);
tmpResult+=amp;
}
......@@ -172,7 +172,7 @@ complex<double> LSDecRefAmps::lsLoopRef(AbsXdecAmp* grandmaAmp, std::string& ref
result+=tmpResult;
}
result*=_preFactor*_isospinCG;
// result*=_preFactor*_isospinCG;
// if(result.real()!=result.real()){
// Alert << "result:\t" << result << endmsg;
// exit(0);
......
......@@ -92,24 +92,21 @@ complex<double> LSOmegaTo3PiDecAmps::XdecAmp(Spin& lamX, EvtData* theData, Spin&
std::vector< std::shared_ptr<const LScomb> >::iterator it;
for (it=_LSs.begin(); it!=_LSs.end(); ++it){
if( fabs(lamX) > _JPCPtr->J ) continue;
double theMag=_currentParamMags[*it];
double thePhi=_currentParamPhis[*it];
complex<double> expi(cos(thePhi), sin(thePhi));
complex<double> amp = theMag*expi*sqrt(2*(*it)->L+1)
complex<double> amp = _currentParamMagExpi.at(*it)*sqrt(2*(*it)->L+1)
*conj( theData->WignerDsStringString.at(_wignerDKey).at(refKey).at(_JPCPtr->J).at(lamX).at(0));
if (_absDyn->isLdependent()) amp*=_absDyn->eval(theData, grandmaAmp, (*it)->L);
result+=amp;
}
result*=sqrt( theData->DoubleStringString[_lambdaDecKey][refKey] );
if (!_absDyn->isLdependent()) result*=_absDyn->eval(theData, grandmaAmp);
if ( _cacheAmps){
theMutex.lock();
_cachedAmpMapNew[evtNo][refKey][_absDyn->grandMaKey(grandmaAmp)][currentSpinIndex]=result;
theMutex.unlock();
}
result*=_absDyn->eval(theData, grandmaAmp);
// result*=_absDyn->eval(theData, grandmaAmp);
return result;
}
......
......@@ -84,6 +84,11 @@ double epemBaseLh::calcEvtIntensity(EvtData* theData, fitParams& theParamVal){
double result=0.;
std::vector< std::shared_ptr<AbsXdecAmp> >::iterator itDecAll;
for (itDecAll=_decAmps.begin(); itDecAll!=_decAmps.end(); ++itDecAll){
(*itDecAll)->calcDynamics(theData);
}
Spin lamSteps=1;
if(_isHighestJaPhoton) lamSteps=2;
......
......@@ -149,6 +149,11 @@ double pbarpBaseLh::calcEvtIntensity(EvtData* theData, fitParams& theParamVal){
double result=0.;
std::vector< std::shared_ptr<AbsXdecAmp> >::iterator itDecAll;
for (itDecAll=_decAmps.begin(); itDecAll!=_decAmps.end(); ++itDecAll){
(*itDecAll)->calcDynamics(theData);
}
std::map <std::shared_ptr<const JPCLS>, std::vector< std::shared_ptr<AbsXdecAmp> >, pawian::Collection::SharedPtrLess >::iterator it;
Spin lamSteps=1;
......
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