Start branch to develop full featured beta-NMR support

src/classes/PMsrHandler.cpp

@@ -632,6 +632,9 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
             case MSR_FITTYPE_MU_MINUS:
               fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << "         (mu minus fit)" << endl ;
               break;
+            case MSR_FITTYPE_BNMR:
+              fout << left << "fittype" << MSR_FITTYPE_BNMR << "         (beta-NMR fit)" << endl ;
+              break;
             case MSR_FITTYPE_NON_MUSR:
               fout << left << "fittype" << MSR_FITTYPE_NON_MUSR << "         (non muSR fit)" << endl ;
               break;
@@ -796,6 +799,9 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
           case MSR_FITTYPE_MU_MINUS:
             fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << "         (mu minus fit)" << endl ;
             break;
+          case MSR_FITTYPE_BNMR:
+            fout << left << "fittype" << MSR_FITTYPE_BNMR << "         (beta-NMR fit)" << endl ;
+            break;
           case MSR_FITTYPE_NON_MUSR:
             fout << left << "fittype" << MSR_FITTYPE_NON_MUSR << "         (non muSR fit)" << endl ;
             break;
@@ -1705,6 +1711,9 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
       case MSR_FITTYPE_MU_MINUS:
         fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << "         (mu minus fit)" << endl ;
         break;
+      case MSR_FITTYPE_BNMR:
+        fout << left << "fittype" << MSR_FITTYPE_BNMR << "         (beta-NMR fit)" << endl ;
+        break;
       case MSR_FITTYPE_NON_MUSR:
         fout << left << "fittype" << MSR_FITTYPE_NON_MUSR << "         (non muSR fit)" << endl ;
         break;
@@ -1896,6 +1905,9 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
       case MSR_FITTYPE_MU_MINUS:
         fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << "         (mu minus fit)" << endl ;
         break;
+      case MSR_FITTYPE_BNMR:
+        fout << left << "fittype" << MSR_FITTYPE_BNMR << "         (beta-NMR fit)" << endl ;
+        break;
       case MSR_FITTYPE_NON_MUSR:
         fout << left << "fittype" << MSR_FITTYPE_NON_MUSR << "         (non muSR fit)" << endl ;
         break;
@@ -2965,6 +2977,7 @@ Bool_t PMsrHandler::HandleGlobalEntry(PMsrLines &lines)
               (fittype == MSR_FITTYPE_ASYM) ||
               (fittype == MSR_FITTYPE_ASYM_RRF) ||
               (fittype == MSR_FITTYPE_MU_MINUS) ||
