include/EvtGenBase/EvtAmpFactory.hh

0001
0002 /***********************************************************************
0003 * Copyright 1998-2020 CERN for the benefit of the EvtGen authors       *
0004 *                                                                      *
0005 * This file is part of EvtGen.                                         *
0006 *                                                                      *
0007 * EvtGen is free software: you can redistribute it and/or modify       *
0008 * it under the terms of the GNU General Public License as published by *
0009 * the Free Software Foundation, either version 3 of the License, or    *
0010 * (at your option) any later version.                                  *
0011 *                                                                      *
0012 * EvtGen is distributed in the hope that it will be useful,            *
0013 * but WITHOUT ANY WARRANTY; without even the implied warranty of       *
0014 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        *
0015 * GNU General Public License for more details.                         *
0016 *                                                                      *
0017 * You should have received a copy of the GNU General Public License    *
0018 * along with EvtGen.  If not, see <https://www.gnu.org/licenses/>.     *
0019 ***********************************************************************/
0020
0021 #ifndef EVT_AMP_FACTORY_HH
0022 #define EVT_AMP_FACTORY_HH
0023
0024 #include "EvtGenBase/EvtAmpPdf.hh"
0025 #include "EvtGenBase/EvtAmplitudeSum.hh"
0026 #include "EvtGenBase/EvtMacros.hh"
0027 #include "EvtGenBase/EvtMultiChannelParser.hh"
0028 #include "EvtGenBase/EvtPdfMax.hh"
0029 #include "EvtGenBase/EvtPdfSum.hh"
0030
0031 #include <stdio.h>
0032 #include <string>
0033 #include <vector>
0034
0035 // Abstract amplitude factory parameterized by a vector of
0036 // strings. Derived classes construct the amplitude, and PDFs for sampling
0037 // points.
0038
0039 template <class T>
0040 class EvtAmpFactory {
0041   public:
0042     EvtAmpFactory() = default;
0043
0044   protected:
0045     EvtAmpFactory( EvtAmpFactory<T>&& ) = default;
0046     EvtAmpFactory( const EvtAmpFactory<T>& other ) :
0047         _amp( other._amp ? other._amp->clone() : nullptr ),
0048         _ampConj( other._ampConj ? other._ampConj->clone() : nullptr ),
0049         _pc( other._pc ? other._pc->clone() : nullptr ),
0050         _names( other._names ),
0051         _dm( other._dm ),
0052         _mixPhase( other._mixPhase ),
0053         _verbose( other._verbose )
0054     {
0055     }
0056
0057   public:
0058     virtual ~EvtAmpFactory() = default;
0059
0060     virtual EvtAmpFactory<T>* clone() const = 0;
0061
0062     virtual void build( const EvtMultiChannelParser& parser, int nItg )
0063     {
0064         _amp = std::make_unique<EvtAmplitudeSum<T>>();
0065         _ampConj = std::make_unique<EvtAmplitudeSum<T>>();
0066         _pc = std::make_unique<EvtPdfSum<T>>();
0067         _dm = parser.dm();
0068         _mixAmpli = parser.mixAmpli();
0069         _mixPhase = parser.mixPhase();
0070
0071         printf( "Amplitude with %d terms\n", parser.getNAmp() );
0072         int i;
0073         for ( i = 0; i < parser.getNAmp(); i++ ) {
0074             std::vector<std::string> v = parser.amp( i );
0075             EvtComplex c = parser.ampCoef( i );
0076             processAmp( c, v );
0077         }
0078
0079         printf( "Conj. amplitude with %d terms\n", parser.getNAmpConj() );
0080         for ( i = 0; i < parser.getNAmpConj(); i++ ) {
0081             std::vector<std::string> v = parser.ampConj( i );
0082             EvtComplex c = parser.ampConjCoef( i );
0083             processAmp( c, v, true );
0084         }
0085
0086         printf( "Calculating pole compensator integrals %d steps\n", nItg );
0087         if ( nItg > 0 )
0088             _pc->getItg( nItg );
0089
0090         printf( "End build\n" );
0091     }
0092
0093     virtual void processAmp( EvtComplex c, std::vector<std::string> v,
0094                              bool conj = false ) = 0;
0095
0096     inline bool isCPModel() const
0097     {
0098         return ( _ampConj->nTerms() > 0 ? true : false );
0099     }
0100     inline double dm() const { return _dm; }
0101     inline double mixAmpli() const { return _mixAmpli; }
0102     inline double mixPhase() const { return _mixPhase; }
0103
0104     void setVerbose() { _verbose = true; }
0105
0106     EvtAmplitudeSum<T>* getAmp() const { return _amp.get(); }
0107     EvtAmplitudeSum<T>* getAmpConj() const { return _ampConj.get(); }
0108     EvtPdfSum<T>* getPC() const { return _pc.get(); }
0109     EvtAmplitude<T>* getAmp( int i ) const { return _amp->getTerm( i ); }
0110     EvtPdf<T>* getPC( int i ) const { return _pc->getPdf( i ); }
0111     const char* compName( int i ) const { return _names[i].c_str(); }
0112
0113     EvtComplex getCoeff( int i ) const { return _amp->c( i ); }
0114
0115     double getTermCoeff( int i ) const { return abs2( _amp->c( i ) ); }
0116     double getTermCoeff( int type, int i, int j ) const
0117     {
0118         switch ( type ) {
0119             case 0:
0120                 return 2 * real( _amp->c( i ) * conj( _amp->c( j ) ) );    //posre
0121             case 1:
0122                 return -2 * real( _amp->c( i ) * conj( _amp->c( j ) ) );    //negre
0123             case 2:
0124                 return -2 * imag( _amp->c( i ) * conj( _amp->c( j ) ) );    //posim
0125             case 3:
0126                 return 2 * imag( _amp->c( i ) * conj( _amp->c( j ) ) );    //negim
0127             default:
0128                 assert( 0 );
0129         }
0130     }
0131
0132   protected:
0133     std::unique_ptr<EvtAmplitudeSum<T>> _amp;        // _owned_ amplitude
0134     std::unique_ptr<EvtAmplitudeSum<T>> _ampConj;    // _owned_ conjugate amplitude
0135     std::unique_ptr<EvtPdfSum<T>> _pc;               // _owned_ pole compensator
0136     std::vector<std::string> _names;    // names of partial amplitudes
0137
0138     double _dm = 0;      // Mass difference for conjugate amplitude
0139     double _mixPhase;    // mixing phase
0140     double _mixAmpli;    // cpv in mixing
0141     bool _verbose = false;
0142 };
0143
0144 #endif