diff --git a/CMakeLists.txt b/CMakeLists.txt
index b6200e8e..e4a7ff9e 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -5,7 +5,7 @@ if (CMAKE_VERSION GREATER_EQUAL 3.12)
cmake_policy(SET CMP0075 NEW)
endif (CMAKE_VERSION GREATER_EQUAL 3.12)
-project(musrfit VERSION 1.7.4 LANGUAGES C CXX)
+project(musrfit VERSION 1.7.5 LANGUAGES C CXX)
#--- musrfit specific options -------------------------------------------------
option(nexus "build optional NeXus support. Needed for ISIS" OFF)
diff --git a/src/classes/PFitter.cpp b/src/classes/PFitter.cpp
index 50633a0d..03142161 100644
--- a/src/classes/PFitter.cpp
+++ b/src/classes/PFitter.cpp
@@ -313,6 +313,9 @@ PFitter::PFitter(PMsrHandler *runInfo, PRunListCollection *runListCollection, Bo
return;
}
+ // create phase bool array
+ GetPhaseParams();
+
// create fit function object
fFitterFcn = new PFitterFcn(runListCollection, fUseChi2);
if (!fFitterFcn) {
@@ -346,6 +349,209 @@ PFitter::~PFitter()
}
}
+//--------------------------------------------------------------------------
+// GetPhaseParams (private)
+//--------------------------------------------------------------------------
+/**
+ *
Checks which parameters are phases. This information is needed to
+ * restrict the phases to the intervall -360 to +360 degrees.
+ */
+void PFitter::GetPhaseParams()
+{
+ fPhase.resize(fRunInfo->GetNoOfParams());
+ for (unsigned int i=0; iGetMsrTheory();
+ TObjArray *tok = nullptr;
+ TObjString *ostr = nullptr;
+ TString str;
+ int pos = -1;
+ for (unsigned int i=0; isize(); i++) {
+ pos = -1;
+ TString line = theo->at(i).fLine;
+ if (line.Contains("TFieldCos") || line.Contains("tf ") ||
+ line.Contains("bessel") || line.Contains("b ") ||
+ line.Contains("skewedGss") || line.Contains("skg ") ||
+ line.Contains("staticNKTF") || line.Contains("snktf ") ||
+ line.Contains("dynamicNKTF") || line.Contains("dnktf ")) { // phase is 1st param
+ pos = 1;
+ }
+ if (line.Contains("internFld") || line.Contains("if ") ||
+ line.Contains("internBsl") || line.Contains("ib ")) { // phase is 2nd param
+ pos = 2;
+ }
+ if (line.Contains("muMinusExpTF") || line.Contains("mmsetf ")) { // phase is 5th param
+ pos = 5;
+ }
+
+ if (pos == -1)
+ continue;
+
+ // extract phase token
+ tok = line.Tokenize(" \t");
+ if (tok == nullptr) {
+ std::cerr << "PFitter::GetPhaseParams(): **ERROR** couldn't tokenize theory line string." << std::endl;
+ return;
+ }
+ if (tok->GetEntries() > pos) {
+ ostr = dynamic_cast(tok->At(pos));
+ str = ostr->GetString();
+ }
+ // clean up
+ delete tok;
+ tok = nullptr;
+
+ // decode phase token. It can be funX, mapX, or a number
+ if (str.Contains("fun")) { // function
+ PIntVector parVec = GetParFromFun(str);
+ for (int i=0; iGetNoOfParams())
+ fPhase[parVec[i]-1] = true;
+ }
+ } else if (str.Contains("map")) { // map
+ PIntVector parVec = GetParFromMap(str);
+ for (int i=0; iGetNoOfParams())
+ fPhase[parVec[i]-1] = true;
+ }
+ } else { // must be a number
+ int idx = str.Atoi();
+ if (idx == 0) { // something went wrong, str is not an integer
+ std::cerr << "PFitter::GetPhaseParams(): **ERROR** str=" << str.View() << " is not an integer!" << std::endl;
+ return;
+ }
+ idx -= 1; // param start at 1, vector at 0
+ if (idx >= fRunInfo->GetNoOfParams()) { // idx is out-of-range
+ std::cerr << "PFitter::GetPhaseParams(): **ERROR** idx=" << idx << " is > #param = " << fRunInfo->GetNoOfParams() << "!" << std::endl;
+ return;
+ }
+ fPhase[idx] = true;
+ }
+ }
+}
+
+//--------------------------------------------------------------------------
+// GetParFromFun (private)
+//--------------------------------------------------------------------------
+/**
+ * Extract from string funX the function number. Base on the function number
+ * the paramter numbers will be collected.
+ *
+ * @param funStr string of the form funX, where X is the function number
+ * @return a vector of all the parameter numbers related to funX
+ */
+PIntVector PFitter::GetParFromFun(const TString funStr)
+{
+ PIntVector parVec;
+
+ PMsrLines *funList = fRunInfo->GetMsrFunctions();
+ TObjArray *tok = nullptr;
+ TObjString *ostr = nullptr;
+ TString str;
+
+ for (int i=0; isize(); i++) {
+ if (funList->at(i).fLine.Contains(funStr)) {
+ // tokenize function string
+ tok = funList->at(i).fLine.Tokenize(" =+-*/");
+ if (tok == nullptr) {
+ std::cerr << "PFitter::GetParFromFun(): **ERROR** couldn't tokenize function string." << std::endl;
+ return parVec;
+ }
+
+ for (int j=1; jGetEntries(); j++) {
+ ostr = dynamic_cast(tok->At(j));
+ str = ostr->GetString();
+ // parse tok for parX
+ if (str.Contains("par")) {
+ // find start idx of par in token
+ Ssiz_t idx = str.Index("par");
+ idx += 3;
+ TString parStr("");
+ do {
+ parStr += str[idx];
+ } while (isdigit(str[idx++]));
+ parVec.push_back(parStr.Atoi());
+ }
+ // parse tok for mapX
+ if (str.Contains("map")) {
+ // find start idx of par in token
+ Ssiz_t idx = str.Index("map");
+ idx += 3;
+ TString mapStr("map");
+ do {
+ mapStr += str[idx];
+ } while (isdigit(str[idx++]));
+ PIntVector mapParVec = GetParFromMap(mapStr);
+ for (int k=0; kExtract from string mapX the map number. Based on the map number the
+ * parameter numbers will be collected.
+ *
+ * @param mapStr string of the form mapX, where X is the map number
+ * @return a vector of all the parameter numbers related to mapX
+ */
+PIntVector PFitter::GetParFromMap(const TString mapStr)
+{
+ PIntVector parVec;
+
+ TString str = mapStr;
+ str.Remove(0,3); // remove map from string
+
+ int idx=str.Atoi();
+ if (idx == 0) {
+ std::cerr << "PFitter::GetParFromMap(): **ERROR** couldn't get propper index from mapX!" << std::endl;
+ return parVec;
+ }
+
+ idx -= 1; // map starts at 1, map vector at 0
+
+ // go through all the runs and collect the parameters from the map vectors
+ PMsrRunList *runList = fRunInfo->GetMsrRunList();
+ if (runList == nullptr) {
+ std::cerr << "PFitter::GetParFromMap(): **ERROR** couldn't get required run list information!" << std::endl;
+ return parVec;
+ }
+
+ PIntVector *map = nullptr;
+ for (int i=0; isize(); i++) {
+ map = runList->at(i).GetMap();
+ if (map == nullptr) {
+ std::cerr << "PFitter::GetParFromMap(): **ERROR** couldn't get required map information (idx=" << i << ")!" << std::endl;
+ parVec.clear();
+ return parVec;
+ }
+ if (idx >= map->size()) {
+ std::cerr << "PFitter::GetParFromMap(): **ERROR** requested map index (idx=" << idx << ") out-of-range (" << map->size() << ")!" << std::endl;
+ parVec.clear();
+ return parVec;
+ }
+ parVec.push_back(map->at(idx));
+ }
+
+ return parVec;
+}
+
//--------------------------------------------------------------------------
// DoFit
//--------------------------------------------------------------------------
@@ -1608,7 +1814,12 @@ Bool_t PFitter::ExecuteMigrad()
// fill run info
for (UInt_t i=0; iSetMsrParamValue(i, min.UserState().Value(i));
+ Double_t dval = min.UserState().Value(i);
+ if (fPhase[i]) {
+ Int_t m = (Int_t)(dval/360.0);
+ dval = dval - m*360.0;
+ }
+ fRunInfo->SetMsrParamValue(i, dval);
fRunInfo->SetMsrParamStep(i, min.UserState().Error(i));
fRunInfo->SetMsrParamPosErrorPresent(i, false);
}
@@ -1684,7 +1895,12 @@ Bool_t PFitter::ExecuteMinimize()
// fill run info
for (UInt_t i=0; iSetMsrParamValue(i, min.UserState().Value(i));
+ Double_t dval = min.UserState().Value(i);
+ if (fPhase[i]) {
+ Int_t m = (Int_t)(dval/360.0);
+ dval = dval - m*360.0;
+ }
+ fRunInfo->SetMsrParamValue(i, dval);
fRunInfo->SetMsrParamStep(i, min.UserState().Error(i));
fRunInfo->SetMsrParamPosErrorPresent(i, false);
}
@@ -2392,7 +2608,12 @@ Bool_t PFitter::ExecuteSimplex()
// fill run info
for (UInt_t i=0; iSetMsrParamValue(i, min.UserState().Value(i));
+ Double_t dval = min.UserState().Value(i);
+ if (fPhase[i]) {
+ Int_t m = (Int_t)(dval/360.0);
+ dval = dval - m*360.0;
+ }
+ fRunInfo->SetMsrParamValue(i, dval);
fRunInfo->SetMsrParamStep(i, min.UserState().Error(i));
fRunInfo->SetMsrParamPosErrorPresent(i, false);
}
diff --git a/src/include/PFitter.h b/src/include/PFitter.h
index 3c513c42..2e0e3e05 100644
--- a/src/include/PFitter.h
+++ b/src/include/PFitter.h
@@ -155,6 +155,13 @@ class PFitter
Bool_t fSectorFlag; ///< sector command present flag
std::vector fSector; ///< stores all chisq/maxLH sector information
+ std::vector fPhase; ///< flag array in which an entry is true if the related parameter value is a phase
+
+ // phase related functions
+ void GetPhaseParams();
+ PIntVector GetParFromFun(const TString funStr);
+ PIntVector GetParFromMap(const TString mapStr);
+
// commands
Bool_t CheckCommands();
Bool_t SetParameters();