+              (fittype == MSR_FITTYPE_BNMR) ||
               (fittype == MSR_FITTYPE_NON_MUSR)) {
             global.SetFitType(fittype);
           } else {
@@ -3306,6 +3319,7 @@ Bool_t PMsrHandler::HandleRunEntry(PMsrLines &lines)
               (fittype == MSR_FITTYPE_ASYM) ||
               (fittype == MSR_FITTYPE_ASYM_RRF) ||
               (fittype == MSR_FITTYPE_MU_MINUS) ||
+              (fittype == MSR_FITTYPE_BNMR) ||
               (fittype == MSR_FITTYPE_NON_MUSR)) {
             param.SetFitType(fittype);
           } else {
@@ -5865,7 +5879,7 @@ Bool_t PMsrHandler::CheckHistoGrouping()
   Bool_t result = true;
 
   for (UInt_t i=0; i<fRuns.size(); i++) {
-    if (fRuns[i].GetFitType() == MSR_FITTYPE_ASYM) {
+    if (fRuns[i].GetFitType() == MSR_FITTYPE_ASYM || fRuns[i].GetFitType() == MSR_FITTYPE_BNMR) {
       if (fRuns[i].GetForwardHistoNoSize() != fRuns[i].GetBackwardHistoNoSize()) {
         cerr << endl << ">> PMsrHandler::CheckHistoGrouping: **ERROR** # of forward histos != # of backward histos.";
         cerr << endl << ">> Run #" << i+1;

src/classes/PMusrCanvas.cpp

@@ -742,6 +742,18 @@ void PMusrCanvas::UpdateDataTheoryPad()
         // handle data
         HandleDataSet(i, runNo, data);
         break;
+      case MSR_FITTYPE_BNMR:
+        data = fRunList->GetAsymmetryBNMR(runNo, PRunListCollection::kRunNo);
+        if (!data) { // something wrong
+          fValid = false;
+          // error message
+          cerr << endl << ">> PMusrCanvas::UpdateDataTheoryPad(): **ERROR** couldn't obtain run no " << runNo << " for a beta-NMR asymmetry plot";
+          cerr << endl;
+          return;
+        }
+        // handle data
+        HandleDataSet(i, runNo, data);
+        break;
       case MSR_FITTYPE_NON_MUSR:
         data = fRunList->GetNonMusr(runNo, PRunListCollection::kRunNo);
         if (!data) { // something wrong
@@ -910,7 +922,8 @@ void PMusrCanvas::UpdateInfoPad()
       tstr += grouping;
       tstr += TString(",");
     } else if ((runs[runNo].GetFitType() == MSR_FITTYPE_ASYM) ||
-               (runs[runNo].GetFitType() == MSR_FITTYPE_ASYM_RRF)) {
+               (runs[runNo].GetFitType() == MSR_FITTYPE_ASYM_RRF) ||
+	       (runs[runNo].GetFitType() == MSR_FITTYPE_BNMR)) {
       tstr += TString("h:");
       TString grouping;
       fMsrHandler->GetGroupingString(runNo, "forward", grouping);

src/classes/PRunListCollection.cpp

@@ -76,6 +76,12 @@ PRunListCollection::~PRunListCollection()
   }
   fRunAsymmetryRRFList.clear();
 
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++) {
+    fRunAsymmetryBNMRList[i]->CleanUp();
+    fRunAsymmetryBNMRList[i]->~PRunAsymmetryBNMR();
+  }
+  fRunAsymmetryBNMRList.clear();
+
   for (UInt_t i=0; i<fRunMuMinusList.size(); i++) {
     fRunMuMinusList[i]->CleanUp();
     fRunMuMinusList[i]->~PRunMuMinus();
@@ -133,6 +139,11 @@ Bool_t PRunListCollection::Add(Int_t runNo, EPMusrHandleTag tag)
       if (!fRunAsymmetryRRFList[fRunAsymmetryRRFList.size()-1]->IsValid())
         success = false;
       break;
+    case PRUN_ASYMMETRY_BNMR:
+      fRunAsymmetryBNMRList.push_back(new PRunAsymmetryBNMR(fMsrInfo, fData, runNo, tag));
+      if (!fRunAsymmetryBNMRList[fRunAsymmetryBNMRList.size()-1]->IsValid())
+        success = false;
+      break;
     case PRUN_MU_MINUS:
       fRunMuMinusList.push_back(new PRunMuMinus(fMsrInfo, fData, runNo, tag));
       if (!fRunMuMinusList[fRunMuMinusList.size()-1]->IsValid())
@@ -175,6 +186,8 @@ void PRunListCollection::SetFitRange(const TString fitRange)
     fRunAsymmetryList[i]->SetFitRangeBin(fitRange);
   for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
     fRunAsymmetryRRFList[i]->SetFitRangeBin(fitRange);
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++)
+    fRunAsymmetryBNMRList[i]->SetFitRangeBin(fitRange);
   for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
     fRunMuMinusList[i]->SetFitRangeBin(fitRange);
   for (UInt_t i=0; i<fRunNonMusrList.size(); i++)
@@ -201,6 +214,8 @@ void PRunListCollection::SetFitRange(const PDoublePairVector fitRange)
     fRunAsymmetryList[i]->SetFitRange(fitRange);
   for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
     fRunAsymmetryRRFList[i]->SetFitRange(fitRange);
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++)
+    fRunAsymmetryBNMRList[i]->SetFitRange(fitRange);
   for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
     fRunMuMinusList[i]->SetFitRange(fitRange);
   for (UInt_t i=0; i<fRunNonMusrList.size(); i++)
@@ -291,6 +306,27 @@ Double_t PRunListCollection::GetAsymmetryRRFChisq(const std::vector<Double_t>& p
   return chisq;
 }
 
+//--------------------------------------------------------------------------
+// GetAsymmetryBNMRChisq (public)
+//--------------------------------------------------------------------------
+/**
+ * <p>Calculates chi-square of <em>all</em> asymmetry BNMR runs of a msr-file.
+ *
+ * <b>return:</b>
+ * - chi-square of all asymmetry BNMR runs of the msr-file
+ *
+ * \param par fit parameter vector
+ */
+Double_t PRunListCollection::GetAsymmetryBNMRChisq(const std::vector<Double_t>& par) const
+{
+  Double_t chisq = 0.0;
+
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++)
+    chisq += fRunAsymmetryBNMRList[i]->CalcChiSquare(par);
+
+  return chisq;
+}
+
 //--------------------------------------------------------------------------
 // GetMuMinusChisq (public)
 //--------------------------------------------------------------------------
@@ -380,6 +416,9 @@ Double_t PRunListCollection::GetSingleHistoChisqExpected(const std::vector<Doubl
   case PRUN_ASYMMETRY_RRF:
     expectedChisq = fRunAsymmetryRRFList[subIdx]->CalcChiSquareExpected(par);
     break;
+  case PRUN_ASYMMETRY_BNMR:
+    expectedChisq = fRunAsymmetryBNMRList[subIdx]->CalcChiSquareExpected(par);
+    break;
   case PRUN_MU_MINUS:
     expectedChisq = fRunMuMinusList[subIdx]->CalcChiSquareExpected(par);
     break;
@@ -441,6 +480,9 @@ Double_t PRunListCollection::GetSingleRunChisq(const std::vector<Double_t>& par,
   case PRUN_ASYMMETRY_RRF:
     chisq = fRunAsymmetryRRFList[subIdx]->CalcChiSquare(par);
     break;
+  case PRUN_ASYMMETRY_BNMR:
+    chisq = fRunAsymmetryBNMRList[subIdx]->CalcChiSquare(par);
+    break;
   case PRUN_MU_MINUS:
     chisq = fRunMuMinusList[subIdx]->CalcChiSquare(par);
     break;
@@ -540,6 +582,28 @@ Double_t PRunListCollection::GetAsymmetryRRFMaximumLikelihood(const std::vector<
   return mlh;
 }
 
+//--------------------------------------------------------------------------
+// GetAsymmetryBNMRMaximumLikelihood (public)
+//--------------------------------------------------------------------------
+/**
+ * <p> Since it is not clear yet how to handle asymmetry fits with max likelihood
+ * the chi square will be used!
+ *
+ * <b>return:</b>
+ * - chi-square of all asymmetry BNMR runs of the msr-file
+ *
+ * \param par fit parameter vector
+ */
+Double_t PRunListCollection::GetAsymmetryBNMRMaximumLikelihood(const std::vector<Double_t>& par) const
+{
+  Double_t mlh = 0.0;
+
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++)
+    mlh += fRunAsymmetryBNMRList[i]->CalcChiSquare(par);
+
+  return mlh;
+}
+
 //--------------------------------------------------------------------------
 // GetMuMinusMaximumLikelihood (public)
 //--------------------------------------------------------------------------
@@ -721,6 +785,9 @@ UInt_t PRunListCollection::GetNoOfBinsFitted(const UInt_t idx) const
   case PRUN_ASYMMETRY_RRF:
     result = fRunAsymmetryRRFList[subIdx]->GetNoOfFitBins();
     break;
+  case PRUN_ASYMMETRY_BNMR:
+    result = fRunAsymmetryBNMRList[subIdx]->GetNoOfFitBins();
+    break;
   case PRUN_MU_MINUS:
     result = fRunMuMinusList[subIdx]->GetNoOfFitBins();
     break;
@@ -760,6 +827,9 @@ UInt_t PRunListCollection::GetTotalNoOfBinsFitted() const
   for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
     counts += fRunAsymmetryRRFList[i]->GetNoOfFitBins();
 
+  for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++)
+    counts += fRunAsymmetryBNMRList[i]->GetNoOfFitBins();
+
   for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
     counts += fRunMuMinusList[i]->GetNoOfFitBins();
 
@@ -898,6 +968,49 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
   return data;
 }
 
+//--------------------------------------------------------------------------
+// GetAsymmetryBNMR (public)
+//--------------------------------------------------------------------------
+/**
+ * <p>Get a processed asymmetry data set.
+ *
+ * <b>return:</b>
+ * - pointer to the run data set (processed data) if data set is found
+ * - null pointer otherwise
+ *
+ * \param index msr-file run index
+ * \param tag kIndex -> data at index, kRunNo -> data of given run no
+ */
+PRunData* PRunListCollection::GetAsymmetryBNMR(UInt_t index, EDataSwitch tag)
+{
+  PRunData *data = 0;
+
+  switch (tag) {
+    case kIndex: // called from musrfit when dumping the data
+      if (index > fRunAsymmetryList.size()) {
+        cerr << endl << ">> PRunListCollection::GetAsymmetryBNMR(): **ERROR** index = " << index << " out of bounds";
+        cerr << endl;
+        return 0;
+      }
+
+      fRunAsymmetryList[index]->CalcTheory();
+      data = fRunAsymmetryList[index]->GetData();
+      break;
+    case kRunNo: // called from PMusrCanvas
+      for (UInt_t i=0; i<fRunAsymmetryList.size(); i++) {
+        if (fRunAsymmetryList[i]->GetRunNo() == index) {
+          data = fRunAsymmetryList[i]->GetData();
+          break;
+        }
+      }
+      break;
+    default: // error
+      break;
+  }
+
+  return data;
+}
+
 //--------------------------------------------------------------------------
 // GetAsymmetryRRF (public)
 //--------------------------------------------------------------------------
@@ -941,6 +1054,49 @@ PRunData* PRunListCollection::GetAsymmetryRRF(UInt_t index, EDataSwitch tag)
   return data;
 }
 
