jpacPhoto/amplitudes/amplitude.hpp

0001 // Abstract class for an amplitude. Used so we can easily build observables

0002 // as the incoherent sum of amplitudes in s, t, and u channels.

0003 //

0004 // Author:       Daniel Winney (2020)

0005 // Affiliation:  Joint Physics Analysis Center (JPAC)

0006 // Email:        dwinney@iu.edu

0007 // ---------------------------------------------------------------------------

0008
0009 #ifndef _AMPLITUDE_
0010 #define _AMPLITUDE_
0011
0012 // ---------------------------------------------------------------------------

0013 // Abstract class to define helicity amplitudes. This will allow multiple different

0014 // classes (for s, t, and u- channels but also multiple contibutions in each channel)

0015 // to be added together and evaluated in observables.

0016 //

0017 // Any generic amplitude needs a reaction_kinematics object

0018 // and a way to evaluate the helicity amplitude for given set of helicities,

0019 // CoM energy and scattering angle.

0020 // ---------------------------------------------------------------------------

0021
0022 #include "reaction_kinematics.hpp"
0023
0024 #include "Math/GSLIntegrator.h"
0025 #include "Math/IntegrationTypes.h"
0026 #include "Math/Functor.h"
0027
0028 #include <string>
0029 #include <algorithm>
0030
0031 namespace jpacPhoto
0032 {
0033     class amplitude
0034     {
0035         public:
0036         // Constructor with nParams for backward compatibility (now depricated)

0037         amplitude(reaction_kinematics * xkinem, std::string id = "", int n = 0)
0038         : _kinematics(xkinem), _identifier(id)
0039         {};
0040
0041         // Kinematics object for thresholds and etc.

0042         reaction_kinematics * _kinematics;
0043
0044         // saved energies and angle

0045         double _s, _t, _theta;
0046
0047         // Some saveable string by which to identify the amplitude

0048         std::string _identifier;
0049
0050         // How the calculate the helicity amplitude

0051         // Must be given a specific implementation in a user derived class

0052         virtual std::complex<double> helicity_amplitude(std::array<int, 4> helicities, double s, double t) = 0;
0053
0054         // ---------------------------------------------------------------------------

0055         // Observables

0056         // Evaluatable in terms of s and t or an event object (see reaction_kinematics.hpp)

0057
0058         // Modulus of the amplitude summed over all helicity combinations

0059         double probability_distribution(double s, double t);
0060
0061         // Differential and total cross-section

0062         double differential_xsection(double s, double t);
0063
0064         // integrated crossection

0065         double integrated_xsection(double s);
0066
0067         // Spin asymmetries

0068         double A_LL(double s, double t); // Beam and target

0069         double K_LL(double s, double t); // Beam and recoil

0070
0071         // Spin density matrix elements

0072         std::complex<double> SDME(int alpha, int lam, int lamp, double s, double t);
0073
0074         // Beam Asymmetries

0075         double beam_asymmetry_y(double s, double t);    // Along the y direction

0076         double beam_asymmetry_4pi(double s, double t);  // integrated over decay angles

0077
0078         // Parity asymmetry

0079         double parity_asymmetry(double s, double t);
0080
0081         // ---------------------------------------------------------------------------

0082         // If helicity amplitudes have already been generated for a value of mV, s, t

0083         // store them

0084         double _cached_mX2 = 0., _cached_s = 0., _cached_t = 0.;
0085         std::vector<std::complex<double>> _cached_helicity_amplitude;
0086
0087         void check_cache(double s, double t);
0088
0089         // ---------------------------------------------------------------------------

0090         // nParams error message

0091         int _nParams = 0;
0092         inline void set_nParams(int N){ _nParams = N; };
0093         inline void check_nParams(std::vector<double> params)
0094         {
0095             if (params.size() != _nParams)
0096             {
0097                 std::cout << "\nWarning! Invalid number of parameters (" << params.size() << ") passed to " << _identifier << ".\n";
0098             }
0099         };
0100
0101         // ---------------------------------------------------------------------------

0102         // Each amplitude must supply a function which returns a vector of allowed 2-tuples {J, P}

0103         virtual std::vector<std::array<int,2>> allowedJP() = 0;
0104
0105         // Allowed JP error message

0106         inline void check_JP(std::array<int,2> JP)
0107         {
0108            std::vector<std::array<int,2>> allowed_JP = allowedJP();
0109             if (std::find(allowed_JP.begin(), allowed_JP.end(), JP) == allowed_JP.end())
0110             {
0111                 std::cout << "Error! Amplitude for spin: " << JP[0] << " and parity " << JP[1] << " for " << _identifier << " unavailable.\n";
0112                 exit(0);
0113             }
0114         };
0115     };
0116 };
0117
0118 #endif