EvtMultiChannelParser.cpp

Go to the documentation of this file.

 
/***********************************************************************
* Copyright 1998-2020 CERN for the benefit of the EvtGen authors       *
*                                                                      *
* This file is part of EvtGen.                                         *
*                                                                      *
* EvtGen is free software: you can redistribute it and/or modify       *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 3 of the License, or    *
* (at your option) any later version.                                  *
*                                                                      *
* EvtGen is distributed in the hope that it will be useful,            *
* but WITHOUT ANY WARRANTY; without even the implied warranty of       *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        *
* GNU General Public License for more details.                         *
*                                                                      *
* You should have received a copy of the GNU General Public License    *
* along with EvtGen.  If not, see <https://www.gnu.org/licenses/>.     *
***********************************************************************/
 
#include "EvtGenBase/EvtMultiChannelParser.hh"
 
#include "EvtGenBase/EvtDecayMode.hh"
#include "EvtGenBase/EvtPDL.hh"
#include "EvtGenBase/EvtParser.hh"
 
#include <assert.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string>
#include <vector>
 
using std::string;
using std::vector;
 
EvtDecayMode EvtMultiChannelParser::getDecayMode( const char* file )
{
    // Open file, read tokens
 
    EvtParser parser;
    parser.read( file );
 
    // Seek Decay
 
    int i = 0;
    int N = parser.getNToken();
    while ( i < N ) {
        std::string tok = parser.getToken( i++ );
        if ( tok == std::string( "Decay" ) )
            break;
    }
 
    // Get mother
 
    string mother = string( parser.getToken( i++ ).c_str() );
    std::string bf = parser.getToken( i++ );
 
    vector<string> dauV;
    // Get daughters
 
    while ( 1 ) {
        std::string d = parser.getToken( i++ );
 
        if ( EvtPDL::getStdHep( EvtPDL::getId( d.c_str() ) ) == 0 )
            break;
 
        dauV.push_back( string( d.c_str() ) );
    }
 
    EvtDecayMode mode( mother, dauV );
    printf( "Decay File defines mode %s\n", mode.mode().c_str() );
 
    return mode;
}
 
void EvtMultiChannelParser::parse( const char* file, const char* model )
{
    // Open file, read tokens
 
    EvtParser parser;
    parser.read( file );
 
    // Get parameters (tokens between the model name and ;)
 
    int i = 0;
    int N = parser.getNToken();
 
    // Seek the model name
 
    while ( i < N ) {
        std::string tok = parser.getToken( i++ );
        if ( tok == std::string( model ) )
            break;
    }
    if ( i == N ) {
        printf( "No model %s found in decay file %s", model, file );
        exit( 0 );
    }
 
    // Add all tokens up to a semicolon to vector
 
    std::vector<std::string> v;
    while ( i < N ) {
        std::string tok = parser.getToken( i++ );
        if ( tok == std::string( ";" ) )
            break;
        else
            v.push_back( tok );
    }
 
    if ( i == N ) {
        printf( "No terminating ; found in decay file %s", file );
        assert( 0 );
    }
 
    parse( v );
}
 