+//--------------------------------------------------------------------------
+// GetAsymmetryBNMR (public)
+//--------------------------------------------------------------------------
+/**
+ * <p>Get a processed asymmetry BNMR data set.
+ *
+ * <b>return:</b>
+ * - pointer to the run data set (processed data) if data set is found
+ * - null pointer otherwise
+ *
+ * \param index msr-file run index
+ * \param tag kIndex -> data at index, kRunNo -> data of given run no
+ */
+PRunData* PRunListCollection::GetAsymmetryBNMR(UInt_t index, EDataSwitch tag)
+{
+  PRunData *data = 0;
+
+  switch (tag) {
+    case kIndex: // called from musrfit when dumping the data
+      if (index > fRunAsymmetryBNMRList.size()) {
+        cerr << endl << ">> PRunListCollection::GetAsymmetryBNMR(): **ERROR** index = " << index << " out of bounds";
+        cerr << endl;
+        return 0;
+      }
+
+      fRunAsymmetryBNMRList[index]->CalcTheory();
+      data = fRunAsymmetryBNMRList[index]->GetData();
+      break;
+    case kRunNo: // called from PMusrCanvas
+      for (UInt_t i=0; i<fRunAsymmetryBNMRList.size(); i++) {
+        if (fRunAsymmetryBNMRList[i]->GetRunNo() == index) {
+          data = fRunAsymmetryBNMRList[i]->GetData();
+          break;
+        }
+      }
+      break;
+    default: // error
+      break;
+  }
+
+  return data;
+}
+
 //--------------------------------------------------------------------------
 // GetMuMinus (public)
 //--------------------------------------------------------------------------

