PMusr.cpp

Go to the documentation of this file.

/***************************************************************************

  PMusr.cpp

  Author: Andreas Suter
  e-mail: andreas.suter@psi.ch

***************************************************************************/

/***************************************************************************
 *   Copyright (C) 2007-2023 by Andreas Suter                              *
 *   andreas.suter@psi.ch                                                  *
 *                                                                         *
 *   This program 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 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 *   This program 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 this program; if not, write to the                         *
 *   Free Software Foundation, Inc.,                                       *
 *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
 ***************************************************************************/

#include <cassert>
#include <iostream>

#include <boost/algorithm/string.hpp>
using namespace boost;

#include "TMath.h"

#include "PMusr.h"

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PRunData
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// Constructor
//--------------------------------------------------------------------------
PRunData::PRunData()
{
  fDataTimeStart   = PMUSR_UNDEFINED;
  fDataTimeStep    = PMUSR_UNDEFINED;
  fX.clear();        // only used in non-muSR
  fValue.clear();
  fError.clear();
  fTheoryTimeStart = PMUSR_UNDEFINED;
  fTheoryTimeStep  = PMUSR_UNDEFINED;
  fXTheory.clear();  // only used in non-muSR
  fTheory.clear();
}

//--------------------------------------------------------------------------
// Destructor
//--------------------------------------------------------------------------
PRunData::~PRunData()
{
  fX.clear();        // only used in non-muSR
  fValue.clear();
  fError.clear();
  fXTheory.clear();  // only used in non-muSR
  fTheory.clear();
}

//--------------------------------------------------------------------------
// SetTheoryValue (public)
//--------------------------------------------------------------------------
void PRunData::SetTheoryValue(UInt_t idx, Double_t dval)
{
  if (idx > fTheory.size())
    fTheory.resize(idx+1);

  fTheory[idx] = dval;
}

//--------------------------------------------------------------------------
// ReplaceTheory (public)
//--------------------------------------------------------------------------
void PRunData::ReplaceTheory(const PDoubleVector &theo)
{
  fTheory = theo;
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PNonMusrRawRunData
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// Constructor
//--------------------------------------------------------------------------
PNonMusrRawRunData::PNonMusrRawRunData()
{
  fFromAscii = true;
  fLabels.clear();
  fDataTags.clear();

  for (UInt_t i=0; i<fData.size(); i++)
    fData[i].clear();
  fData.clear();

  for (UInt_t i=0; i<fErrData.size(); i++)
    fErrData[i].clear();
  fErrData.clear();
}

//--------------------------------------------------------------------------
// Destructor
//--------------------------------------------------------------------------
PNonMusrRawRunData::~PNonMusrRawRunData()
{
  fLabels.clear();
  fDataTags.clear();

  for (UInt_t i=0; i<fData.size(); i++)
    fData[i].clear();
  fData.clear();

  for (UInt_t i=0; i<fErrData.size(); i++)
    fErrData[i].clear();
  fErrData.clear();
}

//--------------------------------------------------------------------------
// SetSize (public)
//--------------------------------------------------------------------------
void PNonMusrRawRunData::SetSize(const UInt_t size)
{
  // first clean up
  for (UInt_t i=0; i<fData.size(); i++) {
    fData[i].clear();
  }
  fData.clear();
  for (UInt_t i=0; i<fErrData.size(); i++) {
    fErrData[i].clear();
  }
  fErrData.clear();

  // set size
  fData.resize(size);
  fErrData.resize(size);
}

//--------------------------------------------------------------------------
// SetLabel (public)
//--------------------------------------------------------------------------
void PNonMusrRawRunData::SetLabel(const UInt_t idx, const TString str)
{
  if (idx >= fLabels.size()) {
    std::cerr << std::endl << ">> PNonMusrRawRunData::SetLabel: **WARNING** idx=" << idx << " is out of range [0," << fLabels.size() << "[.";
    std::cerr << std::endl;
    return;
  }
  fLabels[idx] = str;
}

//--------------------------------------------------------------------------
// AppendSubData (public)
//--------------------------------------------------------------------------
void PNonMusrRawRunData::AppendSubData(const UInt_t idx, const Double_t dval)
{
  if (idx >= fData.size()) {
    std::cerr << std::endl << ">> PNonMusrRawRunData::AppendSubData: **WARNING** idx=" << idx << " is out of range [0," << fData.size() << "[.";
    std::cerr << std::endl;
    return;
  }
  fData[idx].push_back(dval);
}

//--------------------------------------------------------------------------
// AppendSubErrData (public)
//--------------------------------------------------------------------------
void PNonMusrRawRunData::AppendSubErrData(const UInt_t idx, const Double_t dval)
{
  if (idx >= fErrData.size()) {
    std::cerr << std::endl << ">> PNonMusrRawRunData::AppendSubErrData: **WARNING** idx=" << idx << " is out of range [0," << fErrData.size() << "[.";
    std::cerr << std::endl;
    return;
  }
  fErrData[idx].push_back(dval);
}


//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PRawRunDataSet
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// Constructor
//--------------------------------------------------------------------------
PRawRunDataSet::PRawRunDataSet()
{
  Clear();
}

//--------------------------------------------------------------------------
// Clear (public)
//--------------------------------------------------------------------------
void PRawRunDataSet::Clear()
{
  fName = TString("n/a");
  fHistoNo = -1;
  fTimeZeroBin = 0.0;
  fTimeZeroBinEstimated = 0.0;
  fFirstGoodBin = 0;
  fLastGoodBin = 0;
  fFirstBkgBin = 0;
  fLastBkgBin = 0;
  fData.clear();
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PRawRunDataVector
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// IsPresent (public)
//--------------------------------------------------------------------------
Bool_t PRawRunDataVector::IsPresent(UInt_t histoNo)
{
  Bool_t found=false;

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      found = true;
      break;
    }
  }

  return found;
}