void EvtMultiChannelParser::parse( const std::vector<std::string>& v )
{
    // place holder for strtod
    char** tc = nullptr;
 
    // Get PDF maximum or number of points to
    // use in the scan.
 
    if ( v[0] == std::string( "MAXPDF" ) ) {
        m_pdfMax = strtod( v[1].c_str(), tc );
        if ( m_pdfMax <= 0 ) {
            printf( "Bad pdfMax=%f\n", m_pdfMax );
            assert( 0 );
        }
    } else if ( v[0] == std::string( "SCANPDF" ) ) {
        m_nScan = atoi( v[1].c_str() );
    } else {
        printf( "Error parsing decay file\n" );
        assert( 0 );
    }
 
    // Now parse the rest of file for amplitude specifications.
 
    bool conjugate = false;
    size_t i = 2;
    assert( isKeyword( v[2] ) );
 
    while ( i < v.size() ) {
        [[maybe_unused]] size_t i0 = i;
 
        // Switch to conjugate amplitudes after keyword
        if ( v[i] == std::string( "CONJUGATE" ) ) {
            assert( conjugate == false );
            conjugate = true;
            i++;
            m_dm = strtod( v[i++].c_str(), tc );
            m_mixAmpli = strtod( v[i++].c_str(), tc );
            m_mixPhase = strtod( v[i++].c_str(), tc );
        }
 
        if ( i >= v.size() )
            break;
        std::vector<std::string> params;
        EvtComplex c;
        int format;
 
        if ( !conjugate && v[i] == std::string( "AMPLITUDE" ) ) {
            while ( !isKeyword( v[++i] ) )
                params.push_back( v[i] );
            m_amp.push_back( params );
 
            parseComplexCoef( i, v, c, format );
            m_ampCoef.push_back( c );
            m_coefFormat.push_back( format );
            continue;
        } else if ( conjugate && v[i] == std::string( "AMPLITUDE" ) ) {
            while ( !isKeyword( v[++i] ) )
                params.push_back( v[i] );
            m_ampConj.push_back( params );
            parseComplexCoef( i, v, c, format );
            m_ampConjCoef.push_back( c );
            m_coefConjFormat.push_back( format );
            continue;
        } else {
            printf( "Expect keyword, found parameter %s\n", v[i].c_str() );
            assert( 0 );
        }
 
        assert( i > i0 );
    }
 
    printf( "PARSING SUCCESSFUL\n" );
    printf( "%d amplitude terms\n", (int)m_amp.size() );
    printf( "%d conj amplitude terms\n", (int)m_ampConj.size() );
}
 
void EvtMultiChannelParser::parseComplexCoef( size_t& i,
                                              const std::vector<std::string>& v,
                                              EvtComplex& c, int& format )
{
    // place holder for strtod
    char** tc = nullptr;
 
    std::string coefString = v[i++];
    assert( coefString == std::string( "COEFFICIENT" ) );
 
    if ( v[i] == std::string( "POLAR_DEG" ) ) {
        double mag = strtod( v[i + 1].c_str(), tc );
        double phaseRad = strtod( v[i + 2].c_str(), tc ) * EvtConst::pi / 180.0;
        i += 3;
        c = EvtComplex( mag * cos( phaseRad ), mag * sin( phaseRad ) );
        format = POLAR_DEG;
    } else if ( v[i] == std::string( "POLAR_RAD" ) ) {
        double mag = strtod( v[i + 1].c_str(), tc );
        double phaseRad = strtod( v[i + 2].c_str(), tc );
        i += 3;
        c = EvtComplex( mag * cos( phaseRad ), mag * sin( phaseRad ) );
        format = POLAR_RAD;
    } else if ( v[i] == std::string( "CARTESIAN" ) ) {
        double re = strtod( v[i + 1].c_str(), tc );
        double im = strtod( v[i + 2].c_str(), tc );
        i += 3;
        c = EvtComplex( re, im );
        format = CARTESIAN;
    } else {
        printf( "Invalid format %s for complex coefficient\n", v[i].c_str() );
        exit( 0 );
    }
}
 
double EvtMultiChannelParser::parseRealCoef( int& i,
                                             const std::vector<std::string>& v )
{
    // place holder for strtod
    char** tc = nullptr;
    double value = 0;
 
    if ( v[i] == std::string( "COEFFICIENT" ) ) {
        value = strtod( v[i + 1].c_str(), tc );
    } else
        assert( 0 );
 
    i += 2;
 
    assert( value > 0. );
    return value;
}
 
bool EvtMultiChannelParser::isKeyword( const std::string& s )
{
    if ( s == std::string( "AMPLITUDE" ) )
        return true;
    if ( s == std::string( "CONJUGATE" ) )
        return true;
    if ( s == std::string( "COEFFICIENT" ) )
        return true;
    return false;
}