From 35f2aa60276f4297662a8838aa99ad17a953603b Mon Sep 17 00:00:00 2001
From: Bertram Kopf <bertram@pc14.ep1.rub.de>
Date: Thu, 8 Jan 2015 12:14:11 +0100
Subject: [PATCH] improved performance for LS-amp calculations; changed order
 of the lamlam- and LS-loops

---
 Particle/pdtNew.table                    |  1 +
 PwaUtils/AbsDynamics.cc                  |  1 +
 PwaUtils/AbsDynamics.hh                  |  2 +
 PwaUtils/BlattWBarrierDynamics.cc        |  1 +
 PwaUtils/BlattWBarrierTensorDynamics.cc  |  1 +
 PwaUtils/BreitWignerBlattWRelDynamics.cc |  1 +
 PwaUtils/BreitWignerDynamics.cc          |  1 +
 PwaUtils/FlatteDynamics.cc               |  1 +
 PwaUtils/K0star1430LassDynamics.cc       |  1 +
 PwaUtils/KMatrixDynamics.cc              |  1 +
 PwaUtils/KPiSWaveIso12Dynamics.cc        |  1 +
 PwaUtils/KPiSWaveIso32Dynamics.cc        |  1 +
 PwaUtils/LSDecNonRefAmps.cc              | 66 ++++++++++++------------
 PwaUtils/LSDecNonRefAmps.hh              |  2 +
 PwaUtils/PiPiSWaveASDynamics.cc          |  1 +
 PwaUtils/VoigtDynamics.cc                |  1 +
 PwaUtils/WoDynamics.cc                   |  1 +
 17 files changed, 51 insertions(+), 33 deletions(-)