src/include/PMusr.h

@@ -57,6 +57,7 @@ using namespace std;
 #define PRUN_ASYMMETRY        2
 #define PRUN_ASYMMETRY_RRF    3
 #define PRUN_MU_MINUS         4
+#define PRUN_ASYMMETRY_BNMR   5
 #define PRUN_NON_MUSR         8
 
 // muon life time in (us), see PRL99, 032001 (2007)
@@ -99,6 +100,7 @@ using namespace std;
 #define MSR_FITTYPE_ASYM             2
 #define MSR_FITTYPE_ASYM_RRF         3
 #define MSR_FITTYPE_MU_MINUS         4
+#define MSR_FITTYPE_BNMR             5
 #define MSR_FITTYPE_NON_MUSR         8
 
 //-------------------------------------------------------------
@@ -108,6 +110,7 @@ using namespace std;
 #define MSR_PLOT_ASYM             2
 #define MSR_PLOT_ASYM_RRF         3
 #define MSR_PLOT_MU_MINUS         4
+#define MSR_PLOT_BNMR             5
 #define MSR_PLOT_NON_MUSR         8
 
 //-------------------------------------------------------------

src/include/PRunListCollection.h

@@ -40,6 +40,7 @@ using namespace std;
 #include "PRunSingleHistoRRF.h"
 #include "PRunAsymmetry.h"
 #include "PRunAsymmetryRRF.h"