//--------------------------------------------------------------------------
// GetSet (public)
//--------------------------------------------------------------------------
PRawRunDataSet* PRawRunDataVector::GetSet(UInt_t idx)
{
  PRawRunDataSet *result = nullptr;

  if (idx >= fDataVec.size())
    return result;

  return &fDataVec[idx];
}

//--------------------------------------------------------------------------
// Get (public)
//--------------------------------------------------------------------------
PRawRunDataSet* PRawRunDataVector::Get(UInt_t histoNo)
{
  PRawRunDataSet *result = nullptr;

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      result = &fDataVec[i];
      break;
    }
  }

  return result;
}

//--------------------------------------------------------------------------
// operator[] (public)
//--------------------------------------------------------------------------
PRawRunDataSet* PRawRunDataVector::operator[](UInt_t histoNo)
{
  return Get(histoNo);
}

//--------------------------------------------------------------------------
// GetData (public)
//--------------------------------------------------------------------------
PDoubleVector* PRawRunDataVector::GetData(UInt_t histoNo)
{
  PDoubleVector *result = nullptr;

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      result = fDataVec[i].GetData();
      break;
    }
  }

  return result;
}

//--------------------------------------------------------------------------
// GetT0Bin (public)
//--------------------------------------------------------------------------
Double_t PRawRunDataVector::GetT0Bin(UInt_t histoNo)
{
  Double_t result=-1.0; // undefined

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      result = fDataVec[i].GetTimeZeroBin();
      break;
    }
  }

  return result;
}

//--------------------------------------------------------------------------
// GetT0BinEstimated (public)
//--------------------------------------------------------------------------
Double_t PRawRunDataVector::GetT0BinEstimated(UInt_t histoNo)
{
  Double_t result=PMUSR_UNDEFINED;

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      result = fDataVec[i].GetTimeZeroBinEstimated();
      break;
    }
  }

  return result;
}

//--------------------------------------------------------------------------
// GetBkgBin (public)
//--------------------------------------------------------------------------
PIntPair PRawRunDataVector::GetBkgBin(UInt_t histoNo)
{
  PIntPair bkg(-1,-1);

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      bkg.first  = fDataVec[i].GetFirstBkgBin();
      bkg.second = fDataVec[i].GetLastBkgBin();
      break;
    }
  }

  return bkg;
}

//--------------------------------------------------------------------------
// GetGoodDataBin (public)
//--------------------------------------------------------------------------
PIntPair PRawRunDataVector::GetGoodDataBin(UInt_t histoNo)
{
  PIntPair gdb(-1,-1);

  for (UInt_t i=0; i<fDataVec.size(); i++) {
    if (fDataVec[i].GetHistoNo() == static_cast<Int_t>(histoNo)) {
      gdb.first  = fDataVec[i].GetFirstGoodBin();
      gdb.second = fDataVec[i].GetLastGoodBin();
      break;
    }
  }

  return gdb;
}

