tkmeasure.cpp

/********************************************************************************
 *   Copyright (c) : Université de Lyon 1, CNRS/IN2P3, UMR5822,                 *
 *                   IP2I, F-69622 Villeurbanne Cedex, France                   *
 *                   Normandie Université, ENSICAEN, UNICAEN, CNRS/IN2P3,       *
 *                   LPC Caen, F-14000 Caen, France                             *
 *   Contibutor(s) :                                                            *
 *      Jérémie Dudouet jeremie.dudouet@cnrs.fr [2020]                          *
 *      Diego Gruyer    diego.gruyer@cnrs.fr    [2020]                          *
 *                                                                              *
 *   Licensed under the MIT License <http://opensource.org/licenses/MIT>.       *   
 *   SPDX-License-Identifier: MIT                                               *
 ********************************************************************************/
 
#include "tkmeasure.h"
 
#include <algorithm>
#include <limits>
 
namespace tkn {
}
 
using namespace tkn;
 
tkmeasure::tkmeasure(double _value, const tkunit_manager::units_keys &_unit, double _err)
{
    fValue.set_value(_value, _unit);
    fError.set_value(_err, _unit);
    fErrorLow.set_unit(_unit);
    fErrorHigh.set_unit(_unit);
}
 
tkmeasure::tkmeasure(double _value, const tkunit_manager::units_keys &_unit, double _errlow, double _errhigh)
{
    fValue.set_value(_value, _unit);
    fError.set_unit(_unit);
    fErrorLow.set_value(_errlow, _unit);
    fErrorHigh.set_value(_errhigh, _unit);
 
    fAsymError = true;
}
 
void tkmeasure::set_info_tag(const tkstring &_tag)
{
    finfo_tag = _tag;
 
    if(_tag.is_empty()) return;
    if(_tag.contains(";ERR=?")) {
        set_info(kUncertain);
        finfo_tag.remove_all(";ERR=?");
    }
    if(_tag.contains(";ERR=CA")) {
        set_info(kTheo);
        finfo_tag.remove_all(";ERR=CA");
    }
    if(_tag.contains(";ERR=AP")) {
        set_info(kAbout);
        finfo_tag.remove_all(";ERR=AP");
    }
    if(_tag.contains(";ERR=SY")) {
        set_info(kSystematic);
        finfo_tag.remove_all(";ERR=SY");
    }
    if(_tag.contains(";VAL=()")) {
        set_info(kTentative);
        finfo_tag.remove_all(";VAL=()");
    }
    if(_tag.contains(";ERR=LT") || _tag.contains(";ERR=GT") || _tag.contains(";ERR=LE") || _tag.contains(";ERR=GE")) {
        set_info(kLimit);
    }
    if(_tag.contains(";ERR=CONV")) {
        set_converted();
        finfo_tag.remove_all(";ERR=CONV");
    }
}
 
void tkmeasure::set(const double &_val, const tkstring &_unit)
{
    if(_unit != "" ) set_unit(_unit);
    fValue.set_value(_val);
}
 
void tkmeasure::set(const double &_val, const tkunit_manager::units_keys &_unit)
{
    set_unit(_unit);
    fValue.set_value(_val,_unit);
}
 
