FormFactorBase.cc

Go to the documentation of this file.

///
/// @file  FormFactorBase.cc
/// @brief Hammer base form factor class
///

//**** This file is a part of the HAMMER library
//**** Copyright (C) 2016 - 2020 The HAMMER Collaboration
//**** HAMMER is licensed under version 3 of the GPL; see COPYING for details
//**** Please note the MCnet academic guidelines; see GUIDELINES for details

// -*- C++ -*-
#include "Hammer/FormFactorBase.hh"
#include "Hammer/Math/Utils.hh"
#include "Hammer/Tools/Utils.hh"
#include "Hammer/IndexLabels.hh"
#include "Hammer/Math/MultiDim/VectorContainer.hh"
#include "Hammer/Math/MultiDim/SparseContainer.hh"
#include "Hammer/Math/MultiDim/ScalarContainer.hh"
#include "Hammer/Math/Units.hh"

using namespace std;

namespace Hammer {

    namespace MD = MultiDimensional;

    FormFactorBase::FormFactorBase() :
        _errPrefixGroup{"", ""},
        _FFErrLabel{NONE},
        _FFInfoList{},
        _FFErrNames{},
        _tensorList{},
        _units{1.} {
    }

    FormFactorBase::FormFactorBase(const FormFactorBase& other)
        : ParticleData(other),
          SettingsConsumer(other),
          _errPrefixGroup(other._errPrefixGroup),
          _FFErrLabel(other._FFErrLabel),
          _FFInfoList{other._FFInfoList},
          _FFErrNames(other._FFErrNames),
          _tensorList{},
          _units{other._units} {
        for(auto& elem: other._tensorList) {
            addTensor(Tensor{elem});
        }
    }

    // FormFactorBase& FormFactorBase::operator=(const FormFactorBase& other) {
    //     ParticleData::operator=(other);
    //     SettingsConsumer::operator=(other);
    //     _errPrefixGroup = other._errPrefixGroup;
    //     _FFErrLabel = other._FFErrLabel;
    //     _FFErrNames = other._FFErrNames;
    //     _tensorList.clear();
    //     for (auto& elem : other._tensorList) {
    //         addTensor(Tensor{elem});
    //     }
    //     return *this;
    // }

    // FormFactorBase::~FormFactorBase() {
    // }

    vector<double> FormFactorBase::getErrVectorFromDict(const map<std::string, double>& errDict) const {
        vector<double> result(_FFErrNames.size() + 1);
        result[0] = 1.;
        for(size_t pos = 0ul; pos < _FFErrNames.size(); ++pos) {
            auto it = errDict.find(_FFErrNames[pos]);
            if (it != errDict.end()) {
                result[pos+1] = it->second;
            }
        }
        return result;
    }

    vector<double> FormFactorBase::getErrVectorFromSettings(bool useDefault) const {
        vector<double> result(_FFErrNames.size() + 1);
        result[0] = 1.;
        for (size_t pos = 0ul; pos < _FFErrNames.size(); ++pos) {
            if(!useDefault) {
                result[pos + 1] = *getSetting<double>(getFFErrPrefixGroup().get(), _FFErrNames[pos]);
            }
        }
        return result;
    }

    void FormFactorBase::updateFFErrSettings(const vector<double>& values) {
        this->updateVectorOfSettings(values, _FFErrNames, getFFErrPrefixGroup().get(),
                                     WTerm::COMMON);
    }

    void FormFactorBase::updateFFErrSettings(const map<std::string, double>& values) {
        this->updateVectorOfSettings(values, getFFErrPrefixGroup().get(), WTerm::COMMON);
    }

    void FormFactorBase::updateFFErrTensor(vector<double> values, MD::SharedTensorData& data) const {
        ASSERT(values.size() == _FFErrNames.size()+1);
        if (!data || data->rank() == 0) {
            data = MD::SharedTensorData{MD::makeEmptyVector({static_cast<uint16_t>(values.size())}, {_FFInfoList[_signatureIndex].second}).release()};
        }
        for (IndexType i = 0; i < values.size(); ++i) {
            data->element({i}) = values[i];
        }
    }


    const FFPrefixGroup& FormFactorBase::getFFErrPrefixGroup() const {
        return _FFInfoList[_signatureIndex].first;
    }

    std::pair<FFPrefixGroup, IndexLabel> FormFactorBase::getFFErrInfo() const {
        return _FFInfoList[_signatureIndex];
    }

    void FormFactorBase::init() {
        // initSettings();
        ///@todo anything else?
    }

    Tensor& FormFactorBase::getTensor() {
        return _tensorList[_signatureIndex];
    }

    const Tensor& FormFactorBase::getTensor() const {
        return _tensorList[_signatureIndex];
    }

    void FormFactorBase::setGroup(const string& name) {
        bool reInitSettings = (_errPrefixGroup.group != name && _errPrefixGroup.group.size() > 0);
        _errPrefixGroup.group = name;
        for(auto& elem: _FFInfoList){
            elem.first.group = name;  
        }
        if(reInitSettings && getSettingsHandler() != nullptr) {
            initSettings();
        }
    }

    void FormFactorBase::setPrefix(const string& name) {
        _errPrefixGroup.prefix = name;
    }

    const std::string& FormFactorBase::group() const {
        return _FFInfoList[_signatureIndex].first.group;
    }

    void FormFactorBase::setUnits(const string& name) {
        addSetting<string>("Units", name);
    }
    
    void FormFactorBase::calcUnits() {
        string mcunits;
        auto mc = getSetting<string>("Hammer", "Units");
        if (mc != nullptr){
            mcunits = *mc;
        } else {
            throw Error("Hammer units not found. Mars Climate Orbiter would like a word.");
        }
        Units& units = Units::instance();
        _units = units.getUnitsRescalingToMC( mcunits, *(getSetting<string>("Units")) );
    }

    Tensor FormFactorBase::getFFPSIntegrand(const EvaluationGrid& intPoints) {
        auto newdims = getTensor().dims();
        auto newlabs = getTensor().labels();
        auto oldsize = newdims.size();
        newdims.reserve(2 * oldsize);
        newlabs.reserve(2 * oldsize);
        copy_n(newdims.begin(), oldsize, back_inserter(newdims));
        transform_n(newlabs.begin(), oldsize, back_inserter(newlabs),
                    [](IndexLabel l) -> IndexLabel { return static_cast<IndexLabel>(-l); });
        newlabs.push_back(INTEGRATION_INDEX);
        newdims.push_back(static_cast<IndexType>(intPoints.size()));
        Tensor result{"", MD::makeEmptySparse(newdims, newlabs)};
        for (IndexType i = 0; i < intPoints.size(); ++i) {
            evalAtPSPoint(intPoints[i]);
            Tensor t = getTensor();
            t.outerSquare();
            result.addAt(t, INTEGRATION_INDEX, i);
        }
        return result;
    }

    Log& FormFactorBase::getLog() const {
        return Log::getLog("Hammer.FormFactorBase");
    }

    void FormFactorBase::addTensor(Tensor&& tensor) {
        _tensorList.push_back(move(tensor));
    }
    
    void FormFactorBase::addProcessSignature(PdgId parent, const std::vector<PdgId>& daughters) {
        ParticleData::addProcessSignature(parent, daughters);
        _FFInfoList.push_back(pair<FFPrefixGroup, IndexLabel>{_errPrefixGroup, _FFErrLabel});
    }

} // namespace Hammer