diff --git a/Particle/pdtNew.table b/Particle/pdtNew.table
index a92cef84..d10bfeed 100644
--- a/Particle/pdtNew.table
+++ b/Particle/pdtNew.table
@@ -61,6 +61,7 @@ Meson f2(1900) 		f_2(1900)  		0  4  1  1  1  	0  0  0  0 	1.900 	0.1750 	relBW	[
 Meson etaHat(1600) 	\hat{\eta}_(1600)  	0  0 -1 -1  1  	0  0  0  0 	1.600 	0.3000 	relBW	[0^{--} state]
 Meson f0Hat(1600) 	\hat{f}_0(1600)  	0  0  1 -1  1  	0  0  0  0 	1.600 	0.3000 	relBW	[0^{+-} state]
 Meson f2Hat(1800) 	\hat{f}_2(1800)  	0  4  1 -1  1   0  0  0  0 	1.800 	0.3000 	relBW	[2^{+-} state]
+Meson eta1(1400) 	\eta1(1400)  	        0  2 -1  1  1   0  0  0  0 	1.400 	0.300 	relBW	[exotic state]
 #
 # End light unflavoured mesons with Isospin 0    
 #
diff --git a/PwaUtils/AbsDynamics.cc b/PwaUtils/AbsDynamics.cc
index dd841656..0c4a545b 100644
--- a/PwaUtils/AbsDynamics.cc
+++ b/PwaUtils/AbsDynamics.cc
@@ -41,6 +41,7 @@ AbsDynamics::AbsDynamics(std::string& name, std::vector<Particle*>& fsParticles,
   ,_mother(mother)
   ,_dynKey(FunctionUtils::particleListName(fsParticles))
   ,_grandmaKey("default")
+  ,_isLdependent(true)
 {
   if(0!=mother) _massKey=mother->name(); 
 }
diff --git a/PwaUtils/AbsDynamics.hh b/PwaUtils/AbsDynamics.hh
index e80fd8e3..ee80eb49 100644
--- a/PwaUtils/AbsDynamics.hh
+++ b/PwaUtils/AbsDynamics.hh
@@ -56,6 +56,7 @@ public:
   virtual void setMassKey(std::string& theMassKey){_massKey=theMassKey;}
   virtual void addGrandMa(std::shared_ptr<AbsDecay> theDec) {;}
   virtual const std::string& grandMaKey(AbsXdecAmp* grandmaAmp){return _grandmaKey;}
+  virtual bool isLdependent() {return _isLdependent;}
 
 protected:
   std::string _name;
@@ -66,6 +67,7 @@ protected:
   std::mutex theMutex;
   std::string _dynKey;
   std::string _grandmaKey;
+  bool _isLdependent;
 private:
 
 };
diff --git a/PwaUtils/BlattWBarrierDynamics.cc b/PwaUtils/BlattWBarrierDynamics.cc
index d879d2fe..feec0f7f 100644
--- a/PwaUtils/BlattWBarrierDynamics.cc
+++ b/PwaUtils/BlattWBarrierDynamics.cc
@@ -45,6 +45,7 @@ BlattWBarrierDynamics::BlattWBarrierDynamics(std::string& name, std::vector<Part
 {
   Info << "BlattWBarrierDynamics for " << _name <<endmsg;
   if(GlobalEnv::instance()->parser()->fitqRProduction()) _fitqRVals=true;
+  _isLdependent=true;
 }
 
 BlattWBarrierDynamics::~BlattWBarrierDynamics()
diff --git a/PwaUtils/BlattWBarrierTensorDynamics.cc b/PwaUtils/BlattWBarrierTensorDynamics.cc
index feafe06d..824329d5 100644
--- a/PwaUtils/BlattWBarrierTensorDynamics.cc
+++ b/PwaUtils/BlattWBarrierTensorDynamics.cc
@@ -36,6 +36,7 @@
 BlattWBarrierTensorDynamics::BlattWBarrierTensorDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother, const std::string& wignerDKey, double qR) :
   BlattWBarrierDynamics(name, fsParticles, mother, wignerDKey, qR)
 {
+  _isLdependent=true;
 }
 
 BlattWBarrierTensorDynamics::~BlattWBarrierTensorDynamics()
diff --git a/PwaUtils/BreitWignerBlattWRelDynamics.cc b/PwaUtils/BreitWignerBlattWRelDynamics.cc
index fe7a3fe0..f082ccfc 100644
--- a/PwaUtils/BreitWignerBlattWRelDynamics.cc
+++ b/PwaUtils/BreitWignerBlattWRelDynamics.cc
@@ -38,6 +38,7 @@ BreitWignerBlattWRelDynamics::BreitWignerBlattWRelDynamics(std::string& name, st
   BreitWignerRelDynamics(name, fsParticles, mother, massSumDaughter1, massSumDaughter2)
   ,_qR(qR)
 {
+  _isLdependent=true;
 }
 
 BreitWignerBlattWRelDynamics::~BreitWignerBlattWRelDynamics()
diff --git a/PwaUtils/BreitWignerDynamics.cc b/PwaUtils/BreitWignerDynamics.cc
index a89da423..e1f60df5 100644
--- a/PwaUtils/BreitWignerDynamics.cc
+++ b/PwaUtils/BreitWignerDynamics.cc
@@ -37,6 +37,7 @@
 BreitWignerDynamics::BreitWignerDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother) :
   AbsDynamics(name, fsParticles, mother)
 {
+  _isLdependent=false;
 }
 
 BreitWignerDynamics::~BreitWignerDynamics()
diff --git a/PwaUtils/FlatteDynamics.cc b/PwaUtils/FlatteDynamics.cc
index e2580911..ec0074ec 100644
--- a/PwaUtils/FlatteDynamics.cc
+++ b/PwaUtils/FlatteDynamics.cc
@@ -62,6 +62,7 @@ FlatteDynamics::FlatteDynamics(std::string& key, std::vector<Particle*>& fsParti
   _g22Key=_massKey+secondPair1stPartStr+secondPair2ndPartStr;  
 
   _flattePtr=std::shared_ptr<Flatte>(new Flatte(decPair1stChannel, decPair2ndChannel));
+  _isLdependent=false;
 }
 
 FlatteDynamics::~FlatteDynamics()
diff --git a/PwaUtils/K0star1430LassDynamics.cc b/PwaUtils/K0star1430LassDynamics.cc
index 2cc1744b..aaf0eee9 100644
--- a/PwaUtils/K0star1430LassDynamics.cc
+++ b/PwaUtils/K0star1430LassDynamics.cc
@@ -45,6 +45,7 @@ K0star1430LassDynamics::K0star1430LassDynamics(std::string& name, std::vector<Pa
   ,_RLASSKey(_massKey+"RLASS_PosNeg")
   ,_phiRKey(_massKey+"phiR_PosNeg")
 {
+  _isLdependent=false;
 }
 
 K0star1430LassDynamics::~K0star1430LassDynamics()
diff --git a/PwaUtils/KMatrixDynamics.cc b/PwaUtils/KMatrixDynamics.cc
index 0425b895..baf51002 100644
--- a/PwaUtils/KMatrixDynamics.cc
+++ b/PwaUtils/KMatrixDynamics.cc
@@ -56,6 +56,7 @@ KMatrixDynamics::KMatrixDynamics(std::string& name, std::vector<Particle*>& fsPa
   ,_kMatrixParser(new KMatrixParser(pathToConfigParser))
 {
   init();
+  _isLdependent=false;
 }
 
 KMatrixDynamics::~KMatrixDynamics()
diff --git a/PwaUtils/KPiSWaveIso12Dynamics.cc b/PwaUtils/KPiSWaveIso12Dynamics.cc
index 1942642c..838d5164 100644
--- a/PwaUtils/KPiSWaveIso12Dynamics.cc
+++ b/PwaUtils/KPiSWaveIso12Dynamics.cc
@@ -42,6 +42,7 @@ KPiSWaveIso12Dynamics::KPiSWaveIso12Dynamics(std::string& name, std::vector<Part
   AbsDynamics(name, fsParticles, mother)
 {
   _kMatr =  std::shared_ptr<KMatrixKPiSFocus> (new KMatrixKPiSFocus(1));
+  _isLdependent=false;
 }
 
 KPiSWaveIso12Dynamics::~KPiSWaveIso12Dynamics()
diff --git a/PwaUtils/KPiSWaveIso32Dynamics.cc b/PwaUtils/KPiSWaveIso32Dynamics.cc
index 4e47cec2..0f2511c7 100644
--- a/PwaUtils/KPiSWaveIso32Dynamics.cc
+++ b/PwaUtils/KPiSWaveIso32Dynamics.cc
@@ -42,6 +42,7 @@ KPiSWaveIso32Dynamics::KPiSWaveIso32Dynamics(std::string& name, std::vector<Part
   AbsDynamics(name, fsParticles, mother)
 {
   _kMatr =  std::shared_ptr<KMatrixKPiSFocus> (new KMatrixKPiSFocus(3));
+  _isLdependent=false;
 }
 
 KPiSWaveIso32Dynamics::~KPiSWaveIso32Dynamics()
diff --git a/PwaUtils/LSDecNonRefAmps.cc b/PwaUtils/LSDecNonRefAmps.cc
index b228d43e..33c82859 100644
--- a/PwaUtils/LSDecNonRefAmps.cc
+++ b/PwaUtils/LSDecNonRefAmps.cc
@@ -41,12 +41,22 @@
 
 LSDecNonRefAmps::LSDecNonRefAmps(std::shared_ptr<IsobarLSDecay> theDec, ChannelID channelID) :
   AbsLSDecAmps(theDec, channelID)
+  ,_Smax(0)
 {
+  std::vector< std::shared_ptr<const LScomb> >::iterator it;
+  for (it=_LSs.begin(); it!=_LSs.end(); ++it){
+    if( (*it)->S > _Smax ) _Smax=(*it)->S;
+  }
 }
 
 LSDecNonRefAmps::LSDecNonRefAmps(std::shared_ptr<AbsDecay> theDec, ChannelID channelID) :
   AbsLSDecAmps(theDec, channelID)
+  ,_Smax(0)
 {
+  std::vector< std::shared_ptr<const LScomb> >::iterator it;
+  for (it=_LSs.begin(); it!=_LSs.end(); ++it){
+    if( (*it)->S > _Smax ) _Smax=(*it)->S;
+  }
 }
 
 LSDecNonRefAmps::~LSDecNonRefAmps()
@@ -78,6 +88,8 @@ complex<double> LSDecNonRefAmps::XdecPartAmp(Spin& lamX, Spin& lamDec, short fix
   }
 
   result=lsLoop( grandmaAmp, lamX, theData, _lam1Min, _lam1Max, _lam2Min, _lam2Max, false);
+  result*=_preFactor*_isospinCG;
+  if (!_absDyn->isLdependent()) result *=_absDyn->eval(theData, grandmaAmp);
 
   return result;
 }
@@ -111,6 +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);
 
   if ( _cacheAmps){
      theMutex.lock();
@@ -133,44 +147,30 @@ complex<double> LSDecNonRefAmps::lsLoop(AbsXdecAmp* grandmaAmp, Spin& lamX, EvtD
  
   complex<double> result(0.,0.);
 
-  //  map<Spin,complex<double> >& currentWignerDsMap=theData->WignerDsString.at(_wignerDKey).at(_JPCPtr->J).at(lamX);
-  //  Spin currentJ=_JPCPtr->J;
-  //  std::map<Id3StringType, complex<double> >& currentWignerDMap=theData->WignerDStringStringId.at(_wignerDKey).at(refKey);
-
-  std::map<Id3StringType, complex<double> >& currentWignerDMap=theData->WignerDStringId.at(_wignerDKey);
-
   std::vector< std::shared_ptr<const LScomb> >::iterator it;
-  for (it=_LSs.begin(); it!=_LSs.end(); ++it){
+  std::map<Id3StringType, complex<double> >& currentWignerDMap=theData->WignerDStringId.at(_wignerDKey);
 
-    map<Spin,map<Spin, double > >& currentCgFactor=_cgPreFactor.at(*it);
-
-    double theMag=_currentParamMags.at(*it);
-    double thePhi=_currentParamPhis.at(*it);
-    complex<double> expi(cos(thePhi), sin(thePhi));
-
-    complex<double> tmpResult(0.,0.);
-    for(Spin lambda1=lam1Min; lambda1<=lam1Max; ++lambda1){
-      for(Spin lambda2=lam2Min; lambda2<=lam2Max; ++lambda2){
-	Spin lambda = lambda1-lambda2;
-	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 = theMag*expi*currentCgFactor.at(lambda1).at(lambda2)*conj(currentWignerDsMap.at(lambda));
-      	if(withDecs) amp *=daughterAmp(lambda1, lambda2, theData, lamFs);
-	tmpResult+=amp;
+  for(Spin lambda1=lam1Min; lambda1<=lam1Max; ++lambda1){
+    for(Spin lambda2=lam2Min; lambda2<=lam2Max; ++lambda2){
+      Spin lambda = lambda1-lambda2;
+      if( fabs(lambda)>_JPCPtr->J || fabs(lambda)>_Smax) continue;
+      
+      complex<double> amp(0.,0.);     
+      for (it=_LSs.begin(); it!=_LSs.end(); ++it){
+	if( fabs(lambda)>(*it)->S) continue;
+	double theMag=_currentParamMags.at(*it);
+	double thePhi=_currentParamPhis.at(*it);
+	complex<double> expi(cos(thePhi), sin(thePhi));
+	if (_absDyn->isLdependent()) amp+=theMag*expi*_cgPreFactor.at(*it).at(lambda1).at(lambda2)*_absDyn->eval(theData, grandmaAmp, (*it)->L);
+	else amp+=theMag*expi*_cgPreFactor.at(*it).at(lambda1).at(lambda2);
       }
+      Id3StringType IdJLamXLam12=FunctionUtils::spin3Index(_J, lamX, lambda);
+      amp *= conj(currentWignerDMap.at(IdJLamXLam12));
+      if(withDecs) amp *=daughterAmp(lambda1, lambda2, theData, lamFs);
+      result+=amp;    
     }
-
-    tmpResult*=_absDyn->eval(theData, grandmaAmp, (*it)->L);
-
-    result+=tmpResult; 
   }
-
-  result*=_preFactor*_isospinCG;
-  // if(result.real()!=result.real()){
-  //   Alert << "result:\t" << result << endmsg;
-  //   exit(0);
-  // }
+  
   return result;
 }
 
diff --git a/PwaUtils/LSDecNonRefAmps.hh b/PwaUtils/LSDecNonRefAmps.hh
index a00d314f..c6bbd912 100644
--- a/PwaUtils/LSDecNonRefAmps.hh
+++ b/PwaUtils/LSDecNonRefAmps.hh
@@ -63,6 +63,8 @@ public:
 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;
+
 private:
 
 
diff --git a/PwaUtils/PiPiSWaveASDynamics.cc b/PwaUtils/PiPiSWaveASDynamics.cc
index b7f15e17..3b8c4c94 100644
--- a/PwaUtils/PiPiSWaveASDynamics.cc
+++ b/PwaUtils/PiPiSWaveASDynamics.cc
@@ -53,6 +53,7 @@ PiPiSWaveASDynamics::PiPiSWaveASDynamics(std::string& name, std::vector<Particle
   ,_kPoles(_KmatrixPiPiS->kpoles())
 {
   DebugMsg << "projection index is " << _projectionIndex << endmsg;
+  _isLdependent=false;
 }
 
 PiPiSWaveASDynamics::~PiPiSWaveASDynamics()
diff --git a/PwaUtils/VoigtDynamics.cc b/PwaUtils/VoigtDynamics.cc
index 1919e4e8..7d3499cf 100644
--- a/PwaUtils/VoigtDynamics.cc
+++ b/PwaUtils/VoigtDynamics.cc
@@ -40,6 +40,7 @@ VoigtDynamics::VoigtDynamics(std::string& name, std::vector<Particle*>& fsPartic
   _voigtPtr(new Voigtian())
 {
   if(0!=mother) _massSigmaKey=_massKey+"Sigma";
+  _isLdependent=false;
 }
 
 VoigtDynamics::~VoigtDynamics()
diff --git a/PwaUtils/WoDynamics.cc b/PwaUtils/WoDynamics.cc
index e9ccd55d..e8210d4b 100644
--- a/PwaUtils/WoDynamics.cc
+++ b/PwaUtils/WoDynamics.cc
@@ -35,6 +35,7 @@
 WoDynamics::WoDynamics(std::string& name, std::vector<Particle*>& fsParticles, Particle* mother) :
   AbsDynamics(name, fsParticles, mother)
 {
+  _isLdependent=false;
 }
 
 WoDynamics::~WoDynamics()
-- 
GitLab