diff --git a/Examples/Jamfile b/Examples/Jamfile
index f4b29065016a8ea257a58602360152a5631cb4e0..f4fe0fcbbdf052ec885b3e471ea195b4c0dc76db 100644
--- a/Examples/Jamfile
+++ b/Examples/Jamfile
@@ -10,3 +10,5 @@ build-project JpsiGamKsKlKK ;
build-project JpsiGamEtaPiPiNew ;
build-project JpsiToPhiPhiGam ;
build-project D0ToKsPipPim ;
+build-project JpsiToOmegaPhiGam ;
+
diff --git a/Examples/JpsiToOmegaPhiGam/ConfigFiles/JpsiToOmegaPhiDefault.cfg b/Examples/JpsiToOmegaPhiGam/ConfigFiles/JpsiToOmegaPhiDefault.cfg
deleted file mode 100644
index 2246fba5255afba99cae4bfb8187cece9c842729..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/ConfigFiles/JpsiToOmegaPhiDefault.cfg
+++ /dev/null
@@ -1,32 +0,0 @@
-#################################################################################
-# Configuration file for JpsiToEtaPiPi PWA #
-# #
-# Authors: Bertram Kopf (RUB) #
-#################################################################################
-
-errLogMode = debug
-noOfThreads = 32
-#errLogMode = trace
-
-datFile = /data/liema/pascal/bes/marc/FS4Vectors-JpsiGamPhiOm_daten.dat
-mcFile = /data/liema/pascal/bes/marc/FS4Vectors-JpsiGamPhiOm1G.dat
-
-ratioMcToData = 10
-useEventWeight = false
-
-paramFile = ./defaultparams.dat
-
-startHypo = production
-
-
-enableHyp = GammaEta_2981
-
-
-massRangeMin = 1.
-massRangeMax = 4.
-
-mode = qaMode
-# Determines whether information about configuration options should be emitted
-verbose = 1
-
-#mnParFix = J1P-1C-1Lama0Lamb1GammaEta_2981Phi
diff --git a/Examples/JpsiToOmegaPhiGam/Jamfile b/Examples/JpsiToOmegaPhiGam/Jamfile
index 0baeb6e8d4509d36a655c463d568cbf7193c46b2..5ef09b65611832d0b64b57f8f0cced6415160bf7 100644
--- a/Examples/JpsiToOmegaPhiGam/Jamfile
+++ b/Examples/JpsiToOmegaPhiGam/Jamfile
@@ -15,3 +15,4 @@ exe MJpsiToOmegaPhiGamApp : MJpsiToOmegaPhiGamApp.cc JpsiToOmegaPhiGam : ;
+
diff --git a/Examples/JpsiToOmegaPhiGam/Jamfile~ b/Examples/JpsiToOmegaPhiGam/Jamfile~
deleted file mode 100644
index 348fb7c2df4b1e355da749635e926430934ba46e..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/Jamfile~
+++ /dev/null
@@ -1,17 +0,0 @@
-project :
- ;
-
-lib JpsiToPhiPhiGam :
- [ glob *.cc : *App.cc ]
- $(TOP)/PwaUtils//PwaUtils
- $(TOP)/Event//Event
- $(TOP)/qft++//qft++
- :
- :
- : ;
-
-exe MJpsiToOmegaPhiGamApp : MJpsiToOmegaPhiGamApp.cc JpsiToOmegaPhiGam : ;
-
-
-
-
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiGamPipPimPi0KK.root b/Examples/JpsiToOmegaPhiGam/JpsiGamPipPimPi0KK.root
deleted file mode 100644
index 33c4ff0a6b88d32079beab5499a2720fe5014c6b..0000000000000000000000000000000000000000
Binary files a/Examples/JpsiToOmegaPhiGam/JpsiGamPipPimPi0KK.root and /dev/null differ
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc
index 16ddaf9f73f4ab1e7a74e52274a10af050fd7642..c431245a6fa816951632e46f04908832a16defab 100644
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc
+++ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc
@@ -79,27 +79,56 @@ void JpsiToOmegaPhiGamEventList::read4Vecs(EventList& evtList, std::vector<EvtDa
Vector4<float> V4_PipPimPi0_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab);
Vector4<float> V4_KpKm_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_KpKm_Lab);
+ Vector4<float> V4_Kp_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Kp_Lab);
+ Vector4<float> V4_Km_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Km_Lab);
+
Vector4<float> V4_Pip_HeliPipPimPi0=helicityVec(V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, V4_pip_Lab);
Vector4<float> V4_Pim_HeliPipPimPi0=helicityVec(V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, V4_pim_Lab);
Vector4<float> V4_Pi0_HeliPipPimPi0=helicityVec(V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, V4_pi0_Lab);
Vector4<float> V4_Kp_HeliKpKm=helicityVec(V4_PipPimPi0KpKm_Lab, V4_KpKm_Lab, V4_Kp_Lab);
- Vector4<float> V4_Kp_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Kp_Lab);
- Vector4<float> V4_Km_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Km_Lab);
-
+ //
+ // Omega decay plane normal and lambda parameter
+ //
Vector4<float> V4_omegaDecPlaneNormal_HeliPipPimPi0 ( 0.5*( V4_Pip_HeliPipPimPi0.T()+ V4_Pim_HeliPipPimPi0.T()+ V4_Pi0_HeliPipPimPi0.T() ),
V4_Pim_HeliPipPimPi0.Y()*V4_Pip_HeliPipPimPi0.Z() - V4_Pim_HeliPipPimPi0.Z()*V4_Pip_HeliPipPimPi0.Y(),
V4_Pim_HeliPipPimPi0.Z()*V4_Pip_HeliPipPimPi0.X() - V4_Pim_HeliPipPimPi0.X()*V4_Pip_HeliPipPimPi0.Z(),
V4_Pim_HeliPipPimPi0.X()*V4_Pip_HeliPipPimPi0.Y() - V4_Pim_HeliPipPimPi0.Y()*V4_Pip_HeliPipPimPi0.X() );
-
-
-
double theQ = V4_Pip_HeliPipPimPi0.E()- V4_Pip_HeliPipPimPi0.M() + V4_Pim_HeliPipPimPi0.E()- V4_Pim_HeliPipPimPi0.M() + V4_Pi0_HeliPipPimPi0.E()- V4_Pi0_HeliPipPimPi0.M();
double lambdaNorm=theQ*theQ*(theQ*theQ/108.+ V4_Pim_HeliPipPimPi0.M()*theQ/9.+ V4_Pim_HeliPipPimPi0.M()* V4_Pim_HeliPipPimPi0.M()/3.);
double lambda= V4_omegaDecPlaneNormal_HeliPipPimPi0.P()* V4_omegaDecPlaneNormal_HeliPipPimPi0.P()/lambdaNorm;
+ ///////
+
+
+ //
+ // omega and phi decay plane normal in omega phi helicity system
+ //
+
+ //omega
+ Vector4<float> normalpipluspiminusheliomegaTLVlab= helicityVecInverse( V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, V4_omegaDecPlaneNormal_HeliPipPimPi0);
+ Vector4<float> V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm= helicityVec(V4_all_Lab,V4_PipPimPi0KpKm_Lab , normalpipluspiminusheliomegaTLVlab);
+
+ //test
+ Vector4<float> testBoosted = helicityVec( V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, normalpipluspiminusheliomegaTLVlab );
+ //cout << "++++++++++++++++++++++++++++++"<<endl;
+ //cout << V4_omegaDecPlaneNormal_HeliPipPimPi0 << endl;
+ //cout << testBoosted << endl;
+ //cout << "+++++++++++++++++++++++++++++++" << endl;
+ //eot
+
+
+
+
+
+ //phi
+ Vector4<float> V4_phiDecPlaneNormal_HeliPipPimPi0KpKm ( 0.5*( V4_Kp_HeliPipPimPi0KpKm.T()+ V4_Km_HeliPipPimPi0KpKm.T() ),
+ V4_Km_HeliPipPimPi0KpKm.Y()*V4_Kp_HeliPipPimPi0KpKm.Z() - V4_Km_HeliPipPimPi0KpKm.Z()*V4_Kp_HeliPipPimPi0KpKm.Y(),
+ V4_Km_HeliPipPimPi0KpKm.Z()*V4_Kp_HeliPipPimPi0KpKm.X() - V4_Km_HeliPipPimPi0KpKm.X()*V4_Kp_HeliPipPimPi0KpKm.Z(),
+ V4_Km_HeliPipPimPi0KpKm.X()*V4_Kp_HeliPipPimPi0KpKm.Y() - V4_Km_HeliPipPimPi0KpKm.Y()*V4_Kp_HeliPipPimPi0KpKm.X() );
+ //
EvtDataNew* evtData=new EvtDataNew();
evtData->FourVecsProd[enumProd4V::Psi] = V4_psi;
@@ -121,6 +150,8 @@ void JpsiToOmegaPhiGamEventList::read4Vecs(EventList& evtList, std::vector<EvtDa
evtData->FourVecsDec[enumJpsiGamX4V::V4_Pi0_HeliPipPimPi0]=V4_Pi0_HeliPipPimPi0;
evtData->FourVecsDec[enumJpsiGamX4V::V4_omegaDecPlaneNormal_HeliPipPimPi0] =V4_omegaDecPlaneNormal_HeliPipPimPi0;
+ evtData->FourVecsDec[enumJpsiGamX4V::V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm] =V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm;
+ evtData->FourVecsDec[enumJpsiGamX4V::V4_phiDecPlaneNormal_HeliPipPimPi0KpKm] =V4_phiDecPlaneNormal_HeliPipPimPi0KpKm;
evtData->KinematicVariables[enumJpsiGamXKin::OmegaDecLambda]=lambda;
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc~ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc~
deleted file mode 100644
index f53eaae208131d3c2e14958d5b82a747968dd2ec..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.cc~
+++ /dev/null
@@ -1,163 +0,0 @@
-#include <getopt.h>
-#include <fstream>
-#include <string>
-
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh"
-
-#include "Event/EventList.hh"
-#include "PwaUtils/KinUtils.hh"
-
-#include "Event/Event.hh"
-#include "ErrLogger/ErrLogger.hh"
-
-
-JpsiToOmegaPhiGamEventList::JpsiToOmegaPhiGamEventList() :
- EvtDataBaseListNew()
-{
-}
-
-
-JpsiToOmegaPhiGamEventList::~JpsiToOmegaPhiGamEventList()
-{
-}
-
-void JpsiToOmegaPhiGamEventList::read4Vecs(EventList& evtList, std::vector<EvtDataNew*>& theEvtList, double& evtWeightSum, int maxEvts){
- Event* anEvent;
- int evtCount = 0;
- while ((anEvent = evtList.nextEvent())){
- if (evtCount>= maxEvts) break;
- if (evtCount%10000 == 0) Info << "4vec calculation for event " << evtCount ; // << endmsg;
-
- Vector4<float> gam = *(anEvent->p4(0));
- Vector4<float> kp = *(anEvent->p4(1));
- Vector4<float> km = *(anEvent->p4(2));
- Vector4<float> pip = *(anEvent->p4(3));
- Vector4<float> pim = *(anEvent->p4(4));
- Vector4<float> pi0 = *(anEvent->p4(5));
-
- Vector4<float> V4_psi = gam+kp+km+pip+pim+pi0;
-
- if (evtCount%10000 == 0){
- Info << "psi 4vec" << "\n"
- << " px: " << V4_psi.Px() <<"\t"
- << " py: " << V4_psi.Py() <<"\t"
- << " pz: " << V4_psi.Pz() <<"\t"
- << " e : " << V4_psi.E() << "\t"
- << " m : " << V4_psi.M() ; // << endmsg;
- }
-
- Vector4<float> V4_all_Lab( pip+pim+pi0+kp+km+gam );
- Vector4<float> V4_PipPimPi0KpKm_Lab( pip+pim+pi0+kp+km );
- Vector4<float> V4_PipPimPi0_Lab( pip+pim+pi0 );
- Vector4<float> V4_pip_Lab( pip );
- Vector4<float> V4_KpKm_Lab( kp+km );
- Vector4<float> V4_Kp_Lab( kp );
- Vector4<float> V4_Km_Lab( km );
- Vector4<float> V4_pim_Lab( pim );
- Vector4<float> V4_pi0_Lab( pi0 );
-
- Vector4<float> V4_PipPimPi0KpKm_HeliPsi( V4_PipPimPi0KpKm_Lab );
- V4_PipPimPi0KpKm_HeliPsi.Boost(V4_psi);
-
- Vector4<float> V4_gamma_HeliPsi( gam );
- V4_gamma_HeliPsi.Boost(V4_psi);
- Vector4<float> V4_KpKm_HeliPsi( kp+km );
- V4_KpKm_HeliPsi.Boost(V4_psi);
- Vector4<float> V4_PipPimPi0_HeliPsi( pip+pim+pi0 );
- V4_PipPimPi0_HeliPsi.Boost(V4_psi);
-
- Vector4<float> V4_Pip_HeliPsi( pip );
- V4_Pip_HeliPsi.Boost( V4_psi );
- Vector4<float> V4_Pim_HeliPsi( pim );
- V4_Pim_HeliPsi.Boost( V4_psi );
- Vector4<float> V4_Pi0_HeliPsi( pi0 );
- V4_Pi0_HeliPsi.Boost( V4_psi );
- Vector4<float> V4_Kp_HeliPsi( kp );
- V4_Kp_HeliPsi.Boost( V4_psi );
- Vector4<float> V4_Km_HeliPsi( km );
- V4_Km_HeliPsi.Boost( V4_psi );
-
- Vector4<float> V4_PipPimPi0_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab);
- Vector4<float> V4_KpKm_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_KpKm_Lab);
- Vector4<float> V4_Pip_HeliPipPimPi0=helicityVec(V4_PipPimPi0KpKm_Lab, V4_PipPimPi0_Lab, V4_Pip_Lab);
- Vector4<float> V4_Kp_HeliKpKm=helicityVec(V4_PipPimPi0KpKm_Lab, V4_KpKm_Lab, V4_Kp_Lab);
-
- Vector4<float> V4_Kp_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Kp_Lab);
- Vector4<float> V4_Km_HeliPipPimPi0KpKm=helicityVec(V4_all_Lab, V4_PipPimPi0KpKm_Lab, V4_Km_Lab);
-
- // Vector4<float> V4_normKpKmDecHeliKsKlKpKm
- // (0.5*(V4_Kp_HeliKsKlKpKm.T()+V4_Km_HeliKsKlKpKm.T()),
- // V4_Km_HeliKsKlKpKm.Y()*V4_Kp_HeliKsKlKpKm.Z()-V4_Km_HeliKsKlKpKm.Z()*V4_Kp_HeliKsKlKpKm.Y(),
- // V4_Km_HeliKsKlKpKm.Z()*V4_Kp_HeliKsKlKpKm.X()-V4_Km_HeliKsKlKpKm.X()*V4_Kp_HeliKsKlKpKm.Z(),
- // V4_Km_HeliKsKlKpKm.X()*V4_Kp_HeliKsKlKpKm.Y()-V4_Km_HeliKsKlKpKm.Y()*V4_Kp_HeliKsKlKpKm.X());
-
- // Vector4<float> V4_Kl_HeliKsKlKpKm=helicityVec(V4_all_Lab, V4_KsKlKpKm_Lab, V4_Kl_Lab);
- // Vector4<float> V4_Ks_HeliKsKlKpKm=helicityVec(V4_all_Lab, V4_KsKlKpKm_Lab, V4_Ks_Lab);
-
- // Vector4<float> V4_normKsKlDecHeliKsKlKpKm
- // (0.5*(V4_Kl_HeliKsKlKpKm.T()+V4_Ks_HeliKsKlKpKm.T()),
- // V4_Ks_HeliKsKlKpKm.Y()*V4_Kl_HeliKsKlKpKm.Z()-V4_Ks_HeliKsKlKpKm.Z()*V4_Kl_HeliKsKlKpKm.Y(),
- // V4_Ks_HeliKsKlKpKm.Z()*V4_Kl_HeliKsKlKpKm.X()-V4_Ks_HeliKsKlKpKm.X()*V4_Kl_HeliKsKlKpKm.Z(),
- // V4_Ks_HeliKsKlKpKm.X()*V4_Kl_HeliKsKlKpKm.Y()-V4_Ks_HeliKsKlKpKm.Y()*V4_Kl_HeliKsKlKpKm.X());
-
-
- EvtDataNew* evtData=new EvtDataNew();
- evtData->FourVecsProd[enumProd4V::Psi] = V4_psi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0KpKm_HeliPsi]=V4_PipPimPi0KpKm_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPsi]=V4_PipPimPi0_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPsi]=V4_KpKm_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_gamma_HeliPsi]=V4_gamma_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPsi]=V4_Pip_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Pim_HeliPsi]=V4_Pim_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Pi0_HeliPsi]=V4_Pi0_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliPsi]=V4_Kp_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Km_HeliPsi]=V4_Km_HeliPsi;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPipPimPi0KpKm]=V4_PipPimPi0_HeliPipPimPi0KpKm;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm]=V4_KpKm_HeliPipPimPi0KpKm;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Ks_HeliPipPimPi0]=V4_Ks_HeliPipPimPi0;
- evtData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliKpKm]= V4_Kp_HeliKpKm;
-
- // evtData->FourVecsDec[enumJpsiGamX4V::V4_normKpKmDecHeliKsKlKpKm]=V4_normKpKmDecHeliKsKlKpKm;
- // evtData->FourVecsDec[enumJpsiGamX4V::V4_normKsKlDecHeliKsKlKpKm]=V4_normKsKlDecHeliKsKlKpKm;
-
-
- // calculate and store WignerD functions for Psi -> X gamma (JPC =0-+, 0++, 2++)
-
- Spin jPsi=1;
- for (Spin M=-1; M<=1; M=M+2){
- for (Spin lam=-1; lam<=1; lam++){
- evtData->WignerDsProd[enumProdDfunc::Psi][jPsi][M][lam]=Wigner_D(0.,V4_PipPimPi0KpKm_HeliPsi.Theta(),0,jPsi,M,lam);
- }
- }
-
- // //WignerD functions for X -> phi phi
-
- // for (Spin J_X=0; J_X<=2; J_X++){
- // for (Spin lam_X=-J_X; lam_X<=J_X; lam_X++){
- // for (Spin lamPhi1mlamPhi2=-J_X; lamPhi1mlamPhi2<=J_X; lamPhi1mlamPhi2++){
- // evtData->WignerDsDec[enumJpsiGamXDfunc::Df_XToPhiPhi1][J_X][lam_X][lamPhi1mlamPhi2]
- // =Wigner_D(V4_KsKl_HeliKsKlKpKm.Phi(),V4_KsKl_HeliKsKlKpKm.Theta(),0,J_X,lam_X,lamPhi1mlamPhi2);
- // evtData->WignerDsDec[enumJpsiGamXDfunc::Df_XToPhiPhi2][J_X][lam_X][lamPhi1mlamPhi2]
- // =Wigner_D(V4_KpKm_HeliKsKlKpKm.Phi(),V4_KpKm_HeliKsKlKpKm.Theta(),0,J_X,lam_X,lamPhi1mlamPhi2);
- // }
- // }
- // }
-
- //WignerD function for phi -> K+ K- and phi -> KS KL
- Spin phiSpin=1;
- for(Spin M=-phiSpin; M<=phiSpin; M++){
- Spin lam=0;
- evtData->WignerDsDec[enumJpsiGamXDfunc::Df_PhiToKpKm][phiSpin][M][lam] = Wigner_D(V4_Kp_HeliKpKm.Phi(),V4_Kp_HeliKpKm.Theta(), 0,phiSpin,M,lam);
- }
-
- evtData->evtWeight=anEvent->Weight();
- evtData->evtNo=_evtNoAll;
- theEvtList.push_back(evtData);
-
- evtWeightSum += anEvent->Weight();
- ++evtCount;
- ++_evtNoAll;
- }
-}
-
-
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh
index 83bf76731752be38e49461cc2dd5b14e3da4a566..4557386d4acc81876bc7e164cf7086953cf1a0b3 100644
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh
+++ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh
@@ -15,7 +15,8 @@ struct enumJpsiGamX4V{
V4_PipPimPi0_HeliPipPimPi0KpKm, V4_KpKm_HeliPipPimPi0KpKm,
V4_Pip_HeliPipPimPi0, V4_Pim_HeliPipPimPi0, V4_Pi0_HeliPipPimPi0,
V4_Kp_HeliKpKm,
- V4_omegaDecPlaneNormal_HeliPipPimPi0, n4Vecs};
+ V4_omegaDecPlaneNormal_HeliPipPimPi0, V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm, V4_phiDecPlaneNormal_HeliPipPimPi0KpKm,
+ n4Vecs};
static const std::string& name(unsigned int t)
{
@@ -25,7 +26,7 @@ struct enumJpsiGamX4V{
"PipPimPi0_HeliPipPimPi0KpKm", "KpKm_HeliPipPimPi0KpKm",
"Pip_HeliPipPimPi0", "Pim_HeliPipPimPi0", "Pi0_HeliPipPimPi0"
"Kp_HeliKpKm",
- "V4_omegaDecPlaneNormal_HeliPipPimPi0" };
+ "V4_omegaDecPlaneNormal_HeliPipPimPi0","V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm", "V4_phiDecPlaneNormal_HeliPipPimPi0KpKm" };
if (t<0 || t>=enumJpsiGamX4V::n4Vecs) assert(0);
return fitName[t];
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh~ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh~
deleted file mode 100644
index 8a001003a5ab828546f994f5ff4b8d68067fe6ef..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh~
+++ /dev/null
@@ -1,75 +0,0 @@
-#ifndef _JpsiToOmegaPhiGamEventList_H
-#define _JpsiToOmegaPhiGamEventList_H
-
-#include <iostream>
-#include <vector>
-
-#include <cassert>
-// #include <TSystem.h>
-#include "qft++/topincludes/relativistic-quantum-mechanics.hh"
-#include "PwaUtils/EvtDataBaseListNew.hh"
-
-struct enumJpsiGamX4V{
- enum {V4_KsKlKpKm_HeliPsi=0, V4_KsKl_HeliPsi,V4_KpKm_HeliPsi, V4_gamma_HeliPsi,
- V4_Ks_HeliPsi, V4_Kl_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi,
- V4_KsKl_HeliKsKlKpKm,V4_KpKm_HeliKsKlKpKm ,V4_Ks_HeliKsKl,V4_Kp_HeliKpKm,
- V4_normKpKmDecHeliKsKlKpKm, V4_normKsKlDecHeliKsKlKpKm, n4Vecs};
-
- static const std::string& name(unsigned int t)
- {
- static std::string fitName[enumJpsiGamX4V::n4Vecs]
- ={"KsKlKpKm_HeliPsi", "KsKl_HeliPsi,KpKm_HeliPsi", "gamma_HeliPsi",
- "Ks_HeliPsi, Kl_HeliPsi", "Kp_HeliPsi", "Km_HeliPsi",
- "KsKl_HeliKsKlKpKm","KpKm_HeliKsKlKpKm" ,"Ks_HeliKsKl","Kp_HeliKpKm",
- "normKpKmDecHeliKsKlKpKm", "normKsKlDecHeliKsKlKpKm"
- };
-
- if (t<0 || t>=enumJpsiGamX4V::n4Vecs) assert(0);
- return fitName[t];
- }
-};
-
-struct enumJpsiGamXDfunc{
- enum {Df_XToPhiPhi1=0, Df_XToPhiPhi2, Df_PhiToKsKl,Df_PhiToKpKm, nDfuncts};
-
- static const std::string& name(unsigned int t)
- {
- static std::string fitName[enumJpsiGamXDfunc::nDfuncts]
- ={"XToPhiPhi1", "XToPhiPhi2", "PhiToKsKl", "PhiToKpKm"};
- if (t<0 || t>=enumJpsiGamXDfunc::nDfuncts) assert(0);
- return fitName[t];
- }
-};
-
-
-
-class EventList;
-
-class JpsiToOmegaPhiGamEventList : public EvtDataBaseListNew {
-
-public:
-
- // create/copy/destroy:
-
- ///Constructor
- JpsiToOmegaPhiGamEventList();
-
-
-
- /** Destructor */
- virtual ~JpsiToOmegaPhiGamEventList();
-
- // Getters:
-
-
-protected:
-
-
- virtual void read4Vecs(EventList& evtList, std::vector<EvtDataNew*>& theEvtList, double& evtWeightSum, int maxEvts);
-
-private:
-
-
-};
-
-#endif
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc
index e7522fe6fd86875b48600fbf74a76f3b031dfa09..8579868db2e392ce9a9ecb93172d25f61dec5f2f 100644
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc
+++ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc
@@ -334,9 +334,9 @@ void JpsiToOmegaPhiGamHist::initRootStuff()
_phiPhi_OmegaPhiHeliMcHist= new TH1F("_phiPhi_OmegaPhiHeliMcHist", "cos(#Phi_{#phi}) K+ K- Mc", 100, -TMath::Pi(), TMath::Pi());
_phiPhi_OmegaPhiHeliFittedHist= new TH1F("_phiPhi_OmegaPhiHeliFittedHist", "cos(#Phi_{#phi}) K+ K- Fit", 100, -TMath::Pi(), TMath::Pi());
- _chiDataHist= new TH1F("_chiDataHist", "#chi data", 90, 0., 90.);
- _chiMcHist= new TH1F("_chiMcHist", "#chi Mc", 90, 0., 90.);
- _chiFittedHist= new TH1F("_chiFittedHist", "#chi Fit", 90, 0., 90.);
+ _chiDataHist= new TH1F("_chiDataHist", "#chi data", 100, 0., 90.);
+ _chiMcHist= new TH1F("_chiMcHist", "#chi Mc", 100, 0., 90.);
+ _chiFittedHist= new TH1F("_chiFittedHist", "#chi Fit", 100, 0., 90.);
nbins = 50;
xmin=0;
@@ -411,25 +411,20 @@ void JpsiToOmegaPhiGamHist::plotLambda(TH1F* theHisto, EvtDataNew* theData, doub
}
void JpsiToOmegaPhiGamHist::plotChi(TH1F* theChiHisto, EvtDataNew* theData, double weight){
-// Vector4<double>& V4_PipPimPi0KpKm_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0KpKm_HeliPsi] ;
-// Vector4<double>& V4_Pip_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPsi] ;
-// Vector4<double>& V4_Kl_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kl_HeliPsi] ;
-// Vector4<double>& V4_Kp_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliPsi] ;
-// Vector4<double>& V4_Km_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Km_HeliPsi] ;
-
-// double thePhiPhiDecayPlaneAngle = decayAngleChi( V4_PipPimPi0KpKm_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi, V4_Pip_HeliPsi, V4_Kl_HeliPsi );
-
-
- //Vector4<double>& V4_normKpKmDecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normKpKmDecHeliPipPimPi0KpKm];
- //Vector4<double>& V4_normPipPimPi0DecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normPipPimPi0DecHeliPipPimPi0KpKm];
-
- // double cosChi=(V4_normKpKmDecHeliPipPimPi0KpKm.Px()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Px()
- // +V4_normKpKmDecHeliPipPimPi0KpKm.Py()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Py()
- // +V4_normKpKmDecHeliPipPimPi0KpKm.Pz()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Pz())
- // / (V4_normKpKmDecHeliPipPimPi0KpKm.P()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.P());
-
- // double chi=acos(fabs(cosChi));
- // theChiHisto->Fill(chi*180./TMath::Pi(),weight);
+
+ Vector4<double>& V4_normKpKmDecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_phiDecPlaneNormal_HeliPipPimPi0KpKm];
+ Vector4<double>& V4_normPipPimPi0DecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_omegaDecPlaneNormal_HeliPipPimPi0KpKm];
+
+ double cosChi=(V4_normKpKmDecHeliPipPimPi0KpKm.Px()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Px()
+ +V4_normKpKmDecHeliPipPimPi0KpKm.Py()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Py()
+ +V4_normKpKmDecHeliPipPimPi0KpKm.Pz()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Pz())
+ / (V4_normKpKmDecHeliPipPimPi0KpKm.P()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.P());
+
+ double chi=acos(fabs(cosChi));
+ if(chi>0.5*TMath::Pi()){
+ chi=TMath::Pi()-chi;
+ }
+ theChiHisto->Fill(chi*180./TMath::Pi(),weight);
}
void JpsiToOmegaPhiGamHist::fillTuple( TNtuple* theTuple, EvtDataNew* theData, double weight)
diff --git a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc~ b/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc~
deleted file mode 100644
index 0cd9ac3601e2e564badc3f6acc1484684aee1077..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.cc~
+++ /dev/null
@@ -1,672 +0,0 @@
-#include <getopt.h>
-#include <fstream>
-#include <sstream>
-#include <string>
-
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.hh"
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh"
-//#include "Examples/JpsiGamPipPimPi0KK/JpsiGamPipPimPi0KKData.hh"
-//#include "Examples/JpsiGamPipPimPi0KK/JpsiGamPipPimPi0KKProdLh.hh"
-#include "PwaUtils/KinUtils.hh"
-//#include "Examples/JpsiGamPipPimPi0KK/JpsiGamPipPimPi0KKFitParams.hh"
-
-#include "TFile.h"
-#include "TH1F.h"
-#include "TH2F.h"
-#include "TMath.h"
-#include "TNtuple.h"
-#include "TLorentzVector.h"
-#include "ErrLogger/ErrLogger.hh"
-
-
-JpsiToOmegaPhiGamHist::JpsiToOmegaPhiGamHist(boost::shared_ptr<const EvtDataBaseListNew> theEvtList) :
- _theTFile(0),
- _dalitzDataHist(0),
- _dalitzMcHist(0),
- _dalitzFittedHist(0),
- _PhiPhiMassDataHist(0),
- _PhiPhiMassMcHist(0),
- _PhiPhiMassFittedHist(0),
- _KpKmMassDataHist(0),
- _KpKmMassMcHist(0),
- _KpKmMassFittedHist(0),
- _PipPimPi0MassDataHist(0),
- _PipPimPi0MassMcHist(0),
- _PipPimPi0MassFittedHist(0),
- _costPhi_KpKmDataHist(0),
- _costPhi_KpKmMcHist(0),
- _costPhi_KpKmFittedHist(0),
- _phiPhi_KpKmDataHist(0),
- _phiPhi_KpKmMcHist(0),
- _phiPhi_KpKmFittedHist(0),
- _chiDataHist(0),
- _chiMcHist(0),
- _chiFittedHist(0),
- _dataTuple(0),
- _mcTuple(0),
- _massIndepTuple(0),
- _massRange(make_pair(0.,0.) )
-{
- if(0==theEvtList){
- Alert <<"JpsiGamPipPimPi0KKEventList* theEvtList is a 0 pointer !!!!" ; // << endmsg;
- exit(1);
- }
-
- initRootStuff();
-
- const std::vector<EvtDataNew*> dataList=theEvtList->getDataVecs();
-
- std::vector<EvtDataNew*>::const_iterator it=dataList.begin();
- while(it!=dataList.end())
- {
- const double evtWeight=(*it)->evtWeight;
- plotDalitz(_dalitzDataHist, (*it), evtWeight);
- plotPhiPhi(_PhiPhiMassDataHist, (*it), evtWeight );
- plotPipPimPi0( _PipPimPi0MassDataHist, (*it), evtWeight );
- plotKpKm( _KpKmMassDataHist, (*it), evtWeight );
- plotCostPhiPip( _costPip_PipPimPi0HeliDataHist, _phiPip_PipPimPi0HeliDataHist, (*it), evtWeight );
- plotCostPhiKp( _costKp_KpKmHeliDataHist, _phiKp_KpKmHeliDataHist, (*it), evtWeight );
- plotCostGam( _costGamCmDataHist, (*it), evtWeight );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmDataHist, _phiPhi_KpKmDataHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], evtWeight);
- plotChi(_chiDataHist, (*it), evtWeight );
- fillTuple(_dataTuple, (*it), evtWeight);
-
- ++it;
- }
-
- const std::vector<EvtDataNew*> mcList=theEvtList->getMcVecs();
- it=mcList.begin();
- while(it!=mcList.end())
- {
- plotDalitz(_dalitzMcHist, (*it), 1.);
- plotPhiPhi(_PhiPhiMassMcHist, (*it), 1. );
- plotPipPimPi0( _PipPimPi0MassMcHist, (*it), 1. );
- plotKpKm( _KpKmMassMcHist, (*it), 1. );
- plotCostPhiPip( _costPip_PipPimPi0HeliMcHist, _phiPip_PipPimPi0HeliMcHist, (*it), 1. );
- plotCostPhiKp( _costKp_KpKmHeliMcHist, _phiKp_KpKmHeliMcHist,(*it), 1. );
- plotCostGam( _costGamCmMcHist,(*it), 1. );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmMcHist, _phiPhi_KpKmMcHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], 1.);
- plotChi(_chiMcHist, (*it), 1. );
-
- // double evtWeight= _theJpsiGamPipPimPi0KKLh->calcEvtIntensity((*it), _fitParam);
- double evtWeight=1.;
- plotDalitz(_dalitzFittedHist, (*it),evtWeight );
- plotPhiPhi(_PhiPhiMassFittedHist, (*it), evtWeight );
- plotPipPimPi0( _PipPimPi0MassFittedHist, (*it), evtWeight );
- plotKpKm( _KpKmMassFittedHist, (*it), evtWeight );
- plotCostPhiPip( _costPip_PipPimPi0HeliFittedHist, _phiPip_PipPimPi0HeliFittedHist,(*it), evtWeight );
- plotCostPhiKp( _costKp_KpKmHeliFittedHist, _phiKp_KpKmHeliFittedHist,(*it), evtWeight );
- plotCostGam( _costGamCmFittedHist,(*it), evtWeight );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmFittedHist, _phiPhi_KpKmFittedHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], evtWeight);
- plotChi(_chiFittedHist, (*it), evtWeight );
-
- fillTuple(_mcTuple, (*it), evtWeight);
- ++it;
- }
-
- double integralDataWoWeight=(double) theEvtList->getDataVecs().size();
- Info <<"No of data events without weight " << integralDataWoWeight ; // << endmsg;
-
- double integralDataWWeight = theEvtList->NoOfWeightedDataEvts();
- Info <<"No of data events with weight " << integralDataWWeight ; // << endmsg;
-
- double integralMC = theEvtList->NoOfWeightedMcEvts();
- Info <<"No of MC events " << integralMC ; // << endmsg;
-
- Info <<"scaling factor " << integralDataWWeight/integralMC ; // << endmsg;
-
- double scaleFactor = integralDataWWeight/integralMC;
-
- _dalitzFittedHist->Scale(scaleFactor);
- _PhiPhiMassFittedHist->Scale(scaleFactor);
- _PipPimPi0MassFittedHist->Scale(scaleFactor);
- _KpKmMassFittedHist->Scale(scaleFactor);
- _costPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
- _phiPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
- _costKp_KpKmHeliFittedHist->Scale(scaleFactor);
- _phiKp_KpKmHeliFittedHist->Scale(scaleFactor);
- _costGamCmFittedHist->Scale(scaleFactor);
- _costPhi_KpKmFittedHist->Scale(scaleFactor);
- _phiPhi_KpKmFittedHist->Scale(scaleFactor);
- _chiFittedHist->Scale(scaleFactor);
-
-
-}
-
-JpsiToOmegaPhiGamHist::JpsiToOmegaPhiGamHist(boost::shared_ptr<AbsLhNew> theLh, fitParamsNew& theFitParms):
- _theTFile(0),
- _dalitzDataHist(0),
- _dalitzMcHist(0),
- _dalitzFittedHist(0),
- _PhiPhiMassDataHist(0),
- _PhiPhiMassMcHist(0),
- _PhiPhiMassFittedHist(0),
- _KpKmMassDataHist(0),
- _KpKmMassMcHist(0),
- _KpKmMassFittedHist(0),
- _PipPimPi0MassDataHist(0),
- _PipPimPi0MassMcHist(0),
- _PipPimPi0MassFittedHist(0),
- _costPhi_KpKmDataHist(0),
- _costPhi_KpKmMcHist(0),
- _costPhi_KpKmFittedHist(0),
- _phiPhi_KpKmDataHist(0),
- _phiPhi_KpKmMcHist(0),
- _phiPhi_KpKmFittedHist(0),
- _chiDataHist(0),
- _chiMcHist(0),
- _chiFittedHist(0),
- _dataTuple(0),
- _mcTuple(0),
- _massIndepTuple(0),
- _massRange(make_pair(0.,0.) )
-{
- if(0==theLh){
- Alert <<"AbsLh* is a 0 pointer !!!!" ; // << endmsg;
- exit(1);
- }
-
- initRootStuff();
- boost::shared_ptr<const EvtDataBaseListNew> theEvtList=theLh->getEventList();
- const std::vector<EvtDataNew*> dataList=theEvtList->getDataVecs();
-
- std::vector<EvtDataNew*>::const_iterator it=dataList.begin();
- while(it!=dataList.end())
- {
- const double evtWeight=(*it)->evtWeight;
- plotDalitz(_dalitzDataHist, (*it), evtWeight);
- plotPhiPhi(_PhiPhiMassDataHist, (*it), evtWeight );
- plotPipPimPi0( _PipPimPi0MassDataHist, (*it), evtWeight );
- plotKpKm( _KpKmMassDataHist, (*it), evtWeight );
- plotCostPhiPip( _costPip_PipPimPi0HeliDataHist, _phiPip_PipPimPi0HeliDataHist, (*it), evtWeight );
- plotCostPhiKp( _costKp_KpKmHeliDataHist, _phiKp_KpKmHeliDataHist, (*it), evtWeight );
- plotCostGam( _costGamCmDataHist, (*it), evtWeight );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmDataHist, _phiPhi_KpKmDataHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], evtWeight);
- plotChi(_chiDataHist, (*it), evtWeight );
- fillTuple(_dataTuple, (*it), evtWeight);
-
- ++it;
- }
-
- const std::vector<EvtDataNew*> mcList=theEvtList->getMcVecs();
- it=mcList.begin();
- while(it!=mcList.end())
- {
- plotDalitz(_dalitzMcHist, (*it), 1.);
- plotPhiPhi(_PhiPhiMassMcHist, (*it), 1. );
- plotPipPimPi0( _PipPimPi0MassMcHist, (*it), 1. );
- plotKpKm( _KpKmMassMcHist, (*it), 1. );
- plotCostPhiPip( _costPip_PipPimPi0HeliMcHist, _phiPip_PipPimPi0HeliMcHist, (*it), 1. );
- plotCostPhiKp( _costKp_KpKmHeliMcHist, _phiKp_KpKmHeliMcHist,(*it), 1. );
- plotCostGam( _costGamCmMcHist,(*it), 1. );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmMcHist, _phiPhi_KpKmMcHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], 1.);
- plotChi(_chiMcHist, (*it), 1. );
-
- // double evtWeight= _theJpsiGamPipPimPi0KKLh->calcEvtIntensity((*it), _fitParam);
- double evtWeight= theLh->calcEvtIntensity( (*it), theFitParms );
- plotDalitz(_dalitzFittedHist, (*it),evtWeight );
- plotPhiPhi(_PhiPhiMassFittedHist, (*it), evtWeight );
- plotPipPimPi0( _PipPimPi0MassFittedHist, (*it), evtWeight );
- plotKpKm( _KpKmMassFittedHist, (*it), evtWeight );
- plotCostPhiPip( _costPip_PipPimPi0HeliFittedHist, _phiPip_PipPimPi0HeliFittedHist,(*it), evtWeight );
- plotCostPhiKp( _costKp_KpKmHeliFittedHist, _phiKp_KpKmHeliFittedHist,(*it), evtWeight );
- plotCostGam( _costGamCmFittedHist,(*it), evtWeight );
- plotCostPhi_PhiPhiHeli(_costPhi_KpKmFittedHist, _phiPhi_KpKmFittedHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], evtWeight);
- plotChi(_chiFittedHist, (*it), evtWeight );
-
- fillTuple(_mcTuple, (*it), evtWeight);
- ++it;
- }
-
- double integralDataWoWeight=(double) theEvtList->getDataVecs().size();
- Info <<"No of data events without weight " << integralDataWoWeight ; // << endmsg;
-
- double integralDataWWeight=(double) theEvtList->NoOfWeightedDataEvts();
- Info <<"No of data events with weight " << integralDataWWeight ; // << endmsg;
-
-
- double integralMC=(double) theEvtList->NoOfWeightedMcEvts();
- Info <<"No of MC events " << integralMC ; // << endmsg;
-
- Info <<"scaling factor " << integralDataWWeight/integralMC ; // << endmsg;
-
- double scaleFactor = integralDataWWeight/integralMC;
-
- _dalitzFittedHist->Scale(scaleFactor);
- _PhiPhiMassFittedHist->Scale(scaleFactor);
- _PipPimPi0MassFittedHist->Scale(scaleFactor);
- _KpKmMassFittedHist->Scale(scaleFactor);
- _costPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
- _phiPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
- _costKp_KpKmHeliFittedHist->Scale(scaleFactor);
- _phiKp_KpKmHeliFittedHist->Scale(scaleFactor);
- _costGamCmFittedHist->Scale(scaleFactor);
- _costPhi_KpKmFittedHist->Scale(scaleFactor);
- _phiPhi_KpKmFittedHist->Scale(scaleFactor);
- _chiFittedHist->Scale(scaleFactor);
-
- double integralFitted=(double) _costGamCmFittedHist->Integral();
- Info <<"No of fit events " << integralFitted ; // << endmsg;
-}
-
-
-// JpsiToOmegaPhiGamHist::JpsiToOmegaPhiGamHist(JpsiGamPipPimPi0KKProdLh* theJpsiGamPipPimPi0KKLh, fitParams& fitParam, FitParamErrorMatrix* theErrorMatrix) :
-// _theTFile(0),
-// _dalitzDataHist(0),
-// _dalitzMcHist(0),
-// _dalitzFittedHist(0),
-// _PhiPhiMassDataHist(0),
-// _PhiPhiMassMcHist(0),
-// _PhiPhiMassFittedHist(0),
-// _KpKmMassDataHist(0),
-// _KpKmMassMcHist(0),
-// _KpKmMassFittedHist(0),
-// _PipPimPi0MassDataHist(0),
-// _PipPimPi0MassMcHist(0),
-// _PipPimPi0MassFittedHist(0),
-// _costPhi_KpKmDataHist(0),
-// _costPhi_KpKmMcHist(0),
-// _costPhi_KpKmFittedHist(0),
-// _phiPhi_KpKmDataHist(0),
-// _phiPhi_KpKmMcHist(0),
-// _phiPhi_KpKmFittedHist(0),
-// _chiDataHist(0),
-// _chiMcHist(0),
-// _chiFittedHist(0),
-// _dataTuple(0),
-// _mcTuple(0),
-// _massIndepTuple(0),
-// _massRange(make_pair(0.,0.) )
-// {
-// if(0==theJpsiGamPipPimPi0KKLh){
-// Alert <<"JpsiGamPipPimPi0KKLh* theJpsiGamPipPimPi0KKLh is a 0 pointer !!!!" ; // << endmsg;
-// exit(1);
-// }
-
-// initRootStuff();
-// _theJpsiGamPipPimPi0KKLh = theJpsiGamPipPimPi0KKLh;
-// _fitParam = fitParam;
-// std::vector<double> data;
-// _errMatrix = theErrorMatrix;
-
-// }
-
-
-// void JpsiToOmegaPhiGamHist::fill(){
-
-// boost::shared_ptr<const EvtDataBaseList> theEvtList=_theJpsiGamPipPimPi0KKLh->getEventList();
-// const std::vector<EvtData*> dataList=theEvtList->getDataVecs();
-
-// std::vector<EvtDataNew*>::const_iterator it=dataList.begin();
-// while(it!=dataList.end())
-// {
-// plotDalitz(_dalitzDataHist, (*it), 1.);
-// plotPhiPhi(_PhiPhiMassDataHist, (*it), 1. );
-// plotPipPimPi0( _PipPimPi0MassDataHist, (*it), 1. );
-// plotKpKm( _KpKmMassDataHist, (*it), 1. );
-// plotCostPhiPip( _costPip_PipPimPi0HeliDataHist, _phiPip_PipPimPi0HeliDataHist, (*it), 1. );
-// plotCostPhiKp( _costKp_KpKmHeliDataHist, _phiKp_KpKmHeliDataHist, (*it), 1. );
-// plotCostGam( _costGamCmDataHist, (*it), 1. );
-// plotCostPhi_PhiPhiHeli(_costPhi_KpKmDataHist, _phiPhi_KpKmDataHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], 1.);
-// plotChi(_chiDataHist, (*it), 1. );
-// fillTuple(_dataTuple, (*it), 1.);
-
-// ++it;
-// }
-
-// const std::vector<EvtDataNew*> mcList=theEvtList->getMcVecs();
-// it=mcList.begin();
-// while(it!=mcList.end())
-// {
-// plotDalitz(_dalitzMcHist, (*it), 1.);
-// plotPhiPhi(_PhiPhiMassMcHist, (*it), 1. );
-// plotPipPimPi0( _PipPimPi0MassMcHist, (*it), 1. );
-// plotKpKm( _KpKmMassMcHist, (*it), 1. );
-// plotCostPhiPip( _costPip_PipPimPi0HeliMcHist, _phiPip_PipPimPi0HeliMcHist, (*it), 1. );
-// plotCostPhiKp( _costKp_KpKmHeliMcHist, _phiKp_KpKmHeliMcHist,(*it), 1. );
-// plotCostGam( _costGamCmMcHist,(*it), 1. );
-// plotCostPhi_PhiPhiHeli(_costPhi_KpKmMcHist, _phiPhi_KpKmMcHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], 1.);
-// plotChi(_chiMcHist, (*it), 1. );
-
-// double evtWeight= _theJpsiGamPipPimPi0KKLh->calcEvtIntensity((*it), _fitParam);
-// plotDalitz(_dalitzFittedHist, (*it),evtWeight );
-// plotPhiPhi(_PhiPhiMassFittedHist, (*it), evtWeight );
-// plotPipPimPi0( _PipPimPi0MassFittedHist, (*it), evtWeight );
-// plotKpKm( _KpKmMassFittedHist, (*it), evtWeight );
-// plotCostPhiPip( _costPip_PipPimPi0HeliFittedHist, _phiPip_PipPimPi0HeliFittedHist,(*it), evtWeight );
-// plotCostPhiKp( _costKp_KpKmHeliFittedHist, _phiKp_KpKmHeliFittedHist,(*it), evtWeight );
-// plotCostGam( _costGamCmFittedHist,(*it), evtWeight );
-// plotCostPhi_PhiPhiHeli(_costPhi_KpKmFittedHist, _phiPhi_KpKmFittedHist, (*it)->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm], evtWeight);
-// plotChi(_chiFittedHist, (*it), evtWeight );
-
-// fillTuple(_mcTuple, (*it), evtWeight);
-// ++it;
-// }
-
-
-// double integralData=(double) theEvtList->getDataVecs().size();
-// Info <<"No of fit data events " << integralData ; // << endmsg;
-
-// double integralFitted=(double) theEvtList->getMcVecs().size();
-// Info <<"No of fit events " << integralFitted ; // << endmsg;
-
-// Info <<"scaling factor " << integralData/integralFitted ; // << endmsg;
-
-// double scaleFactor = integralData/integralFitted;
-
-// _dalitzFittedHist->Scale(scaleFactor);
-// _PhiPhiMassFittedHist->Scale(scaleFactor);
-// _PipPimPi0MassFittedHist->Scale(scaleFactor);
-// _KpKmMassFittedHist->Scale(scaleFactor);
-// _costPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
-// _phiPip_PipPimPi0HeliFittedHist->Scale(scaleFactor);
-// _costKp_KpKmHeliFittedHist->Scale(scaleFactor);
-// _phiKp_KpKmHeliFittedHist->Scale(scaleFactor);
-// _costGamCmFittedHist->Scale(scaleFactor);
-// _costPhi_KpKmFittedHist->Scale(scaleFactor);
-// _phiPhi_KpKmFittedHist->Scale(scaleFactor);
-// _chiFittedHist->Scale(scaleFactor);
-
-
-// double iEta(0.), iEtaErr(0.), iF0(0.), iF0Err(0.), iEta2(0.), iEta2Err(0.), iF1(0.), iF1Err(0.), iF2(0.), iF2Err(0.);
-// double etaReal(0.), etaImg(0.);
-
-// it=mcList.begin();
-// while(it!=mcList.end()){
-
-// std::pair<double, double> intensityEvent = make_pair(0.,0.);
-// _theJpsiGamPipPimPi0KKLh->calcComponentIntensity((*it), _fitParam, *_errMatrix, "etacHyp", intensityEvent);
-// iEta+= intensityEvent.first*scaleFactor;
-// iEtaErr+= intensityEvent.second*scaleFactor;
-
-// _theJpsiGamPipPimPi0KKLh->calcComponentIntensity((*it), _fitParam, *_errMatrix, "f02020Hyp",intensityEvent);
-// iF0 += intensityEvent.first*scaleFactor;
-// iF0Err += intensityEvent.second*scaleFactor;
-
-// _theJpsiGamPipPimPi0KKLh->calcComponentIntensity((*it), _fitParam, *_errMatrix, "f1Hyp", intensityEvent);
-// iF1 += intensityEvent.first*scaleFactor;
-// iF1Err += intensityEvent.second*scaleFactor;
-
-// _theJpsiGamPipPimPi0KKLh->calcComponentIntensity((*it), _fitParam, *_errMatrix, "f22300Hyp",intensityEvent);
-// iF2 += intensityEvent.first*scaleFactor;
-// iF2Err += intensityEvent.second*scaleFactor;
-
-// _theJpsiGamPipPimPi0KKLh->calcComponentIntensity((*it), _fitParam, *_errMatrix, "eta21870Hyp",intensityEvent);
-// iEta2+= intensityEvent.first*scaleFactor;
-// iEta2Err += intensityEvent.second*scaleFactor;
-
-
-
-
-// it++;
-// }
-// double meanMassRange = _massRange.first + 0.5*(_massRange.second-_massRange.first);
-
-
-// Info << "Events for component eta : " << iEta << " +/- " << iEtaErr ;
-// Info << "Events for component f0: " << iF0 << " +/- " << iF0Err ;
-// Info << "Events for component f1: " << iF1 << " +/- " << iF1Err ;
-// Info << "Events for component f2: " << iF2 << " +/- " << iF2Err ;
-// Info << "Events for component eta2: " << iEta2 << " +/- " << iEta2Err ;
-
-// std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamPhi=_fitParam.Phis[paramEnumJpsiGamPipPimPi0KK::PsiToEtacGamma];
-// std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess > PsiToEtacGamMag=_fitParam.Mags[paramEnumJpsiGamPipPimPi0KK::PsiToEtacGamma];
-// std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator phiIter = PsiToEtacGamPhi.begin();
-// for(; phiIter!=PsiToEtacGamPhi.end(); phiIter++){
-// double thePhase = phiIter->second;
-// double theMag = PsiToEtacGamMag[phiIter->first];
-
-// etaImg = theMag*cos(thePhase) ;
-// etaReal = theMag*sin(thePhase);
-// }
-
-// _massIndepTuple->Fill(meanMassRange, iEta, iEtaErr, iF0, iF0Err, iF1, iF1Err, iF2, iF2Err, iEta2, iEta2Err, etaImg, etaReal );
-
-// }
-
-JpsiToOmegaPhiGamHist::~JpsiToOmegaPhiGamHist()
-{
- _theTFile->Write();
- _theTFile->Close();
-}
-
-void JpsiToOmegaPhiGamHist::initRootStuff()
-{
- std::string rootFileName="./JpsiGamPipPimPi0KK.root";
-
- _theTFile=new TFile(rootFileName.c_str(),"recreate");
- int xbins = 50;
- double xmin=0.8;
- double xmax =10.;
- int ybins=xbins;
- double ymin=xmin;
- double ymax=xmax;
-
- _dalitzDataHist= new TH2F("_dalitzDataHist","Dpl K+K- K+#pi^{0} data",xbins, xmin, xmax, ybins, ymin, ymax );
- _dalitzMcHist= new TH2F("_dalitzMcHist","Dpl K+K- K+#pi^{0} MC",xbins, xmin, xmax, ybins, ymin, ymax);
- _dalitzFittedHist= new TH2F("_dalitzFittedHist","Dpl K+K- K+#pi^{0} fit",xbins, xmin, xmax, ybins, ymin, ymax );
-
-
- int nbins=400;
- xmin=2;
- xmax=4.;
- _PhiPhiMassDataHist = new TH1F("_PhiPhiMassDataHist", "_PhiPhiMassDataHist", nbins, xmin, xmax);
- _PhiPhiMassMcHist = new TH1F("_PhiPhiMassMcHist", "_PhiPhiMassMcHist", nbins, xmin, xmax);
- _PhiPhiMassFittedHist = new TH1F("_PhiPhiMassFittedHist", "_PhiPhiMassFittedHist", nbins, xmin, xmax);
-
- nbins=50;
- xmin=0.8;
- xmax=1.2;
- _PipPimPi0MassDataHist = new TH1F("_PipPimPi0MassDataHist", "_PipPimPi0MassDataHist", nbins, xmin, xmax);
- _PipPimPi0MassMcHist = new TH1F("_PipPimPi0MassMcHist", "_PipPimPi0MassMcHist", nbins, xmin, xmax);
- _PipPimPi0MassFittedHist = new TH1F("_PipPimPi0MassFittedHist", "_PipPimPi0MassFittedHist", nbins, xmin, xmax);
-
- _KpKmMassDataHist = new TH1F("_KpKmMassDataHist", "_KpKmMassDataHist", nbins, xmin, xmax);
- _KpKmMassMcHist = new TH1F("_KpKmMassMcHist", "_KpKmMassMcHist", nbins, xmin, xmax);
- _KpKmMassFittedHist = new TH1F("_KpKmMassFittedHist", "_KpKmMassFittedHist", nbins, xmin, xmax);
-
- _costPip_PipPimPi0HeliDataHist= new TH1F("_costPip_PipPimPi0HeliDataHist", "cos(#Theta_{Pip}) PipPimPi0Heli data", 100, -1., 1.);
- _costPip_PipPimPi0HeliMcHist= new TH1F("_costPip_PipPimPi0HeliMcHist", "cos(#Theta_{Pip}) PipPimPi0Heli Mc", 100, -1., 1.);
- _costPip_PipPimPi0HeliFittedHist= new TH1F("_costPip_PipPimPi0HeliFittedHist", "cos(#Theta_{Pip}) PipPimPi0Heli Fit", 100, -1, 1);
- _phiPip_PipPimPi0HeliDataHist= new TH1F("_phiPip_PipPimPi0HeliDataHist", "#Phi_{Pip} PipPimPi0Heli data", 100, -TMath::Pi(), TMath::Pi());
- _phiPip_PipPimPi0HeliMcHist= new TH1F("_phiPip_PipPimPi0HeliMcHist", "#Phi_{Pip} PipPimPi0Heli Mc", 100, -TMath::Pi(), TMath::Pi());
- _phiPip_PipPimPi0HeliFittedHist= new TH1F("_phiPip_PipPimPi0HeliFittedHist", "#Phi_{Pip} PipPimPi0Heli Fit", 100, -TMath::Pi(), TMath::Pi());
-
- _costKp_KpKmHeliDataHist= new TH1F("_costKp_KpKmHeliDataHist", "cos(#Theta_{K+}) PipPimPi0Heli data", 100, -1., 1.);
- _costKp_KpKmHeliMcHist= new TH1F("_costKp_KpKmHeliMcHist", "cos(#Theta_{K+}) PipPimPi0Heli Mc", 100, -1., 1.);
- _costKp_KpKmHeliFittedHist= new TH1F("_costKp_KpKmHeliFittedHist", "cos(#Theta_{K+}) PipPimPi0Heli Fit", 100, -1, 1);
- _phiKp_KpKmHeliDataHist= new TH1F("_phiKp_KpKmHeliDataHist", "#Phi_{K+} PipPimPi0Heli data", 100, -TMath::Pi(), TMath::Pi());
- _phiKp_KpKmHeliMcHist= new TH1F("_phiKp_KpKmHeliMcHist", "#Phi_{K+} PipPimPi0Heli Mc", 100, -TMath::Pi(), TMath::Pi());
- _phiKp_KpKmHeliFittedHist= new TH1F("_phiKp_KpKmHeliFittedHist", "#Phi_{K+} PipPimPi0Heli Fit", 100, -TMath::Pi(), TMath::Pi());
-
-
- _costGamCmDataHist= new TH1F("_costGamCmDataHist", "cos(#Theta_{#gamma}) CM data", 100, -1., 1.);
- _costGamCmMcHist= new TH1F("_costGamCmMcHist", "cos(#Theta_{#gamma}) CM Mc", 100, -1., 1.);
- _costGamCmFittedHist= new TH1F("_costGamCmFittedHist", "cos(#Theta_{#gamma}) CM Fit", 100, -1, 1);
-
- _costPhi_KpKmDataHist= new TH1F("_costPhi_KpKmDataHist", "cos(#Theta_{#phi}) K+ K- data", 100, -1., 1.);
- _costPhi_KpKmMcHist= new TH1F("_costPhi_KpKmMcHist", "cos(#Theta_{#phi}) K+ K- Mc", 100, -1., 1.);
- _costPhi_KpKmFittedHist= new TH1F("_costPhi_KpKmFittedHist", "cos(#Theta_{#phi}) K+ K- Fit", 100, -1., 1.);
-
- _phiPhi_KpKmDataHist= new TH1F("_phiPhi_KpKmDataHist", "cos(#Phi_{#phi}) K+ K- data", 100, -TMath::Pi(), TMath::Pi());
- _phiPhi_KpKmMcHist= new TH1F("_phiPhi_KpKmMcHist", "cos(#Phi_{#phi}) K+ K- Mc", 100, -TMath::Pi(), TMath::Pi());
- _phiPhi_KpKmFittedHist= new TH1F("_phiPhi_KpKmFittedHist", "cos(#Phi_{#phi}) K+ K- Fit", 100, -TMath::Pi(), TMath::Pi());
-
- _chiDataHist= new TH1F("_chiDataHist", "#chi data", 90, 0., 90.);
- _chiMcHist= new TH1F("_chiMcHist", "#chi Mc", 90, 0., 90.);
- _chiFittedHist= new TH1F("_chiFittedHist", "#chi Fit", 90, 0., 90.);
-
- std::string tupleVariables = "mPipPimPi0KpKm:mPipPimPi0:mKpKm:PipPimPi0KpKmCostTheta:gamCosTheta:PipPimPi0CosTheta:KpKmCosTheta:PipCosTheta:KpCosTheta:decPlaneChi:testHeli:weight";
-
-
- _dataTuple=new TNtuple("_dataTuple", "data ntuple", tupleVariables.data());
- _mcTuple=new TNtuple("_mcTuple", "mc ntuple", tupleVariables.data());
-
- _massIndepTuple = new TNtuple("_massIndepTuple","results from mass indep. fits", ("mass:eta:etaErr:f0:f0Err:f1:f1Err:f2:f2Err:eta2:eta2Err:etaReal:etaImg") );
-}
-
-void JpsiToOmegaPhiGamHist::plotDalitz(TH2F* theHisto, EvtDataNew* theData, double weight)
-{
- Vector4<double>& V4_PipPimPi0_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPsi] ;
- Vector4<double>& V4_KpKm_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPsi] ;
- Vector4<double>& V4_gamma_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_gamma_HeliPsi] ;
-
- double gphi1 = (V4_gamma_HeliPsi + V4_PipPimPi0_HeliPsi ).M();
- double gphi2 = (V4_gamma_HeliPsi + V4_KpKm_HeliPsi ).M();
- theHisto->Fill( gphi1*gphi1, gphi2*gphi2 ,weight);
-}
-
-void JpsiToOmegaPhiGamHist::plotPhiPhi(TH1F* theHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0KpKm_HeliPsi];
- theHisto->Fill( v4.M(), weight );
-}
-
-void JpsiToOmegaPhiGamHist::plotPipPimPi0(TH1F* theHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPsi];
- theHisto->Fill( v4.M(), weight );
-}
-void JpsiToOmegaPhiGamHist::plotKpKm(TH1F* theHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm];
- theHisto->Fill( v4.M(), weight );
-}
-
-void JpsiToOmegaPhiGamHist::plotCostPhiPip(TH1F* theCostHisto, TH1F* thePhiHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPipPimPi0];
- theCostHisto->Fill( v4.CosTheta(), weight );
- thePhiHisto->Fill( v4.Phi(), weight );
-}
-
-void JpsiToOmegaPhiGamHist::plotCostPhiKp(TH1F* theCostHisto, TH1F* thePhiHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliKpKm];
- theCostHisto->Fill( v4.CosTheta(), weight );
- thePhiHisto->Fill( v4.Phi(), weight );
-}
-
-void JpsiToOmegaPhiGamHist::plotCostGam(TH1F* theCostHisto, EvtDataNew* theData, double weight){
- Vector4<double>& v4 = theData->FourVecsDec[enumJpsiGamX4V::V4_gamma_HeliPsi];
- theCostHisto->Fill( v4.CosTheta(), weight );
-}
-
-void JpsiToOmegaPhiGamHist::plotCostPhi_PhiPhiHeli(TH1F* theCostHisto, TH1F* thePhiHisto, const Vector4<double>& the4Vec, double weight){
- theCostHisto->Fill( the4Vec.CosTheta(), weight);
- thePhiHisto->Fill( the4Vec.Phi(), weight);
-}
-
-void JpsiToOmegaPhiGamHist::plotChi(TH1F* theChiHisto, EvtDataNew* theData, double weight){
-// Vector4<double>& V4_PipPimPi0KpKm_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0KpKm_HeliPsi] ;
-// Vector4<double>& V4_Pip_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPsi] ;
-// Vector4<double>& V4_Kl_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kl_HeliPsi] ;
-// Vector4<double>& V4_Kp_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliPsi] ;
-// Vector4<double>& V4_Km_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Km_HeliPsi] ;
-
-// double thePhiPhiDecayPlaneAngle = decayAngleChi( V4_PipPimPi0KpKm_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi, V4_Pip_HeliPsi, V4_Kl_HeliPsi );
-
-
- Vector4<double>& V4_normKpKmDecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normKpKmDecHeliPipPimPi0KpKm];
- Vector4<double>& V4_normPipPimPi0DecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normPipPimPi0DecHeliPipPimPi0KpKm];
-
- double cosChi=(V4_normKpKmDecHeliPipPimPi0KpKm.Px()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Px()
- +V4_normKpKmDecHeliPipPimPi0KpKm.Py()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Py()
- +V4_normKpKmDecHeliPipPimPi0KpKm.Pz()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Pz())
- / (V4_normKpKmDecHeliPipPimPi0KpKm.P()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.P());
-
- double chi=acos(fabs(cosChi));
- theChiHisto->Fill(chi*180./TMath::Pi(),weight);
-}
-
-void JpsiToOmegaPhiGamHist::fillTuple( TNtuple* theTuple, EvtDataNew* theData, double weight)
-{
-
- Vector4<double>& V4_PipPimPi0KpKm_HeliPsi = theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0KpKm_HeliPsi] ;
- Vector4<double>& V4_PipPimPi0_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPsi] ;
- Vector4<double>& V4_KpKm_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPsi] ;
- Vector4<double>& V4_gamma_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_gamma_HeliPsi] ;
-
- Vector4<double>& V4_PipPimPi0_HeliPipPimPi0KpKm= theData->FourVecsDec[enumJpsiGamX4V::V4_PipPimPi0_HeliPipPimPi0KpKm] ;
- Vector4<double>& V4_KpKm_HeliPipPimPi0KpKm= theData->FourVecsDec[enumJpsiGamX4V::V4_KpKm_HeliPipPimPi0KpKm] ;
- Vector4<double>& V4_Pip_HeliPipPimPi0= theData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPipPimPi0] ;
- Vector4<double>& V4_Kp_HeliKpKm= theData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliKpKm] ;
-
- Vector4<double>& V4_Pip_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Pip_HeliPsi] ;
- Vector4<double>& V4_Kl_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kl_HeliPsi] ;
- Vector4<double>& V4_Kp_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Kp_HeliPsi] ;
- Vector4<double>& V4_Km_HeliPsi= theData->FourVecsDec[enumJpsiGamX4V::V4_Km_HeliPsi] ;
-
-
- Vector4<double>& V4_normKpKmDecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normKpKmDecHeliPipPimPi0KpKm];
- Vector4<double>& V4_normPipPimPi0DecHeliPipPimPi0KpKm = theData->FourVecsDec[enumJpsiGamX4V::V4_normPipPimPi0DecHeliPipPimPi0KpKm];
- double cosChi=(V4_normKpKmDecHeliPipPimPi0KpKm.Px()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Px()
- +V4_normKpKmDecHeliPipPimPi0KpKm.Py()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Py()
- +V4_normKpKmDecHeliPipPimPi0KpKm.Pz()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.Pz())
- / (V4_normKpKmDecHeliPipPimPi0KpKm.P()*V4_normPipPimPi0DecHeliPipPimPi0KpKm.P());
-
- double chi=acos(fabs(cosChi));
- double thePhiPhiDecayPlaneAngle = chi/TMath::Pi()*180.;
-
- //double thePhiPhiDecayPlaneAngle = decayAngleChi( V4_PipPimPi0KpKm_HeliPsi, V4_Kp_HeliPsi, V4_Km_HeliPsi, V4_Pip_HeliPsi, V4_Kl_HeliPsi );
- double testHeli = costDecHeli( V4_PipPimPi0KpKm_HeliPsi+V4_gamma_HeliPsi, V4_Pip_HeliPsi+V4_Kl_HeliPsi+V4_Km_HeliPsi+V4_Kp_HeliPsi, V4_Pip_HeliPsi+V4_Kl_HeliPsi );
-
- theTuple->Fill( V4_PipPimPi0KpKm_HeliPsi.M(),
- V4_PipPimPi0_HeliPsi.M(),
- V4_KpKm_HeliPsi.M(),
- V4_PipPimPi0KpKm_HeliPsi.CosTheta(),
- V4_gamma_HeliPsi.CosTheta(),
- V4_PipPimPi0_HeliPipPimPi0KpKm.CosTheta(),
- V4_KpKm_HeliPipPimPi0KpKm.CosTheta(),
- V4_Pip_HeliPipPimPi0.CosTheta(),
- V4_Kp_HeliKpKm.CosTheta(),
- thePhiPhiDecayPlaneAngle,
- testHeli,
- weight);
-}
-
-
-double JpsiToOmegaPhiGamHist::decayAngleChi(const Vector4<double>& v4_p,const Vector4<double>& v4_d1,
- const Vector4<double>& v4_d2,const Vector4<double>& v4_h1,
- const Vector4<double>& v4_h2 ) {
-
- TLorentzVector p4_p( v4_p.Px(), v4_p.Py(), v4_p.Pz(), v4_p.E() );
- TLorentzVector p4_d1p( v4_d1.Px(), v4_d1.Py(), v4_d1.Pz(), v4_d1.E() );
- TLorentzVector p4_d2p( v4_d2.Px(), v4_d2.Py(), v4_d2.Pz(), v4_d2.E() );
- TLorentzVector p4_h1p( v4_h1.Px(), v4_h1.Py(), v4_h1.Pz(), v4_h1.E() );
- TLorentzVector p4_h2p( v4_h2.Px(), v4_h2.Py(), v4_h2.Pz(), v4_h2.E() );
-
-
- // boost all vectors parent restframe
-
- p4_d1p.Boost( -p4_p.BoostVector() );
- p4_d2p.Boost( -p4_p.BoostVector() );
- p4_h1p.Boost( -p4_p.BoostVector() );
- p4_h2p.Boost( -p4_p.BoostVector() );
-
-
- TVector3 d1_perp,d1_prime,h1_perp;
- TVector3 D;
-
- D=(p4_d1p+p4_d2p).Vect();
-
- d1_perp=p4_d1p.Vect()-(D.Dot(p4_d1p.Vect())/D.Dot(D))*D;
- h1_perp=p4_h1p.Vect()-(D.Dot(p4_h1p.Vect())/D.Dot(D))*D;
-
- // orthogonal to both D and d1_perp
-
- d1_prime=D.Cross(d1_perp);
-
- d1_perp= d1_perp* (1./d1_perp.Mag());
- d1_prime= d1_prime * (1./d1_prime.Mag());
-
- double x,y;
-
- x=d1_perp.Dot(h1_perp);
- y=d1_prime.Dot(h1_perp);
-
- double chi=atan2(y,x);
-
- if (chi<0.0) chi+=2.*TMath::Pi();
-
- return chi;
-
-}
-
-
-
diff --git a/Examples/JpsiToOmegaPhiGam/MJpsiToOmegaPhiGamApp.cc~ b/Examples/JpsiToOmegaPhiGam/MJpsiToOmegaPhiGamApp.cc~
deleted file mode 100644
index cdbd4cf14012d55269d2c535457c6830399ae1e5..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/MJpsiToOmegaPhiGamApp.cc~
+++ /dev/null
@@ -1,353 +0,0 @@
-#include <iostream>
-#include <cstring>
-#include <string>
-#include <sstream>
-#include <vector>
-#include <map>
-#include <iterator>
-
-#ifdef _OPENMP
-#include <omp.h>
-#endif
-
-#include <boost/shared_ptr.hpp>
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh"
-#include "PwaUtils/PsiToXGamParser.hh"
-#include "PwaUtils/PsiToXGamReader.hh"
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamHist.hh"
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamStates.hh"
-//#include "Examples/JpsiGamEtaPiPiNew/JpsiGamEtaPiPiEventListNew.hh"
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamLh.hh"
-
-#include "PwaUtils/EvtDataBaseListNew.hh"
-#include "PwaUtils/FitParamsBaseNew.hh"
-#include "PwaUtils/StreamFitParmsBaseNew.hh"
-#include "PwaUtils/PwaFcnBaseNew.hh"
-#include "PwaUtils/AbsLhNew.hh"
-#include "Utils/ErrLogUtils.hh"
-
-#include "Setup/PwaEnv.hh"
-#include "Particle/ParticleTable.hh"
-#include "Particle/Particle.hh"
-#include "Event/EventList.hh"
-#include "Event/Event.hh"
-#include "Particle/PdtParser.hh"
-#include "qft++/topincludes/tensor.hh"
-
-#include "ErrLogger/ErrLogger.hh"
-#include "Minuit2/MnUserParameters.h"
-#include "Minuit2/MnMigrad.h"
-#include "Minuit2/FunctionMinimum.h"
-#include "Minuit2/MnMinos.h"
-#include "Minuit2/MnStrategy.h"
-#include "Minuit2/MnPrint.h"
-#include "Minuit2/MnScan.h"
-
-using namespace ROOT::Minuit2;
-
-
-
-int main(int __argc,char *__argv[]){
- clock_t start, end;
- start= clock();
-
- // Parse the command line
- static PsiToXGamParser theAppParams(__argc, __argv);
-
- // Set the desired error logging mode
- setErrLogMode(theAppParams.getErrLogMode());
-
-#ifdef _OPENMP
- const int noOfThreads=theAppParams.noOfThreads();
- omp_set_num_threads(noOfThreads);
-#endif
-
-
- std::string theCfgFile = theAppParams.getConfigFile();
-
- std::string jobOption = theAppParams.getjobOption();
-
- const std::string datFile=theAppParams.dataFile();
- const std::string mcFile=theAppParams.mcFile();
- Info << "data file: " << datFile ; // << endmsg;
- Info << "mc file: " << mcFile ; // << endmsg;
-
- ParticleTable pTable;
- PdtParser parser;
- std::string theSourcePath=getenv("CMAKE_SOURCE_DIR");
- std::string pdtFile(theSourcePath+"/Particle/pdt.table");
- if (!parser.parse(pdtFile, pTable)) {
- 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);
-
- bool withEvtWeight=theAppParams.useEvtWeight();
- Info << "EvtWeight: " << withEvtWeight << endmsg;
- PsiToXGamReader eventReader(fileNames, 5, 0, withEvtWeight);
- EventList eventsData;
- eventReader.fillMassRange(eventsData, massRange );
-
- if (!eventsData.findParticleTypes(pTable))
- Warning << "could not find all particles" ; // << endmsg;
-
- Info << "\nFile has " << eventsData.size() << " events. Each event has "
- << eventsData.nextEvent()->size() << " final state particles.\n" ; // << endmsg;
- eventsData.rewind();
-
- Event* anEvent;
- int evtCount = 0;
- while ((anEvent = eventsData.nextEvent()) != 0 && evtCount < 10) {
- Info << "\n"
- << *(anEvent->p4(0)) << "\tm = " << anEvent->p4(0)->Mass() << "\n"
- << *(anEvent->p4(1)) << "\tm = " << anEvent->p4(1)->Mass() << "\n"
- << *(anEvent->p4(2)) << "\tm = " << anEvent->p4(2)->Mass() << "\n"
- << *(anEvent->p4(3)) << "\tm = " << anEvent->p4(3)->Mass() << "\n"
- << *(anEvent->p4(4)) << "\tm = " << anEvent->p4(4)->Mass() << "\n"
- ; // << endmsg;
- ++evtCount;
- }
- eventsData.rewind();
-
-
- std::vector<std::string> fileNamesMc;
- fileNamesMc.push_back(mcFile);
- PsiToXGamReader eventReaderMc(fileNamesMc, 5, 0, withEvtWeight);
- EventList eventsMc;
- eventReaderMc.fillMassRange(eventsMc, massRange);
- eventsMc.rewind();
-
- //
- //calculate helicity angles, fill map with D-functions
- //
-
- boost::shared_ptr<JpsiToOmegaPhiGamEventList> theJpsiGamXEventListPtr(new JpsiToOmegaPhiGamEventList());
- theJpsiGamXEventListPtr->ratioMcToData(theAppParams.ratioMcToData());
- theJpsiGamXEventListPtr->read(eventsData, eventsMc);
-
- std::string mode=theAppParams.mode();
- std::cout << "Mode: " << mode << std::endl;
- if (mode=="plotmode"){
- // JpsiToOmegaPhiGamHist theHist(theJpsiGamXEventListPtr,theAppParams.massRange());
- JpsiToOmegaPhiGamHist theHist(theJpsiGamXEventListPtr);
- // theHist.setMassRange(theAppParams.massRange() );
- return 0;
- }
-
- //
- //retrieve hypotheses
- //
-
- boost::shared_ptr<JpsiToOmegaPhiGamStates> jpsiToOmegaPhiGamStatesPtr(new JpsiToOmegaPhiGamStates());
- jpsiToOmegaPhiGamStatesPtr->print(std::cout);
-
- const std::vector<std::string> hypVec=theAppParams.enabledHyps();
- boost::shared_ptr<AbsLhNew> theLhPtr;
-
- std::string startWithHyp=theAppParams.startHypo();
-
- if (startWithHyp=="production"){
- theLhPtr = boost::shared_ptr<AbsLhNew> (new JpsiToOmegaPhiGamLh(theJpsiGamXEventListPtr, hypVec, jpsiToOmegaPhiGamStatesPtr));
- }
- else {
- Alert << "start with hypothesis " << startWithHyp << " not supported!!!!" ; // << endmsg;
- exit(1);
- }
-
- boost::shared_ptr<FitParamsBaseNew> theFitParamBase=boost::shared_ptr<FitParamsBaseNew>(new FitParamsBaseNew());
-
- if (mode=="dumpDefaultParams"){
- fitParamsNew defaultVal;
- fitParamsNew defaultErr;
- theLhPtr->getDefaultParams(defaultVal, defaultErr);
- std::string defaultparamsname = "defaultparams" + jobOption + ".dat";
- std::ofstream theStreamDefault ( defaultparamsname.c_str() );
-
- theFitParamBase->dumpParams(theStreamDefault, defaultVal, defaultErr);
- return 0;
- }
-
- std::string paramStreamerPath=theAppParams.fitParamFile();
- StreamFitParmsBaseNew theParamStreamer(paramStreamerPath, theLhPtr);
- fitParamsNew theStartparams=theParamStreamer.getFitParamVal();
- fitParamsNew theErrorparams=theParamStreamer.getFitParamErr();
-
- PwaFcnBaseNew theFcn(theLhPtr, theFitParamBase);
- MnUserParameters upar;
- theFitParamBase->setMnUsrParams(upar, theStartparams, theErrorparams);
-
- std::cout << "\n\n**************** Minuit Fit parameter **************************" << std::endl;
- for (int i=0; i<int(upar.Params().size()); ++i){
- std::cout << upar.Name(i) << "\t" << upar.Value(i) << "\t" << upar.Error(i) << std::endl;
- }
-
- const std::vector<std::string> fixedParams=theAppParams.fixedParams();
-
- const unsigned int noOfFreeFitParams = upar.Params().size()-fixedParams.size();
- if (mode=="qaMode"){
-
- Info << "\nThe parameter values are: " << "\n" << endmsg;
- theFitParamBase->printParams(theStartparams);
-
- Info << "\nThe parameter errors are: " << "\n" << endmsg;
- theFitParamBase->printParams(theErrorparams);
-
- double theLh=theLhPtr->calcLogLh(theStartparams);
- Info <<"theLh = "<< theLh << endmsg;
-
-
- JpsiToOmegaPhiGamHist theHist(theLhPtr, theStartparams);
- theHist.setMassRange(theAppParams.massRange() );
-
- double evtWeightSumData = theJpsiGamXEventListPtr->NoOfWeightedDataEvts();
- double BICcriterion=2.*theLh+noOfFreeFitParams*log(evtWeightSumData);
- double AICcriterion=2.*theLh+2.*noOfFreeFitParams;
- double AICccriterion=AICcriterion+2.*noOfFreeFitParams*(noOfFreeFitParams+1)/(evtWeightSumData-noOfFreeFitParams-1);
-
- Info << "noOfFreeFitParams:\t" <<noOfFreeFitParams;
- Info << "evtWeightSumData:\t" <<evtWeightSumData;
- Info << "BIC:\t" << BICcriterion << endmsg;
- Info << "AIC:\t" << AICcriterion << endmsg;
- Info << "AICc:\t" << AICccriterion << endmsg;
-
- std::string qaSummaryFileName = "qaSummary" + jobOption + ".dat";
- std::ofstream theQaStream ( qaSummaryFileName.c_str() );
- theQaStream << "BIC\t" << BICcriterion << "\n";
- theQaStream << "AICa\t" << AICcriterion << "\n";
- theQaStream << "AICc\t" << AICccriterion << "\n";
- theQaStream << "logLh\t" << theLh << "\n";
- theQaStream << "free parameter\t" << noOfFreeFitParams << "\n";
- theQaStream.close();
-
- end= clock();
- double cpuTime= (end-start)/ (CLOCKS_PER_SEC);
- Info << "cpuTime:\t" << cpuTime << "\tsec" << endmsg;
- return 0;
- }
-
- if (mode=="pwa"){
- bool cacheAmps = theAppParams.cacheAmps();
- Info << "caching amplitudes enabled / disabled:\t" << cacheAmps << endmsg;
- if (cacheAmps) theLhPtr->cacheAmplitudes();
- std::vector<std::string>::const_iterator itFix;
- for (itFix=fixedParams.begin(); itFix!=fixedParams.end(); ++itFix){
- upar.Fix( (*itFix) );
- }
-
- // // bool prescan=false;
- // // if(prescan){
- // // upar.Fix(0);
- // // MnScan theScan(theFcn, upar);
- // // FunctionMinimum smin = theScan();
- // // MnUserParameterState sState = smin.UserState();
- // // cout << "After scan" << endl;
- // // cout << sState << endl;
-
- // // upar = smin.UserParameters();
- // // upar.Release(0);
- // // }
-
- MnMigrad migrad(theFcn, upar);
- Info <<"start migrad "<< endmsg;
- FunctionMinimum min = migrad();
-
- if(!min.IsValid()) {
- //try with higher strategy
- Info <<"FM is invalid, try with strategy = 2."<< endmsg;
- MnMigrad migrad2(theFcn, min.UserState(), MnStrategy(2));
- min = migrad2();
- }
-
- MnUserParameters finalUsrParameters=min.UserParameters();
- const std::vector<double> finalParamVec=finalUsrParameters.Params();
- fitParamsNew finalFitParams=theStartparams;
- theFitParamBase->getFitParamVal(finalParamVec, finalFitParams);
-
- // //MnUserCovariance theCov = min.UserCovariance() ;
- // //cout << "User vov : "<< endl;
- // //cout << theCov << endl;
-
- theFitParamBase->printParams(finalFitParams);
- double theLh=theLhPtr->calcLogLh(finalFitParams);
- Info <<"theLh = "<< theLh << endmsg;
-
-
- const std::vector<double> finalParamErrorVec=finalUsrParameters.Errors();
- fitParamsNew finalFitErrs=theErrorparams;
- theFitParamBase->getFitParamVal(finalParamErrorVec, finalFitErrs);
-
- std::string finalResultname = "finalResult" + jobOption + ".dat";
- std::ofstream theStream ( finalResultname.c_str() );
- //std::ofstream theStream ( "finalResult.dat");
- theFitParamBase->dumpParams(theStream, finalFitParams, finalFitErrs);
-
- MnUserCovariance theCovMatrix = min.UserCovariance();
- std::cout << min << std::endl;
-
- // // //std::ofstream theCompStream ( "componentIntensity.dat");
- // // //theProdLh->dumpComponentIntensity( theCompStream, finalFitParams, theErrMatrix );
- // JpsiGamEtaPiPiHistNew theHist(theLhPtr, finalFitParams,theAppParams.massRange());
- // //theHist.setMassRange(theAppParams.massRange() );
- // // // theHist.fill();
- // end= clock();
- // double cpuTime= (end-start)/ (CLOCKS_PER_SEC);
- // Info << "cpuTime:\t" << cpuTime << "\tsec" << endmsg;
-
-
-
-
- // std::cout << "Start event number calculation for each wave" << std::endl;
-
- // std::vector<std::string> hypVec_test=theAppParams.enabledHyps();
- // fitParamsNew finalFitParams_test=finalFitParams;
- // std::vector<std::string>::iterator it;
- // int hypnumber=0;
-
- // for (it=hypVec_test.begin(); it!=hypVec_test.end();++it){
- // //std::cout << "hypothesis\t" << (*it) << "\t enabled\n";
- // if ((*it).find("Gamma")==0) hypnumber++;
- // }
- // std::cout << std::endl;
- // std::cout << "Number of hypothesis found: " << hypnumber << std::endl;
-
- // int j;
- // for (int i=1;i<=hypnumber;i++){
- // j=1;
- // for (it=hypVec_test.begin(); it!=hypVec_test.end();++it){
- // // Mark bad hypothesis
- // if ((*it).find("Gamma")==0 && i!=j) {
- // (*It).Insert(0, "#");
- // j++;
- // }
- // else{ if ((*it).find("Gamma")==0) j++; }
- // }
- // std::cout << "Start calulation with following hypothesis:" << std::endl;
- // if (startWithHyp=="production"){
- // theLhPtr = boost::shared_ptr<AbsLhNew> (new JpsiGamEtaPiPiProdLhNew(theJpsiGamEtaPiPiEventListPtr, hypVec_test, jpsiGamEtaPiPiStatesPtr));
- // }
- // else {
- // Alert << "start with hypothesis " << startWithHyp << " not supported!!!!" ; // << endmsg;
- // exit(1);
- // }
- // double theLh=theLhPtr->calcLogLh(finalFitParams_test);
- // Info <<"theLh = "<< theLh << endmsg;
- // JpsiGamEtaPiPiHistNew theHist(theLhPtr, finalFitParams_test ,theAppParams.massRange());
- // // Unmark bad hypothesis
- // for (it=hypVec_test.begin(); it!=hypVec_test.end();++it){
- // if ((*it).find("#")==0) (*it).erase(0,1);
- // }
- // std::cout << std::endl;
- // }
-
- return 0;
- }
-
-
- return 0;
-}
-
diff --git a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc b/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc
index 4e27c8a004344a17922a9bf806d1adf13c6e4b50..08776146b4ce9f3c7b119320e9db4b677feaec4c 100644
--- a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc
+++ b/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc
@@ -22,7 +22,7 @@ XToOmegaPhiDecAmps::XToOmegaPhiDecAmps(const std::string& name, const std::vecto
,_flatteMassFit(false)
,_phiMass( 1.019455)
,_omegaMass( 0.78265)
- ,_phiPhiPair(make_pair(_phiMass, _phiMass))
+ ,_phiPhiPair(make_pair(_omegaMass, _omegaMass))
,_omegaPhiPair(make_pair(_omegaMass, _phiMass))
,_theStatesPtr(theStates)
{
@@ -41,7 +41,7 @@ complex<double> XToOmegaPhiDecAmps::XdecAmp(Spin lamX, EvtDataNew* theData){
complex<double> result(0.,0.);
- result+=XToPhiPhiAmp(lamX, theData);
+ result+=XToPhiOmegaAmp(lamX, theData);
complex<double> dynModel(1.,0.);
if (!_massIndependent){
@@ -70,7 +70,7 @@ complex<double> XToOmegaPhiDecAmps::XdecAmp(Spin lamX, EvtDataNew* theData){
return result;
}
-complex<double> XToOmegaPhiDecAmps::XToPhiPhiAmp(Spin lamX, EvtDataNew* theData){
+complex<double> XToOmegaPhiDecAmps::XToPhiOmegaAmp(Spin lamX, EvtDataNew* theData){
complex<double> result(0.,0.);
std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itXMag;
@@ -82,17 +82,17 @@ complex<double> XToOmegaPhiDecAmps::XToPhiPhiAmp(Spin lamX, EvtDataNew* theData)
double theXPhi=_currentParamPhis[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;
+ for(Spin lambdaOmega=-1; lambdaOmega<=1; lambdaOmega++){
+ for(Spin lambdaPhi=-1; lambdaPhi<=1; lambdaPhi++){
+ Spin lambda = lambdaOmega-lambdaPhi;
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 )
+ *Clebsch(1, lambdaOmega, 1, -lambdaPhi, XState->S, lambda )
*conj( theData->WignerDsDec[enumJpsiGamXDfunc::Df_XToOmegaPhi_KK][XState->J][lamX][lambda]);
- amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 );
+ amp = amp * phiOmegaTo2K3piAmp( theData, lambdaOmega, lambdaPhi );
result +=amp;
}
@@ -102,13 +102,13 @@ complex<double> XToOmegaPhiDecAmps::XToPhiPhiAmp(Spin lamX, EvtDataNew* theData)
return result;
}
-complex<double> XToOmegaPhiDecAmps::phiphiTo4KAmp( EvtDataNew* theData, Spin lambdaPhi1, Spin lambdaPhi2 ){
+complex<double> XToOmegaPhiDecAmps::phiOmegaTo2K3piAmp( EvtDataNew* theData, Spin lambdaOmega, Spin lambdaPhi ){
complex<double> result(0.,0.);
double lambda = theData->KinematicVariables[enumJpsiGamXKin::OmegaDecLambda];
- result = 3. * conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_OmegaToPipPimPi0][1][lambdaPhi1][0]) * sqrt(lambda)
- * 3.* conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_PhiToKpKm][1][lambdaPhi2][0]);
+ result = 3. * conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_OmegaToPipPimPi0][1][lambdaOmega][0]) * sqrt(lambda)
+ * 3.* conj(theData->WignerDsDec[enumJpsiGamXDfunc::Df_PhiToKpKm][1][lambdaPhi][0]);
return result;
}
diff --git a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc~ b/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc~
deleted file mode 100644
index ba065bea5f07424e4f616f156b87f465d7e924d5..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.cc~
+++ /dev/null
@@ -1,265 +0,0 @@
-#include <getopt.h>
-#include <fstream>
-#include <string>
-
-#include "Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.hh"
-#include "qft++/relativistic-quantum-mechanics/Utils.hh"
-#include "ErrLogger/ErrLogger.hh"
-// #include "PwaUtils/EvtDataBaseListNew.hh"
-#include "Examples/JpsiToOmegaPhiGam/JpsiToOmegaPhiGamEventList.hh"
-
-#ifdef _OPENMP
-#include <omp.h>
-#endif
-
-XToOmegaPhiDecAmps::XToOmegaPhiDecAmps(const std::string& name, const std::vector<std::string>& hypVec, boost::shared_ptr<JpsiToOmegaPhiGamStates> theStates, Spin spinX, int parity) :
- AbsXdecAmp(name, hypVec, spinX, parity)
- ,_phiPhiKey(name+"ToPhiPhi")
- ,_xBWKey(name+"BreitWigner")
- ,_xFlatteKey(name+"Flatte")
- ,_massIndependent(true)
- ,_bwMassFit(false)
- ,_flatteMassFit(false)
- ,_phiMass( 1.019455)
- ,_omegaMass( 0.78265)
- ,_phiPhiPair(make_pair(_phiMass, _phiMass))
- ,_omegaPhiPair(make_pair(_omegaMass, _phiMass))
- ,_theStatesPtr(theStates)
-{
- initialize();
-}
-
-XToOmegaPhiDecAmps::~XToOmegaPhiDecAmps()
-{
-}
-
-complex<double> XToOmegaPhiDecAmps::XdecAmp(Spin lamX, EvtDataNew* theData){
- int evtNo=theData->evtNo;
-
- if ( _cacheAmps && !_recalculate) return _cachedAmpMap[evtNo][lamX];
-
-
- complex<double> result(0.,0.);
-
- result+=XToPhiPhiAmp(lamX, theData);
-
- complex<double> dynModel(1.,0.);
- if (!_massIndependent){
- Vector4<double> p4PhiPhi = theData->FourVecsDec[enumJpsiGamX4V::V4_KsKlKpKm_HeliPsi];
- if(_bwMassFit){
- dynModel=BreitWigner(p4PhiPhi, _currentXMass, _currentXWidth);
- }
- else if(_flatteMassFit){
- dynModel=Flatte( p4PhiPhi, _phiPhiPair, _omegaPhiPair, _currentXMass, _currentgFactorPhiPhi, _currentgFactorOmegaPhi );
- }
- }
-
- result *=dynModel;
-
- if ( _cacheAmps){
-#ifdef _OPENMP
-#pragma omp critical
- {
-#endif
- _cachedAmpMap[evtNo][lamX]=result;
-#ifdef _OPENMP
- }
-#endif
- }
-
- return result;
-}
-
-complex<double> XToOmegaPhiDecAmps::XToPhiPhiAmp(Spin lamX, EvtDataNew* theData){
-
- complex<double> result(0.,0.);
- std::map< boost::shared_ptr<const JPCLS>, double, pawian::Collection::SharedPtrLess >::iterator itXMag;
-
-
- for ( itXMag=_currentParamMags.begin(); itXMag!=_currentParamMags.end(); ++itXMag){
- boost::shared_ptr<const JPCLS> XState=itXMag->first;
- double theXMag=itXMag->second;
- double theXPhi=_currentParamPhis[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_XToPhiPhi1][XState->J][lamX][lambda]);
-
- amp = amp * phiphiTo4KAmp( theData, lambdaPhi1, lambdaPhi2 );
-
- result +=amp;
- }
- }
- }
-
- return result;
-}
-
-complex<double> XToOmegaPhiDecAmps::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 XToOmegaPhiDecAmps::getDefaultParams(fitParamsNew& fitVal, fitParamsNew& fitErr){
-
- std::vector< boost::shared_ptr<const JPCLS> > PhiPhiStates;
- if(_J_X==0 && _parity==-1) PhiPhiStates=_theStatesPtr->EtaToPhiPhiStates();
- else if(_J_X==0 && _parity==1) PhiPhiStates=_theStatesPtr->F0ToPhiPhiStates();
- else if(_J_X==1 && _parity==1) PhiPhiStates=_theStatesPtr->F1ToPhiPhiStates();
- else if(_J_X==2 && _parity==-1) PhiPhiStates=_theStatesPtr->Eta2ToPhiPhiStates();
- else if(_J_X==2 && _parity==1) PhiPhiStates=_theStatesPtr->F2ToPhiPhiStates();
- 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;
-
- if(_bwMassFit){
- fitVal.Widths[_name]=0.2;
- fitErr.Widths[_name]=0.02;
- }
- else if( _flatteMassFit){
- std::string phiPhiGFactorKey="gFactorPhiPhi_"+_name;
- std::string omegaPhiGFactorKey="gFactorOmegaPhi_"+_name;
- fitVal.gFactors[phiPhiGFactorKey]=0.3;
- fitVal.gFactors[omegaPhiGFactorKey]=0.3;
- fitErr.gFactors[phiPhiGFactorKey]=0.2;
- fitErr.gFactors[omegaPhiGFactorKey]=0.2;
- }
- }
-}
-
-void XToOmegaPhiDecAmps::print(std::ostream& os) const{
- return; //dummy
-}
-
-void XToOmegaPhiDecAmps::initialize(){
- std::vector<std::string>::const_iterator it;
-
- for (it=_hypVec.begin(); it!=_hypVec.end(); ++it){
- if (it->compare(0, _xBWKey.size(), _xBWKey) ==0){
- _bwMassFit=true;
- _massIndependent=false;
- }
- else if (it->compare(0, _xFlatteKey.size(), _xFlatteKey) ==0){
- _flatteMassFit=true;
- _massIndependent=false;
- }
- }
-
-
-}
-
-bool XToOmegaPhiDecAmps::checkRecalculation(fitParamsNew& theParamVal){
- _recalculate=false;
- 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];
- if ( fabs(theXMag - _currentParamMags[XState]) > 1.e-10 ){
- _recalculate=true;
- }
- if ( fabs(theXPhi - _currentParamPhis[XState]) > 1.e-10 ){
- _recalculate=true;
- }
- }
-
- if (!_massIndependent){
- double xMass=theParamVal.Masses[_name];
- if ( fabs(xMass-_currentXMass) > 1.e-10){
- _recalculate=true;
- }
-
- if(_bwMassFit){
- double xWidth=theParamVal.Widths[_name];
- if ( fabs(xWidth - _currentXWidth) > 1.e-10){
- _recalculate=true;
- }
- }
- else if(_flatteMassFit){
- std::string phiPhiGFactorKey="gFactorPhiPhi_"+_name;
- std::string omegaPhiFactorKey="gFactorOmegaPhi_"+_name;
- double gFactorPhiPhi=theParamVal.gFactors[phiPhiGFactorKey];
- double gFactorOmegaPhi=theParamVal.gFactors[omegaPhiFactorKey];
- if ( fabs(gFactorPhiPhi - _currentgFactorPhiPhi) > 1.e-10 ){
- _recalculate=true;
- }
- if ( fabs(gFactorOmegaPhi - _currentgFactorOmegaPhi) > 1.e-10 ){
- _recalculate=true;
- }
- }
- }
- if (_recalculate) Info << "Recalculate amplitude:\t" << _name << endmsg;
-
- return _recalculate;
-}
-
-
-void XToOmegaPhiDecAmps::updateFitParams(fitParamsNew& theParamVal){
- 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];
- _currentParamMags[XState]=theXMag;
- _currentParamPhis[XState]=theXPhi;
- }
-
- if (!_massIndependent){
- double xMass=theParamVal.Masses[_name];
- _currentXMass= xMass;
-
- if(_bwMassFit){
- double xWidth=theParamVal.Widths[_name];
- _currentXWidth=xWidth;
- }
- else if(_flatteMassFit){
- std::string phiPhiGFactorKey="gFactorPhiPhi_"+_name;
- std::string omegaPhiFactorKey="gFactorOmegaPhi_"+_name;
- double gFactorPhiPhi=theParamVal.gFactors[phiPhiGFactorKey];
- double gFactorOmegaPhi=theParamVal.gFactors[omegaPhiFactorKey];
- _currentgFactorPhiPhi=gFactorPhiPhi;
- _currentgFactorOmegaPhi=gFactorOmegaPhi;
- }
- }
-}
diff --git a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.hh b/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.hh
index 2d88afdc2cdb1810c6cf50ed7020198baed36400..6174638e3bf3127554c858758d65f5ab7a12912a 100644
--- a/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.hh
+++ b/Examples/JpsiToOmegaPhiGam/XToOmegaPhiDecAmps.hh
@@ -51,8 +51,8 @@ protected:
boost::shared_ptr<JpsiToOmegaPhiGamStates> _theStatesPtr;
- complex<double> XToPhiPhiAmp(Spin lamX, EvtDataNew* theData);
- complex<double> phiphiTo4KAmp( EvtDataNew* theData, Spin lambdaPhi1, Spin lambdaPhi2 );
+ complex<double> XToPhiOmegaAmp(Spin lamX, EvtDataNew* theData);
+ complex<double> phiOmegaTo2K3piAmp( EvtDataNew* theData, Spin lambdaOmega, Spin lambdaPhi );
virtual void initialize();
diff --git a/Examples/JpsiToOmegaPhiGam/defaultparams.dat b/Examples/JpsiToOmegaPhiGam/defaultparams.dat
deleted file mode 100644
index 1945a9dc6639777f9fdaa26fd28987e46e22be89..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/defaultparams.dat
+++ /dev/null
@@ -1,4 +0,0 @@
-J1P-1C-1Lama0Lamb1GammaEta_2981Mag 0.7 0.6
-J1P-1C-1Lama0Lamb1GammaEta_2981Phi 0 0.3
-J0P-1C1L1S1Eta_2981ToPhiPhiMag 0.2 0.8
-J0P-1C1L1S1Eta_2981ToPhiPhiPhi 0 0.3
diff --git a/Examples/JpsiToOmegaPhiGam/qaSummary.dat b/Examples/JpsiToOmegaPhiGam/qaSummary.dat
deleted file mode 100644
index 082ca6c949578e52a59d51607869be9a399c3187..0000000000000000000000000000000000000000
--- a/Examples/JpsiToOmegaPhiGam/qaSummary.dat
+++ /dev/null
@@ -1,5 +0,0 @@
-BIC 12195.7
-AICa 12170.2
-AICc 12170.2
-logLh 6081.08
-free parameter 4
diff --git a/Examples/JpsiToPhiPhiGam/ConfigFiles/startParamsEta.cfg b/Examples/JpsiToPhiPhiGam/ConfigFiles/startParamsEta.cfg
index 37e6e50975e7724103fec7d3a6221b0393bdcdcb..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/Examples/JpsiToPhiPhiGam/ConfigFiles/startParamsEta.cfg
+++ b/Examples/JpsiToPhiPhiGam/ConfigFiles/startParamsEta.cfg
@@ -1,4 +0,0 @@
-J1P-1C-1Lama0Lamb1GammaEta_2000Mag 0.7 0.6
-J1P-1C-1Lama0Lamb1GammaEta_2000Phi 0 0.3
-J0P-1C1L1S1Eta_2000ToPhiPhiMag 0.2 0.8
-J0P-1C1L1S1Eta_2000ToPhiPhiPhi 0 0.3