//--------------------------------------------------------------------------
// Set (public)
//--------------------------------------------------------------------------
void PRawRunDataVector::Set(PRawRunDataSet dataSet, Int_t idx)
{
  if (idx == -1) { // data set to be appended
    fDataVec.push_back(dataSet);
  } else {
    if (idx >= static_cast<Int_t>(fDataVec.size()))
      fDataVec.resize(idx+1);
    fDataVec[idx] = dataSet;
  }
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PRawRunData
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// Constructor
//--------------------------------------------------------------------------
PRawRunData::PRawRunData()
{
  fVersion = TString("n/a");
  fGenericValidatorURL = TString("n/a");
  fSpecificValidatorURL = TString("n/a");
  fGenerator = TString("n/a");
  fComment = TString("n/a");
  fFileName = TString("n/a");
  fLaboratory = TString("n/a");
  fBeamline = TString("n/a");
  fInstrument = TString("n/a");
  fMuonSource = TString("n/a");
  fMuonSpecies = TString("n/a");
  fMuonBeamMomentum = PMUSR_UNDEFINED;
  fMuonSpinAngle = PMUSR_UNDEFINED;
  fRunName = TString("n/a");
  fRunNumber = -1;
  fRunTitle = TString("n/a");
  fSetup = TString("n/a");
  fStartTime = TString("n/a");
  fStartDate = TString("n/a");
  fStartDateTimeSec = 0;
  fStopTime = TString("n/a");
  fStopDate = TString("n/a");
  fStopDateTimeSec = 0;
  fCryo = TString("n/a");
  fSample = TString("n/a");
  fOrientation = TString("n/a");
  fMagnet = TString("n/a");
  fField = PMUSR_UNDEFINED;
  fEnergy = PMUSR_UNDEFINED;
  fTransport = PMUSR_UNDEFINED;
  fTimeResolution = PMUSR_UNDEFINED;
  fRedGreenOffset.push_back(0);
}

//--------------------------------------------------------------------------
// Destructor
//--------------------------------------------------------------------------
PRawRunData::~PRawRunData()
{
  fTemp.clear();
  fRingAnode.clear();
  fRedGreenOffset.clear();
}

//--------------------------------------------------------------------------
// GetTemperature (public)
//--------------------------------------------------------------------------
const Double_t PRawRunData::GetTemperature(const UInt_t idx)
{
  if (idx >= fTemp.size()) {
    std::cerr << std::endl << ">> PRawRunData::GetTemperature: **WARNING** idx=" << idx << " is out of range [0," << fTemp.size() << "[.";
    std::cerr << std::endl;
    return PMUSR_UNDEFINED;
  }
  return fTemp[idx].first;
}

//--------------------------------------------------------------------------
// GetTempError (public)
//--------------------------------------------------------------------------
const Double_t PRawRunData::GetTempError(const UInt_t idx)
{
  if (idx >= fTemp.size()) {
    std::cerr << std::endl << ">> PRawRunData::GetTempError: **WARNING** idx=" << idx << " is out of range [0," << fTemp.size() << "[.";
    std::cerr << std::endl;
    return PMUSR_UNDEFINED;
  }
  return fTemp[idx].second;
}

//--------------------------------------------------------------------------
// GetRingAnode (public)
//--------------------------------------------------------------------------
const Double_t PRawRunData::GetRingAnode(const UInt_t idx)
{
  if (idx >= fRingAnode.size()) {
    std::cerr << std::endl << ">> PRawRunData::GetRingAnode: **WARNING** idx=" << idx << " is out of range [0," << fRingAnode.size() << "[.";
    std::cerr << std::endl;
    return PMUSR_UNDEFINED;
  }
  return fRingAnode[idx];
}

//--------------------------------------------------------------------------
// GetDataSet (public)
//--------------------------------------------------------------------------
PRawRunDataSet* PRawRunData::GetDataSet(const UInt_t idx, Bool_t wantHistoNo)
{
  if (wantHistoNo)
    return fData[idx];
  else
    return fData.GetSet(idx);
}


//--------------------------------------------------------------------------
// SetRingAnode (public)
//--------------------------------------------------------------------------
void PRawRunData::SetRingAnode(const UInt_t idx, const Double_t dval)
{
  if (idx >= fRingAnode.size())
    fRingAnode.resize(idx+1);
  fRingAnode[idx] = dval;
}

//--------------------------------------------------------------------------
// SetTemperature (public)
//--------------------------------------------------------------------------
void PRawRunData::SetTemperature(const UInt_t idx, const Double_t temp, const Double_t errTemp)
{
  if (idx >= fTemp.size()) {
    fTemp.resize(idx+1);
  }
  fTemp[idx].first  = temp;
  fTemp[idx].second = errTemp;
}

//--------------------------------------------------------------------------
// SetTempError (public)
//--------------------------------------------------------------------------
void PRawRunData::SetTempError(const UInt_t idx, const Double_t errTemp)
{
  if (idx >= fTemp.size()) {
    fTemp.resize(idx+1);
  }
  fTemp[idx].first  = PMUSR_UNDEFINED;
  fTemp[idx].second = errTemp;
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PMsrGlobalBlock
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// PMsrGlobalBlock
//--------------------------------------------------------------------------
PMsrGlobalBlock::PMsrGlobalBlock()
{
  fGlobalPresent = false;
  fRRFFreq = RRF_FREQ_UNDEF; // rotating reference frequency in units given by fRRFUnitTag. Only needed for fittype 1
  fRRFUnitTag = RRF_UNIT_UNDEF; // RRF unit tag. Default: undefined
  fRRFPhase = 0.0;
  fRRFPacking = -1; // undefined RRF packing/rebinning
  fFitType = -1; // undefined fit type
  for (UInt_t i=0; i<4; i++) {
    fDataRange[i] = -1; // undefined data bin range
  }
  fFitRangeInBins = false; // default is that fit range is given in time NOT bins
  fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range
  fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range
  fFitRangeOffset[0] = -1; // undefined start fit range offset
  fFitRangeOffset[1] = -1; // undefined end fit range offset
  fPacking = -1; // undefined packing/rebinning
}

//--------------------------------------------------------------------------
// GetRRFFreq (public)
//--------------------------------------------------------------------------
Double_t PMsrGlobalBlock::GetRRFFreq(const char *unit)
{
  Double_t freq = 0.0;

  // check that the units given make sense
  TString unitStr = unit;
  Int_t unitTag = RRF_UNIT_UNDEF;
  if (!unitStr.CompareTo("MHz", TString::kIgnoreCase))
    unitTag = RRF_UNIT_MHz;
  else if (!unitStr.CompareTo("Mc", TString::kIgnoreCase))
    unitTag = RRF_UNIT_Mcs;
  else if (!unitStr.CompareTo("T", TString::kIgnoreCase))
    unitTag = RRF_UNIT_T;
  else
    return RRF_FREQ_UNDEF;

  // calc the conversion factor
  if (unitTag == fRRFUnitTag)
    freq = fRRFFreq;
  else if ((unitTag == RRF_UNIT_MHz) && (fRRFUnitTag == RRF_UNIT_Mcs))
    freq = fRRFFreq/TMath::TwoPi();
  else if ((unitTag == RRF_UNIT_MHz) && (fRRFUnitTag == RRF_UNIT_T))
    freq = fRRFFreq*1e4*GAMMA_BAR_MUON; // 1e4 need for T -> G since GAMMA_BAR_MUON is given in MHz/G
  else if ((unitTag == RRF_UNIT_Mcs) && (fRRFUnitTag == RRF_UNIT_MHz))
    freq = fRRFFreq*TMath::TwoPi();
  else if ((unitTag == RRF_UNIT_Mcs) && (fRRFUnitTag == RRF_UNIT_T))
    freq = fRRFFreq*1e4*TMath::TwoPi()*GAMMA_BAR_MUON; // 1e4 need for T -> G since GAMMA_BAR_MUON is given in MHz/G
  else if ((unitTag == RRF_UNIT_T) && (fRRFUnitTag == RRF_UNIT_MHz))
    freq = fRRFFreq/GAMMA_BAR_MUON*1e-4; // 1e-4 need for G -> T since GAMMA_BAR_MUON is given in MHz/G
  else if ((unitTag == RRF_UNIT_T) && (fRRFUnitTag == RRF_UNIT_Mcs))
    freq = fRRFFreq/(TMath::TwoPi()*GAMMA_BAR_MUON)*1e-4; // 1e-4 need for G -> T since GAMMA_BAR_MUON is given in MHz/G

  return freq;
}

//--------------------------------------------------------------------------
// SetRRFFreq (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetRRFFreq(Double_t freq, const char *unit)
{
  // check that the units given make sense
  TString unitStr = unit;
  Int_t unitTag = RRF_UNIT_UNDEF;
  if (!unitStr.CompareTo("MHz", TString::kIgnoreCase))
    unitTag = RRF_UNIT_MHz;
  else if (!unitStr.CompareTo("Mc", TString::kIgnoreCase))
    unitTag = RRF_UNIT_Mcs;
  else if (!unitStr.CompareTo("T", TString::kIgnoreCase))
    unitTag = RRF_UNIT_T;
  else {
    std::cerr << std::endl << ">> PMsrGlobalBlock::SetRRFFreq: **ERROR** found undefined RRF unit '" << unit << "'!";
    std::cerr << std::endl << ">>                              Will set RRF frequency to 0.0." << std::endl;
    fRRFFreq = 0.0;
    fRRFUnitTag = RRF_UNIT_UNDEF;
  }

  fRRFFreq = freq;
  fRRFUnitTag = unitTag;
}

//--------------------------------------------------------------------------
// GetRRFUnit (public)
//--------------------------------------------------------------------------
TString PMsrGlobalBlock::GetRRFUnit()
{
  TString unit;

  switch (fRRFUnitTag) {
    case RRF_UNIT_UNDEF:
      unit = TString("??");
      break;
    case RRF_UNIT_kHz:
      unit = TString("kHz");
      break;
    case RRF_UNIT_MHz:
      unit = TString("MHz");
      break;
    case RRF_UNIT_Mcs:
      unit = TString("Mc");
      break;
    case RRF_UNIT_G:
      unit = TString("G");
      break;
    case RRF_UNIT_T:
      unit = TString("T");
      break;
    default:
      unit = TString("??");
      break;
  }

  return unit;
}

//--------------------------------------------------------------------------
// SetRRFPacking (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetRRFPacking(Int_t pack)
{
  if (pack <= 0) {
    std::cerr << std::endl << "PMsrGlobalBlock::SetRRFPacking: **WARNING** found RRF packing <= 0. Likely doesn't make any sense." << std::endl;
    fRRFPacking = -1; // set to undefined
  }

  fRRFPacking = pack;
}

//--------------------------------------------------------------------------
// GetDataRange (public)
//--------------------------------------------------------------------------
Int_t PMsrGlobalBlock::GetDataRange(UInt_t idx)
{
  if (idx >= 4)
    return -1;

  return fDataRange[idx];
}

//--------------------------------------------------------------------------
// SetDataRange (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetDataRange(Int_t ival, Int_t idx)
{
  if (idx >= 4) {
    std::cerr << std::endl << ">> PMsrGlobalBlock::SetDataRange: **WARNING** idx=" << idx << ", only idx=0..3 are sensible.";
    std::cerr << std::endl;
    return;
  }

  fDataRange[idx] = ival;
}

//--------------------------------------------------------------------------
// GetT0Bin (public)
//--------------------------------------------------------------------------
Double_t PMsrGlobalBlock::GetT0Bin(UInt_t idx)
{
  if (idx >= fT0.size())
    return -1;

  return fT0[idx];
}

//--------------------------------------------------------------------------
// SetT0Bin (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetT0Bin(Double_t dval, Int_t idx)
{
  if (idx == -1) {
    fT0.push_back(dval);
    return;
  }

  if (idx >= static_cast<Int_t>(fT0.size()))
    fT0.resize(idx+1);

  fT0[idx] = dval;
}

//--------------------------------------------------------------------------
// GetAddT0BinSize (public)
//--------------------------------------------------------------------------
Int_t PMsrGlobalBlock::GetAddT0BinSize(UInt_t addRunIdx)
{
  if (fAddT0.empty())
    return -1;

  if (addRunIdx >= fAddT0.size())
    return -1;

  return fAddT0[addRunIdx].size();
}

//--------------------------------------------------------------------------
// GetAddT0Bin (public)
//--------------------------------------------------------------------------
Double_t PMsrGlobalBlock::GetAddT0Bin(UInt_t addRunIdx, UInt_t histoIdx)
{
  if (fAddT0.empty())
    return -1.0;

  if (addRunIdx >= fAddT0.size())
    return -1.0;

  if (fAddT0[addRunIdx].empty())
    return -1.0;

  if (histoIdx >= fAddT0[addRunIdx].size())
    return -1.0;

  return fAddT0[addRunIdx][histoIdx];
}

//--------------------------------------------------------------------------
// SetAddT0Bin (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetAddT0Bin(Double_t dval, UInt_t addRunIdx, UInt_t histoNoIdx)
{
  if (addRunIdx >= fAddT0.size())
    fAddT0.resize(addRunIdx+1);

  if (histoNoIdx >= fAddT0[addRunIdx].size())
    fAddT0[addRunIdx].resize(histoNoIdx+1);

  fAddT0[addRunIdx][histoNoIdx] = dval;
}

//--------------------------------------------------------------------------
// GetFitRange (public)
//--------------------------------------------------------------------------
Double_t PMsrGlobalBlock::GetFitRange(UInt_t idx)
{
  if (idx >= 2)
    return PMUSR_UNDEFINED;

  return fFitRange[idx];
}

//--------------------------------------------------------------------------
// SetFitRange (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetFitRange(Double_t dval, UInt_t idx)
{
  if (idx >= 2)
    return;

  fFitRange[idx] = dval;
}

//--------------------------------------------------------------------------
// GetFitRangeOffset (public)
//--------------------------------------------------------------------------
Int_t PMsrGlobalBlock::GetFitRangeOffset(UInt_t idx)
{
  if (idx >= 2)
    return -1;

  return fFitRangeOffset[idx];
}

//--------------------------------------------------------------------------
// SetFitRangeOffset (public)
//--------------------------------------------------------------------------
void PMsrGlobalBlock::SetFitRangeOffset(Int_t ival, UInt_t idx)
{
  if (idx >= 2)
    return;

  fFitRangeOffset[idx] = ival;
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PMsrRunBlock
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// PMsrRunBlock
//--------------------------------------------------------------------------
PMsrRunBlock::PMsrRunBlock()
{
  fFitType = -1; // undefined fit type  
  fAlphaParamNo = -1; // undefined alpha parameter number
  fBetaParamNo = -1; // undefined beta parameter number
  fNormParamNo = -1; // undefined norm parameter number
  fBkgFitParamNo = -1; // undefined background parameter number
  fLifetimeParamNo = -1; // undefined lifetime parameter number
  fLifetimeCorrection = false; // lifetime correction == false by default (used in single histogram musrview)
  for (UInt_t i=0; i<2; i++) {
    fBkgEstimated[i] = PMUSR_UNDEFINED;
    fBkgFix[i] = PMUSR_UNDEFINED;
  }
  for (UInt_t i=0; i<4; i++) {
    fBkgRange[i] = -1; // undefined start background range
    fDataRange[i] = -1; // undefined start data range
  }
  fFitRangeInBins = false; // default is that fit range is given in time NOT bins
  fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range
  fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range
  fFitRangeOffset[0] = -1; // undefined start fit range offset
  fFitRangeOffset[1] = -1; // undefined end fit range offset
  fPacking = -1; // undefined packing
  fXYDataIndex[0] = -1; // undefined x data index (NonMusr)
  fXYDataIndex[1] = -1; // undefined y data index (NonMusr)
  fXYDataLabel[0] = TString(""); // undefined x data label (NonMusr)
  fXYDataLabel[1] = TString(""); // undefined y data label (NonMusr)
}

//--------------------------------------------------------------------------
// ~PMsrRunBlock
//--------------------------------------------------------------------------
PMsrRunBlock::~PMsrRunBlock()
{
  fRunName.clear();
  fBeamline.clear();
  fInstitute.clear();
  fFileFormat.clear();
  fForwardHistoNo.clear();
  fBackwardHistoNo.clear();
  fMap.clear();
  fT0.clear();
  fParGlobal.clear();
  fMapGlobal.clear();
}

//--------------------------------------------------------------------------
// CleanUp (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::CleanUp()
{
  fFitType = -1; // undefined fit type
  fAlphaParamNo = -1; // undefined alpha parameter number
  fBetaParamNo = -1; // undefined beta parameter number
  fNormParamNo = -1; // undefined norm parameter number
  fBkgFitParamNo = -1; // undefined background parameter number
  fLifetimeParamNo = -1; // undefined lifetime parameter number
  fLifetimeCorrection = false; // lifetime correction == false by default (used in single histogram musrview)
  fBkgFix[0] = PMUSR_UNDEFINED; // undefined fixed background for forward
  fBkgFix[1] = PMUSR_UNDEFINED; // undefined fixed background for backward
  for (UInt_t i=0; i<4; i++) {
    fBkgRange[i] = -1; // undefined background range
    fDataRange[i] = -1; // undefined data range
  }
  fFitRangeInBins = false; // default is that fit range is given in time NOT bins
  fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range
  fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range
  fFitRangeOffset[0] = -1; // undefined start fit range offset
  fFitRangeOffset[1] = -1; // undefined end fit range offset
  fPacking = -1; // undefined packing
  fXYDataIndex[0] = -1; // undefined x data index (NonMusr)
  fXYDataIndex[1] = -1; // undefined y data index (NonMusr)
  fXYDataLabel[0] = TString(""); // undefined x data label (NonMusr)
  fXYDataLabel[1] = TString(""); // undefined y data label (NonMusr)

  fRunName.clear();
  fBeamline.clear();
  fInstitute.clear();
  fFileFormat.clear();
  fForwardHistoNo.clear();
  fBackwardHistoNo.clear();
  fMap.clear();
  fT0.clear();
  for (UInt_t i=0; i<fAddT0.size(); i++)
    fAddT0[i].clear();
  fAddT0.clear();

  fParGlobal.clear();
  fMapGlobal.clear();
}

//--------------------------------------------------------------------------
// GetRunName (public)
//--------------------------------------------------------------------------
TString* PMsrRunBlock::GetRunName(UInt_t idx)
{
  if (idx>fRunName.size())
    return nullptr;

  return &fRunName[idx];
}

//--------------------------------------------------------------------------
// SetRunName (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetRunName(TString &str, Int_t idx)
{
  if (idx == -1) {
    fRunName.push_back(str);
    return;
  }

  if (idx >= static_cast<Int_t>(fRunName.size()))
    fRunName.resize(idx+1);

  fRunName[idx] = str;
}

//--------------------------------------------------------------------------
// GetBeamline (public)
//--------------------------------------------------------------------------
TString* PMsrRunBlock::GetBeamline(UInt_t idx)
{
  if (idx>fBeamline.size())
    return nullptr;

  return &fBeamline[idx];
}

//--------------------------------------------------------------------------
// SetBeamline (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetBeamline(TString &str, Int_t idx)
{
  if (idx == -1) {
    fBeamline.push_back(str);
    return;
  }

  if (idx >= static_cast<Int_t>(fBeamline.size()))
    fBeamline.resize(idx+1);

  fBeamline[idx] = str;
}

//--------------------------------------------------------------------------
// GetInstitute (public)
//--------------------------------------------------------------------------
TString* PMsrRunBlock::GetInstitute(UInt_t idx)
{
  if (idx>fInstitute.size())
    return nullptr;

  return &fInstitute[idx];
}

//--------------------------------------------------------------------------
// SetInstitute (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetInstitute(TString &str, Int_t idx)
{
  if (idx == -1) {
    fInstitute.push_back(str);
    return;
  }

  if (idx >= static_cast<Int_t>(fInstitute.size()))
    fInstitute.resize(idx+1);

  fInstitute[idx] = str;
}

//--------------------------------------------------------------------------
// GetFileFormat (public)
//--------------------------------------------------------------------------
TString* PMsrRunBlock::GetFileFormat(UInt_t idx)
{
  if (idx>fFileFormat.size())
    return nullptr;

  return &fFileFormat[idx];
}

//--------------------------------------------------------------------------
// SetFileFormat (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetFileFormat(TString &str, Int_t idx)
{
  if (idx == -1) {
    fFileFormat.push_back(str);
    return;
  }

  if (idx >= static_cast<Int_t>(fFileFormat.size()))
    fFileFormat.resize(idx+1);

  fFileFormat[idx] = str;
}

//--------------------------------------------------------------------------
// GetForwardHistoNo (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetForwardHistoNo(UInt_t idx)
{
  if (fForwardHistoNo.empty())
    return -1;

  if (idx>fForwardHistoNo.size())
    return -1;

  return fForwardHistoNo[idx];
}

//--------------------------------------------------------------------------
// SetForwardHistoNo (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetForwardHistoNo(Int_t histoNo, Int_t idx)
{
  if (idx == -1) { // i.e. append forward histo no
    fForwardHistoNo.push_back(histoNo);
  } else {
    if (idx >= static_cast<Int_t>(fForwardHistoNo.size()))
      fForwardHistoNo.resize(idx+1);
    fForwardHistoNo[idx] = histoNo;
  }
}

//--------------------------------------------------------------------------
// GetBackwardHistoNo (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetBackwardHistoNo(UInt_t idx)
{
  if (fBackwardHistoNo.empty())
    return -1;

  if (idx>fBackwardHistoNo.size())
    return -1;

  return fBackwardHistoNo[idx];
}

//--------------------------------------------------------------------------
// SetBackwardHistoNo (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetBackwardHistoNo(Int_t histoNo, Int_t idx)
{
  if (idx == -1) { // i.e. append forward histo no
    fBackwardHistoNo.push_back(histoNo);
  } else {
    if (idx >= static_cast<Int_t>(fBackwardHistoNo.size()))
      fBackwardHistoNo.resize(idx+1);
    fBackwardHistoNo[idx] = histoNo;
  }
}

//--------------------------------------------------------------------------
// GetMap (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetMap(UInt_t idx)
{
  if (idx>fMap.size())
    return -1;

  return fMap[idx];
}

//--------------------------------------------------------------------------
// SetMap (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetMap(Int_t mapVal, Int_t idx)
{
  if (idx == -1) {
    fMap.push_back(mapVal);
    return;
  }

  if (idx >= static_cast<Int_t>(fMap.size()))
    fMap.resize(idx+1);

  fMap[idx] = mapVal;
}

//--------------------------------------------------------------------------
// GetBkgEstimated (public)
//--------------------------------------------------------------------------
Double_t PMsrRunBlock::GetBkgEstimated(UInt_t idx)
{
  if (idx >= 2)
    return PMUSR_UNDEFINED;

  return fBkgEstimated[idx];
}


//--------------------------------------------------------------------------
// SetBkgEstimated (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetBkgEstimated(Double_t dval, Int_t idx)
{
  if (idx >= 2) {
    std::cerr << std::endl << ">> PMsrRunBlock::SetBkgEstimated: **WARNING** idx=" << idx << ", only idx=0,1 are sensible.";
    std::cerr << std::endl;
    return;
  }

  fBkgEstimated[idx] = dval;
}

//--------------------------------------------------------------------------
// GetBkgFix (public)
//--------------------------------------------------------------------------
Double_t PMsrRunBlock::GetBkgFix(UInt_t idx)
{
  if (idx >= 2)
    return PMUSR_UNDEFINED;

  return fBkgFix[idx];
}

//--------------------------------------------------------------------------
// SetBkgFix (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetBkgFix(Double_t dval, Int_t idx)
{
  if (idx >= 2) {
    std::cerr << std::endl << ">> PMsrRunBlock::SetBkgFix: **WARNING** idx=" << idx << ", only idx=0,1 are sensible.";
    std::cerr << std::endl;
    return;
  }

  fBkgFix[idx] = dval;
}

//--------------------------------------------------------------------------
// GetBkgRange (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetBkgRange(UInt_t idx)
{
  if (idx >= 4) {
    return -1;
  }

  return fBkgRange[idx];
}

//--------------------------------------------------------------------------
// SetBkgRange (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetBkgRange(Int_t ival, Int_t idx)
{
  if (idx >= 4) {
    std::cerr << std::endl << ">> PMsrRunBlock::SetBkgRange: **WARNING** idx=" << idx << ", only idx=0..3 are sensible.";
    std::cerr << std::endl;
    return;
  }

  fBkgRange[idx] = ival;
}


//--------------------------------------------------------------------------
// GetDataRange (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetDataRange(UInt_t idx)
{
  if (idx >= 4) {
    return -1;
  }

  return fDataRange[idx];
}

//--------------------------------------------------------------------------
// SetDataRange (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetDataRange(Int_t ival, Int_t idx)
{
  if (idx >= 4) {
    std::cerr << std::endl << ">> PMsrRunBlock::SetDataRange: **WARNING** idx=" << idx << ", only idx=0..3 are sensible.";
    std::cerr << std::endl;
    return;
  }

  fDataRange[idx] = ival;
}

//--------------------------------------------------------------------------
// GetT0Bin (public)
//--------------------------------------------------------------------------
Double_t PMsrRunBlock::GetT0Bin(UInt_t idx)
{
  if (idx >= fT0.size())
    return -1;

  return fT0[idx];
}

//--------------------------------------------------------------------------
// SetT0Bin (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetT0Bin(Double_t dval, Int_t idx)
{
  if (idx == -1) {
    fT0.push_back(dval);
    return;
  }

  if (idx >= static_cast<Int_t>(fT0.size()))
    fT0.resize(idx+1);

  fT0[idx] = dval;
}

//--------------------------------------------------------------------------
// GetAddT0BinSize (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetAddT0BinSize(UInt_t addRunIdx)
{
  if (fAddT0.empty())
    return -1;

  if (addRunIdx >= fAddT0.size())
    return -1;

  return fAddT0[addRunIdx].size();
}

//--------------------------------------------------------------------------
// GetAddT0Bin (public)
//--------------------------------------------------------------------------
Double_t PMsrRunBlock::GetAddT0Bin(UInt_t addRunIdx, UInt_t histoIdx)
{
  if (fAddT0.empty())
    return -1.0;

  if (addRunIdx >= fAddT0.size())
    return -1.0;

  if (fAddT0[addRunIdx].empty())
    return -1.0;

  if (histoIdx >= fAddT0[addRunIdx].size())
    return -1.0;

  return fAddT0[addRunIdx][histoIdx];
}

//--------------------------------------------------------------------------
// SetAddT0Bin (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetAddT0Bin(Double_t dval, UInt_t addRunIdx, UInt_t histoNoIdx)
{
  if (addRunIdx >= fAddT0.size())
    fAddT0.resize(addRunIdx+1);

  if (histoNoIdx >= fAddT0[addRunIdx].size())
    fAddT0[addRunIdx].resize(histoNoIdx+1);

  fAddT0[addRunIdx][histoNoIdx] = dval;
}

//--------------------------------------------------------------------------
// GetFitRange (public)
//--------------------------------------------------------------------------
Double_t PMsrRunBlock::GetFitRange(UInt_t idx)
{
  if (idx >= 2)
    return PMUSR_UNDEFINED;

  return fFitRange[idx];
}

//--------------------------------------------------------------------------
// SetFitRange (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetFitRange(Double_t dval, UInt_t idx)
{
  if (idx >= 2)
    return;

  fFitRange[idx] = dval;
}

//--------------------------------------------------------------------------
// GetFitRangeOffset (public)
//--------------------------------------------------------------------------
Int_t PMsrRunBlock::GetFitRangeOffset(UInt_t idx)
{
  if (idx >= 2)
    return -1;

  return fFitRangeOffset[idx];
}

//--------------------------------------------------------------------------
// SetFitRangeOffset (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetFitRangeOffset(Int_t ival, UInt_t idx)
{
  if (idx >= 2)
    return;

  fFitRangeOffset[idx] = ival;
}

//--------------------------------------------------------------------------
// SetParGlobal (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetParGlobal(const TString &str, Int_t ival)
{
  fParGlobal[str] = ival; // will either create a new entry or overwrite an old one if the key "str" is present
  return;
}

//--------------------------------------------------------------------------
// SetMapGlobal (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetMapGlobal(UInt_t idx, Int_t ival)
{
  if (fMapGlobal.size() != fMap.size())
    fMapGlobal.resize(fMap.size(), -1);
  if (idx < fMap.size() && fMap[idx] > 0)
    fMapGlobal[idx] = ival;
  // else do nothing at the moment
  return;
}

//--------------------------------------------------------------------------
// SetEstimatedAlpha (public)
//--------------------------------------------------------------------------
void PMsrRunBlock::SetEstimatedAlpha(Double_t dval)
{
  fAlpha = dval;
}

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PStringNumberList
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

//--------------------------------------------------------------------------
// Parse (public)
//--------------------------------------------------------------------------
bool PStringNumberList::Parse(std::string &errorMsg, bool ignoreFirstToken)
{
  bool result=true;
  std::vector<std::string> splitVec;
  int ival;

  // before checking tokens, remove 'forbidden' " - " and " : "
  StripSpaces();

  // split string into space separated tokens
  split(splitVec, fString, is_any_of(" "), token_compress_on);

  unsigned int start=0;
  if (ignoreFirstToken)
    start=1;

  for (unsigned int i=start; i<splitVec.size(); i++) {
    if (splitVec[i].length() != 0) { // ignore empty tokens
      if (splitVec[i].find("-") != std::string::npos) { // check for potential range
        std::vector<std::string> subSplitVec;
        // split potential nS-nE token
        split(subSplitVec, splitVec[i], is_any_of("-"),  token_compress_on);

        int start=-1, end=-1;
        unsigned int count=0;
        for (unsigned int j=0; j<subSplitVec.size(); j++) {
          if (subSplitVec[j].length() != 0) { // ignore empty tokens
            if (!IsNumber(subSplitVec[j])) {
              result = false;
            } else {
              count++;
              if (count == 1)
                start = atoi(subSplitVec[j].c_str());
              else if (count == 2)
                end = atoi(subSplitVec[j].c_str());
              else
                result = false;
            }
          }
        }
        if ((start < 0) || (end < 0)) { // check that there is a vaild start and end
          errorMsg = "**ERROR** start or end of a range is not valid";
          result = false;
        }
        if (result) { // no error, hence check start and end
          if (start > end) {
            int swap = end;
            std::cerr << "**WARNING** start=" << start << " > end=" << end << ", hence I will swap them" << std::endl;
            end = start;
            start = swap;
          }
          for (int j=start; j<=end; j++)
            fList.push_back(j);
        }
      } else if (splitVec[i].find(":") != std::string::npos) { // check for potential sequence
        std::vector<std::string> subSplitVec;
        // split potential rStart:rEnd:rStep token
        split(subSplitVec, splitVec[i], is_any_of(":"),  token_compress_on);

        int start=-1, end=-1, step=-1;
        unsigned int count=0;
        for (unsigned int j=0; j<subSplitVec.size(); j++) {
          if (subSplitVec[j].length() != 0) { // ignore empty tokens
            if (!IsNumber(subSplitVec[j])) {
              result = false;
            } else {
              count++;
              if (count == 1)
                start = atoi(subSplitVec[j].c_str());
              else if (count == 2)
                end = atoi(subSplitVec[j].c_str());
              else if (count == 3)
                step = atoi(subSplitVec[j].c_str());
              else
                result = false;
            }
          }
        }
        if ((start < 0) || (end < 0) || (step < 0)) { // check that there is a vaild start and end
          errorMsg = "**ERROR** start, end, or step of a sequence is not valid";
          result = false;
        }
        if (result) { // no error, hence check start and end
          if (start > end) {
            int swap = end;
            std::cerr << "**WARNING** start=" << start << " > end=" << end << ", hence I will swap them" << std::endl;
            end = start;
            start = swap;
          }
          for (int j=start; j<=end; j+=step)
            fList.push_back(j);
        }
      } else if (IsNumber(splitVec[i])) {
        ival = atoi(splitVec[i].c_str());
        fList.push_back(ival);
      } else {
        errorMsg = "**ERROR** invalid token: " + splitVec[i];
        result = false;
      }
    }
  }

  return result;
}

//--------------------------------------------------------------------------
// StripSpaces (private)
//--------------------------------------------------------------------------
void PStringNumberList::StripSpaces()
{
  std::string str=fString;
  int pos=-1;

  // backward scan
  for (int i=str.size(); i>=0; --i) { // check if first space is found
    if ((str[i] == ' ') && (pos == -1)) {
      pos = i;
    } else if ((str[i] == '-') || (str[i] == ':')) { // check for '-' or ':'
      if (pos != -1) {
        str.erase(i+1, pos-i);
      }
    } else if (str[i] != ' ') { // anything but different than a space leads to a reset of the pos counter
      pos = -1;
    }
  }
  // forward scan
  for (unsigned int i=0; i<str.size(); i++) { // check if first space is found
    if ((str[i] == ' ') && (pos == -1)) {
      pos = i;
    } else if ((str[i] == '-') || (str[i] == ':')) { // check for '-' or ':'
      if (pos != -1) {
        str.erase(pos, i-pos);
        i = pos;
      }
    } else if (str[i] != ' ') { // anything but different than a space leads to a reset of the pos counter
      pos = -1;
    }
  }

  fString = str;
}