bool tkmeasure::set_unit(const tkunit_manager::units_keys &_unit)
{
    const auto old_unit = fValue.get_unit_key();
    if (old_unit == _unit) return true;
 
    const auto old_type = gunits->get_type(old_unit);
    const auto new_type = gunits->get_type(_unit);
    const bool cross_type = old_type != gunits->kUndef_type && old_type != new_type;
 
    if (!cross_type) {
        if (!fValue.set_unit(_unit)) return false;
        fError.set_unit(_unit);
        fErrorLow.set_unit(_unit);
        fErrorHigh.set_unit(_unit);
        return true;
    }
 
    // Energy/width <-> lifetime conversion is reciprocal, so uncertainty
    // endpoints must be converted independently.  Work on local copies and
    // commit only after every conversion has succeeded.
    tkunit converted_value = fValue;
    if (!converted_value.set_unit(_unit)) return false;
 
    const double error_low = fAsymError ? fErrorLow.get_value() : fError.get_value();
    const double error_high = fAsymError ? fErrorHigh.get_value() : fError.get_value();
    auto reverse_limit = [this, new_type] {
        if (ftype == "half-life" || ftype == "decay width")
            ftype = new_type == gunits->kEnergy ? "decay width" : "half-life";
        const auto tag = get_info_tag();
        if (tag.contains(";ERR=LT")) set_info_tag(tag.copy().replace_all(";ERR=LT", ";ERR=GT"));
        else if (tag.contains(";ERR=LE")) set_info_tag(tag.copy().replace_all(";ERR=LE", ";ERR=GE"));
        else if (tag.contains(";ERR=GT")) set_info_tag(tag.copy().replace_all(";ERR=GT", ";ERR=LT"));
        else if (tag.contains(";ERR=GE")) set_info_tag(tag.copy().replace_all(";ERR=GE", ";ERR=LE"));
    };
 
    if (error_low < 0. || error_high < 0.) {
        fValue = converted_value;
        fError.set_value(-1., _unit);
        fErrorLow.set_value(-1., _unit);
        fErrorHigh.set_value(-1., _unit);
        fAsymError = false;
        reverse_limit();
        return true;
    }
 
    const double value = fValue.get_value();
    tkunit endpoint_low(value - error_low, old_unit);
    tkunit endpoint_high(value + error_high, old_unit);
    if (!endpoint_low.set_unit(_unit) || !endpoint_high.set_unit(_unit)) return false;
 
    const double converted = converted_value.get_value();
    const double converted_min = std::min(endpoint_low.get_value(), endpoint_high.get_value());
    const double converted_max = std::max(endpoint_low.get_value(), endpoint_high.get_value());
    const double converted_error_low = converted - converted_min;
    const double converted_error_high = converted_max - converted;
 
    fValue = converted_value;
    if (converted_error_low == converted_error_high) {
        fError.set_value(converted_error_low, _unit);
        fErrorLow.set_value(-1., _unit);
        fErrorHigh.set_value(-1., _unit);
        fAsymError = false;
    } else {
        fError.set_value(-1., _unit);
        fErrorLow.set_value(converted_error_low, _unit);
        fErrorHigh.set_value(converted_error_high, _unit);
        fAsymError = true;
    }
    reverse_limit();
    return true;
}
 
double tkmeasure::get_value(const tkunit_manager::units_keys &_unit) const
{
    if (get_unit_key() == _unit) return get_value();
    tkmeasure converted = *this;
    return converted.set_unit(_unit) ? converted.get_value() : std::numeric_limits<double>::quiet_NaN();
}
 
double tkmeasure::get_error() const
{
    if(fAsymError) {
        glog << warning_o << "The measure contains assymetric errors, the return value is the largest one" << do_endl;
        return std::max(fErrorLow.get_value(),fErrorHigh.get_value());
    }
    return fError.get_value();
}
 
double tkmeasure::get_error(const tkunit_manager::units_keys &_unit) const
{
    if (get_unit_key() == _unit) return get_error();
    tkmeasure converted = *this;
    return converted.set_unit(_unit) ? converted.get_error() : std::numeric_limits<double>::quiet_NaN();
}
 
double tkmeasure::get_error_low() const
{
    if(!fAsymError) {
        return fError.get_value();
    }
    return fErrorLow.get_value();
}
 
double tkmeasure::get_error_low(const tkunit_manager::units_keys &_unit) const
{
    if (get_unit_key() == _unit) return get_error_low();
    tkmeasure converted = *this;
    return converted.set_unit(_unit) ? converted.get_error_low() : std::numeric_limits<double>::quiet_NaN();
}
 
double tkmeasure::get_error_high() const
{
    if(!fAsymError) {
        return fError.get_value();
    }
    return fErrorHigh.get_value();
}
 
double tkmeasure::get_error_high(const tkunit_manager::units_keys &_unit) const
{
    if (get_unit_key() == _unit) return get_error_high();
    tkmeasure converted = *this;
    return converted.set_unit(_unit) ? converted.get_error_high() : std::numeric_limits<double>::quiet_NaN();
}
 
#ifdef HAS_ROOT
ClassImp(tkmeasure);
#endif