+#include "PRunAsymmetryBNMR.h"
 #include "PRunMuMinus.h"
 #include "PRunNonMusr.h"
 
@@ -63,6 +64,7 @@ class PRunListCollection
     virtual Double_t GetSingleHistoRRFChisq(const std::vector<Double_t>& par) const;
     virtual Double_t GetAsymmetryChisq(const std::vector<Double_t>& par) const;
     virtual Double_t GetAsymmetryRRFChisq(const std::vector<Double_t>& par) const;
+    virtual Double_t GetAsymmetryBNMRChisq(const std::vector<Double_t>& par) const;
     virtual Double_t GetMuMinusChisq(const std::vector<Double_t>& par) const;
     virtual Double_t GetNonMusrChisq(const std::vector<Double_t>& par) const;
 
@@ -73,6 +75,7 @@ class PRunListCollection
     virtual Double_t GetSingleHistoRRFMaximumLikelihood(const std::vector<Double_t>& par) const;
     virtual Double_t GetAsymmetryMaximumLikelihood(const std::vector<Double_t>& par) const;
     virtual Double_t GetAsymmetryRRFMaximumLikelihood(const std::vector<Double_t>& par) const;
+    virtual Double_t GetAsymmetryBNMRMaximumLikelihood(const std::vector<Double_t>& par) const;
     virtual Double_t GetMuMinusMaximumLikelihood(const std::vector<Double_t>& par) const;
     virtual Double_t GetNonMusrMaximumLikelihood(const std::vector<Double_t>& par) const;
 
@@ -86,6 +89,7 @@ class PRunListCollection
     virtual UInt_t GetNoOfSingleHistoRRF() const { return fRunSingleHistoRRFList.size(); } ///< returns the number of single histogram RRF data sets present in the msr-file
     virtual UInt_t GetNoOfAsymmetry() const { return fRunAsymmetryList.size(); }     ///< returns the number of asymmetry data sets present in the msr-file
     virtual UInt_t GetNoOfAsymmetryRRF() const { return fRunAsymmetryRRFList.size(); }  ///< returns the number of asymmetry RRF data sets present in the msr-file
+    virtual UInt_t GetNoOfAsymmetryBNMR() const { return fRunAsymmetryBNMRList.size(); }  ///< returns the number of asymmetry BNMR data sets present in the msr-file
     virtual UInt_t GetNoOfMuMinus() const { return fRunMuMinusList.size(); }         ///< returns the number of mu minus data sets present in the msr-file
     virtual UInt_t GetNoOfNonMusr() const { return fRunNonMusrList.size(); }         ///< returns the number of non-muSR data sets present in the msr-file
 
@@ -93,6 +97,7 @@ class PRunListCollection
     virtual PRunData* GetSingleHistoRRF(UInt_t index, EDataSwitch tag=kIndex);
     virtual PRunData* GetAsymmetry(UInt_t index, EDataSwitch tag=kIndex);
     virtual PRunData* GetAsymmetryRRF(UInt_t index, EDataSwitch tag=kIndex);
+    virtual PRunData* GetAsymmetryBNMR(UInt_t index, EDataSwitch tag=kIndex);
     virtual PRunData* GetMuMinus(UInt_t index, EDataSwitch tag=kIndex);
     virtual PRunData* GetNonMusr(UInt_t index, EDataSwitch tag=kIndex);
 
@@ -111,6 +116,7 @@ class PRunListCollection
     vector<PRunSingleHistoRRF*> fRunSingleHistoRRFList; ///< stores all processed single histogram RRF data
     vector<PRunAsymmetry*>      fRunAsymmetryList;      ///< stores all processed asymmetry data
     vector<PRunAsymmetryRRF*>   fRunAsymmetryRRFList;   ///< stores all processed asymmetry RRF data
+    vector<PRunAsymmetryBNMR*>  fRunAsymmetryBNMRList;  ///< stores all processed asymmetry BNMR data
     vector<PRunMuMinus*>        fRunMuMinusList;        ///< stores all processed mu-minus data
     vector<PRunNonMusr*>        fRunNonMusrList;        ///< stores all processed non-muSR data
 };