From ae5bdd189737911f4f29905719d9afbdd93a71a3 Mon Sep 17 00:00:00 2001 From: amfodajlab Date: Wed, 29 Jul 2020 11:39:55 -0600 Subject: [PATCH 1/4] include omegapi_amplitude --- src/programs/Simulation/gen_omegapi/gen_omegapi.cc | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/programs/Simulation/gen_omegapi/gen_omegapi.cc b/src/programs/Simulation/gen_omegapi/gen_omegapi.cc index ac18245fb..c0b0eaea5 100644 --- a/src/programs/Simulation/gen_omegapi/gen_omegapi.cc +++ b/src/programs/Simulation/gen_omegapi/gen_omegapi.cc @@ -16,6 +16,7 @@ #include "AMPTOOLS_DATAIO/HDDMDataWriter.h" #include "AMPTOOLS_DATAIO/ASCIIDataWriter.h" +#include "AMPTOOLS_AMPS/omegapi_amplitude.h" #include "AMPTOOLS_AMPS/omegapiAngAmp.h" #include "AMPTOOLS_AMPS/omegapiAngles.h" #include "AMPTOOLS_AMPS/BreitWigner.h" @@ -223,6 +224,7 @@ int main( int argc, char* argv[] ){ cout << "TRandom3 Seed : " << seed << endl; // setup AmpToolsInterface + AmpToolsInterface::registerAmplitude( omegapi_amplitude() ); AmpToolsInterface::registerAmplitude( omegapiAngAmp() ); AmpToolsInterface::registerAmplitude( BreitWigner() ); AmpToolsInterface::registerAmplitude( Uniform() ); From 4e79568a0c9952f7e7975f48134cffbffd385533 Mon Sep 17 00:00:00 2001 From: amfodajlab Date: Wed, 29 Jul 2020 11:43:57 -0600 Subject: [PATCH 2/4] adding omegapi_amplitude --- .../AMPTOOLS_AMPS/omegapi_amplitude.cc | 219 ++++++++++++++++++ .../AMPTOOLS_AMPS/omegapi_amplitude.h | 100 ++++++++ 2 files changed, 319 insertions(+) create mode 100644 src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc create mode 100644 src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.h diff --git a/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc new file mode 100644 index 000000000..95aa9ffdf --- /dev/null +++ b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc @@ -0,0 +1,219 @@ +//Jan 18th 2020, Based on model by Adam Szczepaniak & Vincent Mathieu +#include +#include +#include + +#include +#include +#include +#include +//#include "UTILITIES/CobremsGeneration.hh" +//#include "UTILITIES/BeamProperties.h" + +#include "TLorentzVector.h" +#include "TLorentzRotation.h" + +#include "IUAmpTools/AmpParameter.h" +#include "omegapi_amplitude.h" +#include "barrierFactor.h" +#include "clebschGordan.h" +#include "wignerD.h" +#include "breakupMomentum.h" +#include "omegapiAngles.h" + +#include +#include +#include +#include "TMath.h" + + +omegapi_amplitude::omegapi_amplitude( const vector< string >& args ): + UserAmplitude< omegapi_amplitude >( args ) +{ + assert( args.size() == (7+4+2) || args.size() == (7+4+3) ); + + if(args.size() == (7+4+3)){ + polAngle = atof(args[7+4+1].c_str() ); // azimuthal angle of the photon polarization vector in the lab measured in degrees. + polFraction = AmpParameter( args[7+4+2] ); // polarization fraction + std::cout << "Fixed polarization fraction =" << polFraction << " and pol.angle= " << polAngle << " degrees." << std::endl; + } +/* else if (args.size() == (7+4+2)){//beam properties requires halld_sim + // BeamProperties configuration file + TString beamConfigFile = args[7+4+1].c_str(); + BeamProperties beamProp(beamConfigFile); + polFrac_vs_E = (TH1D*)beamProp.GetPolFrac(); + polAngle = beamProp.GetPolAngle(); + std::cout << "Polarisation angle of " << polAngle << " from BeamProperties." << std::endl; + if(polAngle == -1) + std::cout << "This is an amorphous run. Set beam polarisation to 0." << std::endl; + for(Int_t i=0; iGetXaxis()->GetNbins()+2; i++){ + //cout << polFrac_vs_E->GetBinContent(i) << endl; + } + }*/ + else + assert(0); + + sign = atoi(args[0].c_str() ); + lambda_gamma = atoi(args[1].c_str() ); + spin = atoi(args[2].c_str() ); + parity = atoi(args[3].c_str() ); + spin_proj = atoi(args[4].c_str() ); + + c_0 = AmpParameter(args[5]); + registerParameter(c_0); + + c_1 = AmpParameter(args[6]); + registerParameter(c_1); + + c_2 = AmpParameter(args[7]); + registerParameter(c_2); + + //Dalitz Parameters + dalitz_alpha = AmpParameter(args[7+1]); + dalitz_beta = AmpParameter(args[7+2]); + dalitz_gamma = AmpParameter(args[7+3]); + dalitz_delta = AmpParameter(args[7+4]); + + registerParameter(dalitz_alpha); + registerParameter(dalitz_beta); + registerParameter(dalitz_gamma); + registerParameter(dalitz_delta); + +} +////////////////////////////////////////////////// User Vars ////////////////////////////////// +void +omegapi_amplitude::calcUserVars( GDouble** pKin, GDouble* userVars ) const +{ + + TLorentzVector beam (pKin[0][1], pKin[0][2], pKin[0][3], pKin[0][0]); + TLorentzVector recoil(pKin[1][1], pKin[1][2], pKin[1][3], pKin[1][0]); + + TLorentzVector rhos_pip(pKin[4][1], pKin[4][2], pKin[4][3], pKin[4][0]); + TLorentzVector rhos_pim(pKin[5][1], pKin[5][2], pKin[5][3], pKin[5][0]); + TLorentzVector rho = rhos_pip + rhos_pim; + + TLorentzVector omegas_pi(pKin[3][1], pKin[3][2], pKin[3][3], pKin[3][0]); + TLorentzVector omega = rho + omegas_pi; + + TLorentzVector Xs_pi(pKin[2][1], pKin[2][2], pKin[2][3], pKin[2][0]); + + TLorentzVector X = omega + Xs_pi; + + //////////////////////// Boost Particles and Get Angles////////////////////////////////// + + TLorentzVector target(0,0,0,0.938); + //Helicity coordinate system + TLorentzVector Gammap = beam + target; + +// polarization BeamProperties + GDouble Pgamma=polFraction;//fixed beam polarization fraction + if(polAngle == -1) + Pgamma = 0.;//if beam is amorphous set polarization fraction to 0 +/* else if(polFrac_vs_E!=NULL){ + //This part causes seg fault with 34 amplitudes or more with gen_amp and gen_omegapi. + //Not needed for fixed beam pol angle and frac. + int bin = polFrac_vs_E->GetXaxis()->FindBin(beam.E()); + + if (bin == 0 || bin > polFrac_vs_E->GetXaxis()->GetNbins()){ + Pgamma = 0.; + } + else + Pgamma = polFrac_vs_E->GetBinContent(bin); + }*/ + + //Calculate decay angles in helicity frame + vector locthetaphi = getomegapiAngles(polAngle, omega, X, beam, Gammap); + + vector locthetaphih = getomegapiAngles(rhos_pip, omega, X, Gammap, rhos_pim); + + userVars[uv_cosTheta] = TMath::Cos(locthetaphi[0]); + userVars[uv_Phi] = locthetaphi[1]; + + userVars[uv_cosThetaH] = TMath::Cos(locthetaphih[0]); + userVars[uv_PhiH] = locthetaphih[1]; + + userVars[uv_prod_angle] = locthetaphi[2]; + + userVars[uv_Pgamma] = Pgamma; + +///////////////////////////////////////////// Dalitz Parameters /////////////////////////////// + double dalitz_s = rho.M2();//s=M2(pip pim) + double dalitz_t = (rhos_pip+omegas_pi).M2();//t=M2(pip pi0) + double dalitz_u = (rhos_pim+omegas_pi).M2();//u=M2(pim pi0) + double m3pi = (2*139.57018)+134.9766; + double dalitz_d = 2*omega.M()*( omega.M() - m3pi); + double dalitz_sc = (1/3)*( omega.M2() - rhos_pip.M2() - rhos_pim.M2() - omegas_pi.M2()); + double dalitzx = sqrt(3)*(dalitz_t - dalitz_u)/dalitz_d; + double dalitzy = 3*(dalitz_sc - dalitz_s)/dalitz_d; + double dalitz_z = dalitzx*dalitzx + dalitzy*dalitzy; + double dalitz_sin3theta = TMath::Sin(3 * TMath::ASin( (dalitzy/sqrt(dalitz_z) )) ); + + userVars[uv_dalitz_z] = dalitz_z; + userVars[uv_dalitz_sin3theta] = dalitz_sin3theta; + +} + +////////////////////////////////////////////////// Amplitude Calculation ////////////////////////////////// + +complex< GDouble > +omegapi_amplitude::calcAmplitude( GDouble** pKin, GDouble* userVars ) const +{ + + complex COne(1,0); + complex CZero(0,0); + + GDouble cosTheta = userVars[uv_cosTheta]; + GDouble Phi = userVars[uv_Phi]; + GDouble cosThetaH = userVars[uv_cosThetaH]; + GDouble PhiH = userVars[uv_PhiH]; + GDouble prod_angle = userVars[uv_prod_angle]; + GDouble polfrac = userVars[uv_Pgamma]; + GDouble dalitz_z = userVars[uv_dalitz_z]; + GDouble dalitz_sin3theta = userVars[uv_dalitz_sin3theta]; + + GDouble G = sqrt(1 + 2 * dalitz_alpha * dalitz_z + 2 * dalitz_beta * pow(dalitz_z,3/2) * dalitz_sin3theta + + 2 * dalitz_gamma * pow(dalitz_z,2) + 2 * dalitz_delta * pow(dalitz_z,5/2) * dalitz_sin3theta ); + + GDouble hel_c[3] = { c_0, c_1, c_2}; + + complex amplitude = CZero; + + for (int lambda = -1; lambda <= 1; lambda++)//omega helicity + { + GDouble hel_amp = 0.0; + + for(int l = 0; l <= 2; l++)//partial waves (l). + {//if ( (parity == -1 && l% 2 == 0) || (parity == 1 && l%2 != 0) ) continue; + + hel_amp += hel_c[l] * clebschGordan(l, 1, 0, lambda, spin, lambda); + }//loop over l + + amplitude += wignerD( spin, spin_proj, lambda, cosTheta, Phi ) * hel_amp * wignerD( 1, lambda, 0, cosThetaH, PhiH ) * G; + }//loop over lambda + + + complex prefactor ( cos( lambda_gamma * prod_angle ), sin( lambda_gamma * prod_angle )); + + if (sign == -1 && lambda_gamma == -1){amplitude *= -1*sqrt( ( 1 + (sign * polfrac) )/2 ) * prefactor;} + else{amplitude *= sqrt( ( 1 + (sign * polfrac) )/2 ) * prefactor;} + +//Instead of the vertices "scale" in configuration file. +// if( (parity == +1 && ((spin+spin_proj)%2 == 0) ) || (parity == -1 && ((spin+spin_proj)%2 != 0) ) ){amplitude *= -1;} + +return amplitude; +} + +void omegapi_amplitude::updatePar( const AmpParameter& par ){ + + // could do expensive calculations here on parameter updates +} + +#ifdef GPU_ACCELERATION +void +omegapi_amplitude::launchGPUKernel( dim3 dimGrid, dim3 dimBlock, GPU_AMP_PROTO ) const { + + GPUomegapi_amplitude_exec( dimGrid, dimBlock, GPU_AMP_ARGS, + sign, lambda_gamma, spin, parity, spin_proj, c_0, c_1, c_2, dalitz_alpha, dalitz_beta, dalitz_gamma, dalitz_delta, polAngle, polFraction); +} +#endif //GPU_ACCELERATION + diff --git a/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.h b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.h new file mode 100644 index 000000000..dbb29efe6 --- /dev/null +++ b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.h @@ -0,0 +1,100 @@ +//Jan 18th 2020, Based on model by Adam Szczepaniak & Vincent Mathieu +#if !defined(OMEGAPI_AMPLITUDE) +#define OMEGAPI_AMPLITUDE + +#include "IUAmpTools/Amplitude.h" +#include "IUAmpTools/AmpParameter.h" +#include "IUAmpTools/UserAmplitude.h" +#include "GPUManager/GPUCustomTypes.h" + +#include +#include +#include + +#include "TLorentzVector.h" +#include "TH1D.h" +#include "TFile.h" + +#ifdef GPU_ACCELERATION +void GPUomegapi_amplitude_exec( dim3 dimGrid, dim3 dimBlock, GPU_AMP_PROTO, + int sign, int lambda_gamma, int spin, int parity, int spin_proj, GDouble c_0, GDouble c_1, GDouble c_2, + GDouble dalitz_alpha, GDouble dalitz_beta, GDouble dalitz_gamma, GDouble dalitz_delta, GDouble polAngle, GDouble polFraction); + +#endif // GPU_ACCELERATION + +using std::complex; +using namespace std; + +class Kinematics; + +class omegapi_amplitude : public UserAmplitude< omegapi_amplitude > +{ + +public: + + omegapi_amplitude() : UserAmplitude< omegapi_amplitude >() { } + omegapi_amplitude( const vector< string >& args ); + ~omegapi_amplitude(){} + + string name() const { return "omegapi_amplitude"; } + +complex< GDouble > calcAmplitude( GDouble** pKin, GDouble* userVars ) const; +//complex< GDouble > calcAmplitude( GDouble** pKin ) const; + + // ********************** + // The following lines are optional and can be used to precalcualte + // user-defined data that the amplitudes depend on. + + // Use this for indexing a user-defined data array and notifying + // the framework of the number of user-defined variables. + + enum UserVars { uv_cosTheta = 0, uv_Phi = 1, uv_cosThetaH = 2, uv_PhiH = 3, uv_prod_angle = 4, uv_Pgamma = 5, uv_dalitz_z = 6, uv_dalitz_sin3theta = 7, kNumUserVars }; + unsigned int numUserVars() const { return kNumUserVars; } + + // This function needs to be defined -- see comments and discussion + // in the .cc file. + void calcUserVars( GDouble** pKin, GDouble* userVars ) const; + + // This is an optional addition if the calcAmplitude routine + // can run with only the user-defined data and not the original + // four-vectors. It is used to optimize memory usage in GPU + // based fits. + bool needsUserVarsOnly() const { return true; } + // ** end of optional lines ** + + void updatePar( const AmpParameter& par ); + +#ifdef GPU_ACCELERATION + + void launchGPUKernel( dim3 dimGrid, dim3 dimBlock, GPU_AMP_PROTO ) const; + + bool isGPUEnabled() const { return true; } + +#endif // GPU_ACCELERATION + +private: + + int sign; + int lambda_gamma; + int spin; + int parity; + int spin_proj; + + AmpParameter c_0; + AmpParameter c_1; + AmpParameter c_2; + + AmpParameter dalitz_alpha; + AmpParameter dalitz_beta; + AmpParameter dalitz_gamma; + AmpParameter dalitz_delta; + + double polAngle, polFraction; + + TH1D *totalFlux_vs_E; + TH1D *polFlux_vs_E; + TH1D *polFrac_vs_E; + +}; + +#endif From b696563037e4cf54065d675d1a8244162de860b5 Mon Sep 17 00:00:00 2001 From: amfodajlab Date: Wed, 29 Jul 2020 11:45:28 -0600 Subject: [PATCH 3/4] Update omegapi_amplitude.cc --- src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc index 95aa9ffdf..a499a0618 100644 --- a/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc +++ b/src/libraries/AMPTOOLS_AMPS/omegapi_amplitude.cc @@ -7,8 +7,8 @@ #include #include #include -//#include "UTILITIES/CobremsGeneration.hh" -//#include "UTILITIES/BeamProperties.h" +#include "UTILITIES/CobremsGeneration.hh" +#include "UTILITIES/BeamProperties.h" #include "TLorentzVector.h" #include "TLorentzRotation.h" @@ -37,7 +37,7 @@ omegapi_amplitude::omegapi_amplitude( const vector< string >& args ): polFraction = AmpParameter( args[7+4+2] ); // polarization fraction std::cout << "Fixed polarization fraction =" << polFraction << " and pol.angle= " << polAngle << " degrees." << std::endl; } -/* else if (args.size() == (7+4+2)){//beam properties requires halld_sim + else if (args.size() == (7+4+2)){//beam properties requires halld_sim // BeamProperties configuration file TString beamConfigFile = args[7+4+1].c_str(); BeamProperties beamProp(beamConfigFile); @@ -49,7 +49,7 @@ omegapi_amplitude::omegapi_amplitude( const vector< string >& args ): for(Int_t i=0; iGetXaxis()->GetNbins()+2; i++){ //cout << polFrac_vs_E->GetBinContent(i) << endl; } - }*/ + } else assert(0); @@ -109,7 +109,7 @@ omegapi_amplitude::calcUserVars( GDouble** pKin, GDouble* userVars ) const GDouble Pgamma=polFraction;//fixed beam polarization fraction if(polAngle == -1) Pgamma = 0.;//if beam is amorphous set polarization fraction to 0 -/* else if(polFrac_vs_E!=NULL){ + else if(polFrac_vs_E!=NULL){ //This part causes seg fault with 34 amplitudes or more with gen_amp and gen_omegapi. //Not needed for fixed beam pol angle and frac. int bin = polFrac_vs_E->GetXaxis()->FindBin(beam.E()); @@ -119,7 +119,7 @@ omegapi_amplitude::calcUserVars( GDouble** pKin, GDouble* userVars ) const } else Pgamma = polFrac_vs_E->GetBinContent(bin); - }*/ + } //Calculate decay angles in helicity frame vector locthetaphi = getomegapiAngles(polAngle, omega, X, beam, Gammap); From 6df924d231cf47636301bd926a132d5b0326303f Mon Sep 17 00:00:00 2001 From: amfodajlab Date: Wed, 29 Jul 2020 11:52:09 -0600 Subject: [PATCH 4/4] example omegapi_amplitude config file --- .../gen_omegapi/gen_omegapi_amplitude_all.cfg | 375 ++++++++++++++++++ 1 file changed, 375 insertions(+) create mode 100644 src/programs/Simulation/gen_omegapi/gen_omegapi_amplitude_all.cfg diff --git a/src/programs/Simulation/gen_omegapi/gen_omegapi_amplitude_all.cfg b/src/programs/Simulation/gen_omegapi/gen_omegapi_amplitude_all.cfg new file mode 100644 index 000000000..c08d9b618 --- /dev/null +++ b/src/programs/Simulation/gen_omegapi/gen_omegapi_amplitude_all.cfg @@ -0,0 +1,375 @@ +# This config file generates/fits omega pi0 for up tp spin 2 states based on amplitude analysis expansion by model by Adam Szczepaniak & Vincent Mathieu +##################################### +#### THIS IS A CONFIG FILE #### +##################################### +## +## Blank lines or lines beginning with a "#" are ignored. +## +## Double colons (::) are treated like a space. +## This is sometimes useful for grouping (for example, +## grouping strings like "reaction::sum::amplitudeName") +## +## All non-comment lines must begin with one of the following keywords. +## +## (note: means necessary +## (word) means optional) +## +## include +## define (defn1) (defn2) (defn3) ... +## fit +## keyword +## reaction (particle3) ... +## data (arg1) (arg2) (arg3) ... +## genmc (arg1) (arg2) (arg3) ... +## accmc (arg1) (arg2) (arg3) ... +## normintfile +## sum (sum2) (sum3) ... +## amplitude (arg1) (arg2) ([par]) ... +## initialize <"events"/"polar"/"cartesian"> +## ("fixed"/"real") +## scale +## constrain ... +## permute ... +## parameter ("fixed"/"bounded"/"gaussian") +## (lower/central) (upper/error) +## DEPRECATED: +## datafile (file2) (file3) ... +## genmcfile (file2) (file3) ... +## accmcfile (file2) (file3) ... +## +##################################### + +# useful masses and widths +define b1 1.235 0.142 + +parameter hel_c_0_m_1 1.0 + +parameter hel_c_1_m_1 1.0 + +parameter hel_c_1_p_0 1.0 +parameter hel_c_1_p_2 0.27 + +parameter hel_c_2_m_1 1.0 +parameter hel_c_2_p_2 1.0 + +parameter dalitz_alpha 94.0 +#94 84 80 83 #136 125 113 114 + +parameter dalitz_beta 0.0 fixed +#0 28 27 22 #0 30 27 24 + +parameter dalitz_gamma 0.0 fixed +#0 0 8 1 #0 0 24 20 + +parameter dalitz_delta 0.0 fixed +#0 0 0 14 #0 0 0 6 + + +define l_0m 0.0 [hel_c_0_m_1] 0.0 +define l_1m 0.0 [hel_c_1_m_1] 0.0 +define l_1p [hel_c_1_p_0] 0.0 [hel_c_1_p_2] +define l_2m 0.0 [hel_c_2_m_1] 0.0 +define l_2p 0.0 0.0 [hel_c_2_p_2] + +define dalitz [dalitz_alpha] [dalitz_beta] [dalitz_gamma] [dalitz_delta] + +fit omegapi + +define trange 0.1 0.3 +define erange 8.3 9.3 +define mrange 1.165 1.3 + +# some definitions for adjusting the beam polarization +#define beamconfig TEMPBEAMCONFIG +define angle 0.0 +define fraction 0.4 + +#When generating use this reaction to account for indistinguishable pi0s +reaction omegapi Beam Proton Pi0 Pi0 Pi+ Pi- + +#with fit use the following reaction to avoid indistinguishable pi0s +#reaction omegapi Beam Proton Pi01 Pi02 Pi+ Pi- + + +# Sums +sum omegapi term1 term2 + +genmc omegapi ROOTDataReader phasespace.root +accmc omegapi ROOTDataReader phasespace.root +data omegapi ROOTDataReader zerom.root + +############################################################# +#amplitude omegapi::term1::amp omegapi_amplitude sign lambda_gamma spin parity spin_proj l_1p angle fraction + +############################ spin 0 parity -1 ################################## +amplitude omegapi::term1::pp0m0 omegapi_amplitude 1 1 0 -1 0 l_0m dalitz angle fraction +amplitude omegapi::term1::pm0m0 omegapi_amplitude 1 -1 0 -1 0 l_0m dalitz angle fraction + +amplitude omegapi::term2::mp0m0 omegapi_amplitude -1 1 0 -1 0 l_0m dalitz angle fraction +amplitude omegapi::term2::mm0m0 omegapi_amplitude -1 -1 0 -1 0 l_0m dalitz angle fraction + +initialize omegapi::term1::pp0m0 cartesian 1.0 0.0 +initialize omegapi::term1::pm0m0 cartesian 1.0 0.0 + +initialize omegapi::term2::mp0m0 cartesian 1.0 0.0 +initialize omegapi::term2::mm0m0 cartesian 1.0 0.0 + +constrain omegapi term1 pp0m0 omegapi term1 pm0m0 +constrain omegapi term2 mp0m0 omegapi term2 mm0m0 + +constrain omegapi term1 pp0m0 omegapi term2 mp0m0 + +############################ spin 1 parity -1 ################################## +amplitude omegapi::term1::pp1mp omegapi_amplitude 1 1 1 -1 1 l_1m dalitz angle fraction +amplitude omegapi::term1::pp1m0 omegapi_amplitude 1 1 1 -1 0 l_1m dalitz angle fraction +amplitude omegapi::term1::pp1mm omegapi_amplitude 1 1 1 -1 -1 l_1m dalitz angle fraction + +amplitude omegapi::term1::pm1mp omegapi_amplitude 1 -1 1 -1 1 l_1m dalitz angle fraction +amplitude omegapi::term1::pm1m0 omegapi_amplitude 1 -1 1 -1 0 l_1m dalitz angle fraction +amplitude omegapi::term1::pm1mm omegapi_amplitude 1 -1 1 -1 -1 l_1m dalitz angle fraction + +amplitude omegapi::term2::mp1mp omegapi_amplitude -1 1 1 -1 1 l_1m dalitz angle fraction +amplitude omegapi::term2::mp1m0 omegapi_amplitude -1 1 1 -1 0 l_1m dalitz angle fraction +amplitude omegapi::term2::mp1mm omegapi_amplitude -1 1 1 -1 -1 l_1m dalitz angle fraction + +amplitude omegapi::term2::mm1mp omegapi_amplitude -1 -1 1 -1 1 l_1m dalitz angle fraction +amplitude omegapi::term2::mm1m0 omegapi_amplitude -1 -1 1 -1 0 l_1m dalitz angle fraction +amplitude omegapi::term2::mm1mm omegapi_amplitude -1 -1 1 -1 -1 l_1m dalitz angle fraction + +initialize omegapi::term1::pp1mp cartesian 1.0 0.0 +initialize omegapi::term1::pp1m0 cartesian 1.0 0.0 +initialize omegapi::term1::pp1mm cartesian 1.0 0.0 + +initialize omegapi::term1::pm1mp cartesian 1.0 0.0 +initialize omegapi::term1::pm1m0 cartesian 1.0 0.0 +initialize omegapi::term1::pm1mm cartesian 1.0 0.0 + +initialize omegapi::term2::mp1mp cartesian 1.0 0.0 +initialize omegapi::term2::mp1m0 cartesian 1.0 0.0 +initialize omegapi::term2::mp1mm cartesian 1.0 0.0 + +initialize omegapi::term2::mm1mp cartesian 1.0 0.0 +initialize omegapi::term2::mm1m0 cartesian 1.0 0.0 +initialize omegapi::term2::mm1mm cartesian 1.0 0.0 + + +constrain omegapi term1 pp1mp omegapi term1 pm1mm +constrain omegapi term1 pp1m0 omegapi term1 pm1m0 +constrain omegapi term1 pp1mm omegapi term1 pm1mp + +constrain omegapi term2 mp1mp omegapi term2 mm1mm +constrain omegapi term2 mp1m0 omegapi term2 mm1m0 +constrain omegapi term2 mp1mm omegapi term2 mm1mp + +constrain omegapi term1 pp1mp omegapi term2 mp1mp +constrain omegapi term1 pp1m0 omegapi term2 mp1m0 +constrain omegapi term1 pp1mm omegapi term2 mp1mm + +scale omegapi term1 pm1m0 -1.0 + +scale omegapi term2 mm1m0 -1.0 + +############################ spin 1 parity +1 ################################## +amplitude omegapi::term1::pp1pp omegapi_amplitude 1 1 1 1 1 l_1p dalitz angle fraction +amplitude omegapi::term1::pp1p0 omegapi_amplitude 1 1 1 1 0 l_1p dalitz angle fraction +amplitude omegapi::term1::pp1pm omegapi_amplitude 1 1 1 1 -1 l_1p dalitz angle fraction + +amplitude omegapi::term1::pm1pp omegapi_amplitude 1 -1 1 1 1 l_1p dalitz angle fraction +amplitude omegapi::term1::pm1p0 omegapi_amplitude 1 -1 1 1 0 l_1p dalitz angle fraction +amplitude omegapi::term1::pm1pm omegapi_amplitude 1 -1 1 1 -1 l_1p dalitz angle fraction + +amplitude omegapi::term2::mp1pp omegapi_amplitude -1 1 1 1 1 l_1p dalitz angle fraction +amplitude omegapi::term2::mp1p0 omegapi_amplitude -1 1 1 1 0 l_1p dalitz angle fraction +amplitude omegapi::term2::mp1pm omegapi_amplitude -1 1 1 1 -1 l_1p dalitz angle fraction + +amplitude omegapi::term2::mm1pp omegapi_amplitude -1 -1 1 1 1 l_1p dalitz angle fraction +amplitude omegapi::term2::mm1p0 omegapi_amplitude -1 -1 1 1 0 l_1p dalitz angle fraction +amplitude omegapi::term2::mm1pm omegapi_amplitude -1 -1 1 1 -1 l_1p dalitz angle fraction + + +initialize omegapi::term1::pp1pp cartesian 1.0 0.0 +initialize omegapi::term1::pp1p0 cartesian 1.0 0.0 +initialize omegapi::term1::pp1pm cartesian 1.0 0.0 + +initialize omegapi::term1::pm1pp cartesian 1.0 0.0 +initialize omegapi::term1::pm1p0 cartesian 1.0 0.0 +initialize omegapi::term1::pm1pm cartesian 1.0 0.0 + +initialize omegapi::term2::mp1pp cartesian 1.0 0.0 +initialize omegapi::term2::mp1p0 cartesian 1.0 0.0 +initialize omegapi::term2::mp1pm cartesian 1.0 0.0 + +initialize omegapi::term2::mm1pp cartesian 1.0 0.0 +initialize omegapi::term2::mm1p0 cartesian 1.0 0.0 +initialize omegapi::term2::mm1pm cartesian 1.0 0.0 + + +constrain omegapi term1 pp1pp omegapi term1 pm1pm +constrain omegapi term1 pp1p0 omegapi term1 pm1p0 +constrain omegapi term1 pp1pm omegapi term1 pm1pp + +constrain omegapi term2 mp1pp omegapi term2 mm1pm +constrain omegapi term2 mp1p0 omegapi term2 mm1p0 +constrain omegapi term2 mp1pm omegapi term2 mm1pp + +constrain omegapi term1 pp1pp omegapi term2 mp1pp +constrain omegapi term1 pp1p0 omegapi term2 mp1p0 +constrain omegapi term1 pp1pm omegapi term2 mp1pm + +scale omegapi term1 pm1pm -1.0 +scale omegapi term1 pm1pp -1.0 + +scale omegapi term2 mm1pm -1.0 +scale omegapi term2 mm1pp -1.0 + +############################ spin 2 parity -1 ################################## +amplitude omegapi::term1::pp2mp2 omegapi_amplitude 1 1 2 -1 2 l_2m dalitz angle fraction +amplitude omegapi::term1::pp2mp omegapi_amplitude 1 1 2 -1 1 l_2m dalitz angle fraction +amplitude omegapi::term1::pp2m0 omegapi_amplitude 1 1 2 -1 0 l_2m dalitz angle fraction +amplitude omegapi::term1::pp2mm omegapi_amplitude 1 1 2 -1 -1 l_2m dalitz angle fraction +amplitude omegapi::term1::pp2mm2 omegapi_amplitude 1 1 2 -1 -2 l_2m dalitz angle fraction + +amplitude omegapi::term1::pm2mp2 omegapi_amplitude 1 -1 2 -1 2 l_2m dalitz angle fraction +amplitude omegapi::term1::pm2mp omegapi_amplitude 1 -1 2 -1 1 l_2m dalitz angle fraction +amplitude omegapi::term1::pm2m0 omegapi_amplitude 1 -1 2 -1 0 l_2m dalitz angle fraction +amplitude omegapi::term1::pm2mm omegapi_amplitude 1 -1 2 -1 -1 l_2m dalitz angle fraction +amplitude omegapi::term1::pm2mm2 omegapi_amplitude 1 -1 2 -1 -2 l_2m dalitz angle fraction + +amplitude omegapi::term2::mp2mp2 omegapi_amplitude -1 1 2 -1 2 l_2m dalitz angle fraction +amplitude omegapi::term2::mp2mp omegapi_amplitude -1 1 2 -1 1 l_2m dalitz angle fraction +amplitude omegapi::term2::mp2m0 omegapi_amplitude -1 1 2 -1 0 l_2m dalitz angle fraction +amplitude omegapi::term2::mp2mm omegapi_amplitude -1 1 2 -1 -1 l_2m dalitz angle fraction +amplitude omegapi::term2::mp2mm2 omegapi_amplitude -1 1 2 -1 -2 l_2m dalitz angle fraction + +amplitude omegapi::term2::mm2mp2 omegapi_amplitude -1 -1 2 -1 2 l_2m dalitz angle fraction +amplitude omegapi::term2::mm2mp omegapi_amplitude -1 -1 2 -1 1 l_2m dalitz angle fraction +amplitude omegapi::term2::mm2m0 omegapi_amplitude -1 -1 2 -1 0 l_2m dalitz angle fraction +amplitude omegapi::term2::mm2mm omegapi_amplitude -1 -1 2 -1 -1 l_2m dalitz angle fraction +amplitude omegapi::term2::mm2mm2 omegapi_amplitude -1 -1 2 -1 -2 l_2m dalitz angle fraction + + +initialize omegapi::term1::pp2mp2 cartesian 1.0 0.0 +initialize omegapi::term1::pp2mp cartesian 1.0 0.0 +initialize omegapi::term1::pp2m0 cartesian 1.0 0.0 +initialize omegapi::term1::pp2mm cartesian 1.0 0.0 +initialize omegapi::term1::pp2mm2 cartesian 1.0 0.0 + +initialize omegapi::term1::pm2mp2 cartesian 1.0 0.0 +initialize omegapi::term1::pm2mp cartesian 1.0 0.0 +initialize omegapi::term1::pm2m0 cartesian 1.0 0.0 +initialize omegapi::term1::pm2mm cartesian 1.0 0.0 +initialize omegapi::term1::pm2mm2 cartesian 1.0 0.0 + +initialize omegapi::term2::mp2mp2 cartesian 1.0 0.0 +initialize omegapi::term2::mp2mp cartesian 1.0 0.0 +initialize omegapi::term2::mp2m0 cartesian 1.0 0.0 +initialize omegapi::term2::mp2mm cartesian 1.0 0.0 +initialize omegapi::term2::mp2mm2 cartesian 1.0 0.0 + +initialize omegapi::term2::mm2mp2 cartesian 1.0 0.0 +initialize omegapi::term2::mm2mp cartesian 1.0 0.0 +initialize omegapi::term2::mm2m0 cartesian 1.0 0.0 +initialize omegapi::term2::mm2mm cartesian 1.0 0.0 +initialize omegapi::term2::mm2mm2 cartesian 1.0 0.0 + + +constrain omegapi term1 pp2mp2 omegapi term1 pm2mm2 +constrain omegapi term1 pp2mp omegapi term1 pm2mm +constrain omegapi term1 pp2m0 omegapi term1 pm2m0 +constrain omegapi term1 pp2mm omegapi term1 pm2mp +constrain omegapi term1 pp2mm2 omegapi term1 pm2mp2 + +constrain omegapi term2 mp2mp2 omegapi term2 mm2mm2 +constrain omegapi term2 mp2mp omegapi term2 mm2mm +constrain omegapi term2 mp2m0 omegapi term2 mm2m0 +constrain omegapi term2 mp2mm omegapi term2 mm2mp +constrain omegapi term2 mp2mm2 omegapi term2 mm2mp2 + +constrain omegapi term1 pp2mp2 omegapi term2 mp2mp2 +constrain omegapi term1 pp2mp omegapi term2 mp2mp +constrain omegapi term1 pp2m0 omegapi term2 mp2m0 +constrain omegapi term1 pp2mm omegapi term2 mp2mm +constrain omegapi term1 pp2mm2 omegapi term2 mp2mm2 + +scale omegapi term1 pm2mm -1.0 +scale omegapi term1 pm2mp -1.0 + +scale omegapi term2 mm2mm -1.0 +scale omegapi term2 mm2mp -1.0 + +############################ spin 2 parity +1 ################################## +amplitude omegapi::term1::pp2pp2 omegapi_amplitude 1 1 2 1 2 l_2p dalitz angle fraction +amplitude omegapi::term1::pp2pp omegapi_amplitude 1 1 2 1 1 l_2p dalitz angle fraction +amplitude omegapi::term1::pp2p0 omegapi_amplitude 1 1 2 1 0 l_2p dalitz angle fraction +amplitude omegapi::term1::pp2pm omegapi_amplitude 1 1 2 1 -1 l_2p dalitz angle fraction +amplitude omegapi::term1::pp2pm2 omegapi_amplitude 1 1 2 1 -2 l_2p dalitz angle fraction + +amplitude omegapi::term1::pm2pp2 omegapi_amplitude 1 -1 2 1 2 l_2p dalitz angle fraction +amplitude omegapi::term1::pm2pp omegapi_amplitude 1 -1 2 1 1 l_2p dalitz angle fraction +amplitude omegapi::term1::pm2p0 omegapi_amplitude 1 -1 2 1 0 l_2p dalitz angle fraction +amplitude omegapi::term1::pm2pm omegapi_amplitude 1 -1 2 1 -1 l_2p dalitz angle fraction +amplitude omegapi::term1::pm2pm2 omegapi_amplitude 1 -1 2 1 -2 l_2p dalitz angle fraction + +amplitude omegapi::term2::mp2pp2 omegapi_amplitude -1 1 2 1 2 l_2p dalitz angle fraction +amplitude omegapi::term2::mp2pp omegapi_amplitude -1 1 2 1 1 l_2p dalitz angle fraction +amplitude omegapi::term2::mp2p0 omegapi_amplitude -1 1 2 1 0 l_2p dalitz angle fraction +amplitude omegapi::term2::mp2pm omegapi_amplitude -1 1 2 1 -1 l_2p dalitz angle fraction +amplitude omegapi::term2::mp2pm2 omegapi_amplitude -1 1 2 1 -2 l_2p dalitz angle fraction + +amplitude omegapi::term2::mm2pp2 omegapi_amplitude -1 -1 2 1 2 l_2p dalitz angle fraction +amplitude omegapi::term2::mm2pp omegapi_amplitude -1 -1 2 1 1 l_2p dalitz angle fraction +amplitude omegapi::term2::mm2p0 omegapi_amplitude -1 -1 2 1 0 l_2p dalitz angle fraction +amplitude omegapi::term2::mm2pm omegapi_amplitude -1 -1 2 1 -1 l_2p dalitz angle fraction +amplitude omegapi::term2::mm2pm2 omegapi_amplitude -1 -1 2 1 -2 l_2p dalitz angle fraction + + +initialize omegapi::term1::pp2pp2 cartesian 1.0 0.0 +initialize omegapi::term1::pp2pp cartesian 1.0 0.0 +initialize omegapi::term1::pp2p0 cartesian 1.0 0.0 +initialize omegapi::term1::pp2pm cartesian 1.0 0.0 +initialize omegapi::term1::pp2pm2 cartesian 1.0 0.0 + +initialize omegapi::term1::pm2pp2 cartesian 1.0 0.0 +initialize omegapi::term1::pm2pp cartesian 1.0 0.0 +initialize omegapi::term1::pm2p0 cartesian 1.0 0.0 +initialize omegapi::term1::pm2pm cartesian 1.0 0.0 +initialize omegapi::term1::pm2pm2 cartesian 1.0 0.0 + +initialize omegapi::term2::mp2pp2 cartesian 1.0 0.0 +initialize omegapi::term2::mp2pp cartesian 1.0 0.0 +initialize omegapi::term2::mp2p0 cartesian 1.0 0.0 +initialize omegapi::term2::mp2pm cartesian 1.0 0.0 +initialize omegapi::term2::mp2pm2 cartesian 1.0 0.0 + +initialize omegapi::term2::mm2pp2 cartesian 1.0 0.0 +initialize omegapi::term2::mm2pp cartesian 1.0 0.0 +initialize omegapi::term2::mm2p0 cartesian 1.0 0.0 +initialize omegapi::term2::mm2pm cartesian 1.0 0.0 +initialize omegapi::term2::mm2pm2 cartesian 1.0 0.0 + + +constrain omegapi term1 pp2pp2 omegapi term1 pm2pm2 +constrain omegapi term1 pp2pp omegapi term1 pm2pm +constrain omegapi term1 pp2p0 omegapi term1 pm2p0 +constrain omegapi term1 pp2pm omegapi term1 pm2pp +constrain omegapi term1 pp2pm2 omegapi term1 pm2pp2 + +constrain omegapi term2 mp2pp2 omegapi term2 mm2pm2 +constrain omegapi term2 mp2pp omegapi term2 mm2pm +constrain omegapi term2 mp2p0 omegapi term2 mm2p0 +constrain omegapi term2 mp2pm omegapi term2 mm2pp +constrain omegapi term2 mp2pm2 omegapi term2 mm2pp2 + +constrain omegapi term1 pp2pp2 omegapi term2 mp2pp2 +constrain omegapi term1 pp2pp omegapi term2 mp2pp +constrain omegapi term1 pp2p0 omegapi term2 mp2p0 +constrain omegapi term1 pp2pm omegapi term2 mp2pm +constrain omegapi term1 pp2pm2 omegapi term2 mp2pm2 + +scale omegapi term1 pm2pm2 -1.0 +scale omegapi term1 pm2p0 -1.0 +scale omegapi term1 pm2pp2 -1.0 + +scale omegapi term2 mm2pm2 -1.0 +scale omegapi term2 mm2p0 -1.0 +scale omegapi term2 mm2pp2 -1.0