LMU/musrfit
5
0
Fork 0
Code Issues 5 Pull Requests Actions Packages Projects Releases 4 Activity

merged root6

Browse Source
This commit is contained in:
suter_a
2019-09-11 12:12:13 +02:00
parent 5953786611 a1924c2756
commit 8b416b45e0
99 changed files with 6983 additions and 19284 deletions
Download Patch File Download Diff File Expand all files Collapse all files

69
src/classes/CMakeLists.txt
View File

@@ -2,35 +2,55 @@
#--- generate necessary dictionaries ------------------------------------------
set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include)
# ROOT requires that the dictonary header files are found at configuration time.
# Hence, target_include_directories cannot be used here because, targets are
# setup only afterwards.
include_directories(${MUSRFIT_INC})
root_generate_dictionary(
PFourierCanvasDict
-I${FFTW3_INCLUDE_DIR} ${MUSRFIT_INC}/PFourierCanvas.h
PFourierCanvas.h
OPTIONS
-I${MUSRFIT_INC}
-I${FFTW3_INCLUDE_DIR}
-inlineInputHeader
LINKDEF ${MUSRFIT_INC}/PFourierCanvasLinkDef.h
MODULE PFourierCanvas
)
root_generate_dictionary(
PMusrCanvasDict
-I${Boost_INCLUDE_DIR} -I${FFTW3_INCLUDE_DIR} ${MUSRFIT_INC}/PMusrCanvas.h
PMusrCanvasDict
PMusrCanvas.h
LINKDEF ${MUSRFIT_INC}/PMusrCanvasLinkDef.h
OPTIONS
-I${Boost_INCLUDE_DIR} -I${FFTW3_INCLUDE_DIR} -I${MUSRFIT_INC}
-inlineInputHeader
MODULE PMusrCanvas
)
root_generate_dictionary(
PMusrT0Dict
-I${Boost_INCLUDE_DIR} -I${FFTW3_INCLUDE_DIR} ${MUSRFIT_INC}/PMusrT0.h
PMusrT0.h
LINKDEF ${MUSRFIT_INC}/PMusrT0LinkDef.h
OPTIONS
-I${Boost_INCLUDE_DIR} -I${FFTW3_INCLUDE_DIR} -I${MUSRFIT_INC}
-inlineInputHeader
MODULE PMusrT0
)
root_generate_dictionary(
PStartupHandlerDict
-I${FFTW3_INCLUDE_DIR} ${MUSRFIT_INC}/PStartupHandler.h
PStartupHandler.h
LINKDEF ${MUSRFIT_INC}/PStartupHandlerLinkDef.h
OPTIONS
-I${FFTW3_INCLUDE_DIR} -I${MUSRFIT_INC}
-inlineInputHeader
MODULE PStartupHandler
)
root_generate_dictionary(
PUserFcnBaseDict
-I${FFTW3_INCLUDE_DIR} ${MUSRFIT_INC}/PUserFcnBase.h
PUserFcnBaseDict
PUserFcnBase.h
LINKDEF ${MUSRFIT_INC}/PUserFcnBaseLinkDef.h
OPTIONS
-inlineInputHeader
-I${FFTW3_INCLUDE_DIR} -I${MUSRFIT_INC}
MODULE PUserFcnBase
)
@@ -39,7 +59,7 @@ set(prefix "${CMAKE_INSTALL_PREFIX}")
set(exec_prefix "\$\{prefix\}")
set(libdir "\$\{exec_prefix\}/lib")
set(includedir "\$\{prefix\}/include")
set(MUSR_VERSION "1.2.0")
set(MUSR_VERSION "1.3.0")
set(MUSR_LIBRARY_NAME "PMusr")
configure_file("PMusr.pc.in" "PMusr.pc" @ONLY)
set(USERFCN_LIBRARY_NAME "PUserFcnBase")
@@ -66,6 +86,7 @@ add_library(PMusr SHARED
PMusrT0Dict.cxx
PPrepFourier.cpp
PRunAsymmetry.cpp
PRunAsymmetryBNMR.cpp
PRunAsymmetryRRF.cpp
PRunBase.cpp
PRunDataHandler.cpp
@@ -78,13 +99,6 @@ add_library(PMusr SHARED
PStartupHandlerDict.cxx
PTheory.cpp
)
#--- set target properties, e.g. version --------------------------------------
set_target_properties(PMusr
PROPERTIES
VERSION ${MUSR_VERSION}
)
#--- make sure that the include directory is found ----------------------------
target_include_directories(
PMusr BEFORE PRIVATE
@@ -106,11 +120,14 @@ add_library(PUserFcnBase SHARED
)
#--- set target properties, e.g. version --------------------------------------
set_target_properties(PMusr
PROPERTIES
VERSION ${MUSR_VERSION}
)
set_target_properties(PUserFcnBase
PROPERTIES
VERSION ${MUSR_VERSION}
)
#--- make sure that the include directory is found ----------------------------
target_include_directories(
PUserFcnBase BEFORE PRIVATE $<BUILD_INTERFACE:${MUSRFIT_INC}>
@@ -174,15 +191,31 @@ install(
)
#--- install headers ----------------------------------------------------------
#--- only install headers needed for classes with dictionary
install(
FILES ${MUSRFIT_INC}/PFourierCanvas.h
FILES ${MUSRFIT_INC}/PFitterFcn.h
${MUSRFIT_INC}/PFitter.h
${MUSRFIT_INC}/PFourierCanvas.h
${MUSRFIT_INC}/PFourier.h
${MUSRFIT_INC}/PFunctionGrammar.h
${MUSRFIT_INC}/PFunction.h
${MUSRFIT_INC}/PFunctionHandler.h
${MUSRFIT_INC}/PMsr2Data.h
${MUSRFIT_INC}/PMsrHandler.h
${MUSRFIT_INC}/PMusrCanvas.h
${MUSRFIT_INC}/PMusr.h
${MUSRFIT_INC}/PMusrT0.h
${MUSRFIT_INC}/PPrepFourier.h
${MUSRFIT_INC}/PRunAsymmetry.h
${MUSRFIT_INC}/PRunAsymmetryRRF.h
${MUSRFIT_INC}/PRunBase.h
${MUSRFIT_INC}/PRunDataHandler.h
${MUSRFIT_INC}/PRunListCollection.h
${MUSRFIT_INC}/PRunMuMinus.h
${MUSRFIT_INC}/PRunNonMusr.h
${MUSRFIT_INC}/PRunSingleHisto.h
${MUSRFIT_INC}/PRunSingleHistoRRF.h
${MUSRFIT_INC}/PStartupHandler.h
${MUSRFIT_INC}/PTheory.h
${MUSRFIT_INC}/PUserFcnBase.h
DESTINATION include
)

4
src/classes/PFitter.cpp
View File

@@ -45,8 +45,8 @@
#include <limits>
#include <cmath>
#include <cstdlib>
#include <cstring>
#include <sys/time.h>
#include "Minuit2/FunctionMinimum.h"
#include "Minuit2/MnContours.h"

2
src/classes/PFitterFcn.cpp
View File

@@ -75,6 +75,7 @@ Double_t PFitterFcn::operator()(const std::vector<Double_t>& par) const
value += fRunListCollection->GetSingleHistoRRFChisq(par);
value += fRunListCollection->GetAsymmetryChisq(par);
value += fRunListCollection->GetAsymmetryRRFChisq(par);
value += fRunListCollection->GetAsymmetryBNMRChisq(par);
value += fRunListCollection->GetMuMinusChisq(par);
value += fRunListCollection->GetNonMusrChisq(par);
} else { // max likelihood
@@ -82,6 +83,7 @@ Double_t PFitterFcn::operator()(const std::vector<Double_t>& par) const
value += fRunListCollection->GetSingleHistoRRFMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryRRFMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryBNMRMaximumLikelihood(par);
value += fRunListCollection->GetMuMinusMaximumLikelihood(par);
value += fRunListCollection->GetNonMusrMaximumLikelihood(par);
}

93
src/classes/PMsrHandler.cpp
View File

@@ -68,8 +68,6 @@ PMsrHandler::PMsrHandler(const Char_t *fileName, PStartupOptions *startupOptions
fFuncHandler = nullptr;
fNoOfMaps = -1;
// check if the file name given is a path-file-name, and if yes, split it into path and file name.
if (fFileName.Contains("/")) {
Int_t idx = -1;
@@ -164,9 +162,9 @@ Int_t PMsrHandler::ReadMsrFile()
continue;
}
// remove leading spaces
line.Remove(TString::kLeading, ' ');
// remove leading spaces
line.Remove(TString::kLeading, ' ');
if (!line.IsWhitespace()) { // if not an empty line, handle it
// check for a msr block
if (line_no == 1) { // title
@@ -642,6 +640,9 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_FITTYPE_MU_MINUS:
fout << std::left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << std::endl ;
break;
case MSR_FITTYPE_BNMR:
fout << std::left << "fittype" << MSR_FITTYPE_BNMR << " (beta-NMR fit)" << std::endl ;
break;
case MSR_FITTYPE_NON_MUSR:
fout << std::left << "fittype" << MSR_FITTYPE_NON_MUSR << " (non muSR fit)" << std::endl ;
break;
@@ -806,6 +807,9 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_FITTYPE_MU_MINUS:
fout << std::left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << std::endl ;
break;
case MSR_FITTYPE_BNMR:
fout << std::left << "fittype" << MSR_FITTYPE_BNMR << " (beta-NMR fit)" << std::endl ;
break;
case MSR_FITTYPE_NON_MUSR:
fout << std::left << "fittype" << MSR_FITTYPE_NON_MUSR << " (non muSR fit)" << std::endl ;
break;
@@ -1169,6 +1173,9 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_PLOT_MU_MINUS:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (mu minus plot)" << std::endl;
break;
case MSR_PLOT_BNMR:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (beta-NMR asymmetry plot)" << std::endl;
break;
case MSR_PLOT_NON_MUSR:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (non muSR plot)" << std::endl;
break;
@@ -1718,6 +1725,9 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, std::map<UInt_t, TString
case MSR_FITTYPE_MU_MINUS:
fout << std::left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << std::endl ;
break;
case MSR_FITTYPE_BNMR:
fout << std::left << "fittype" << MSR_FITTYPE_BNMR << " (beta-NMR fit)" << std::endl ;
break;
case MSR_FITTYPE_NON_MUSR:
fout << std::left << "fittype" << MSR_FITTYPE_NON_MUSR << " (non muSR fit)" << std::endl ;
break;
@@ -1909,6 +1919,9 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, std::map<UInt_t, TString
case MSR_FITTYPE_MU_MINUS:
fout << std::left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << std::endl ;
break;
case MSR_FITTYPE_BNMR:
fout << std::left << "fittype" << MSR_FITTYPE_BNMR << " (beta-NMR fit)" << std::endl ;
break;
case MSR_FITTYPE_NON_MUSR:
fout << std::left << "fittype" << MSR_FITTYPE_NON_MUSR << " (non muSR fit)" << std::endl ;
break;
@@ -2265,6 +2278,9 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, std::map<UInt_t, TString
case MSR_PLOT_MU_MINUS:
fout << "PLOT " << fPlots[i].fPlotType << " (mu minus plot)" << std::endl;
break;
case MSR_PLOT_BNMR:
fout << "PLOT " << fPlots[i].fPlotType << " (beta-NMR asymmetry plot)" << std::endl;
break;
case MSR_PLOT_NON_MUSR:
fout << "PLOT " << fPlots[i].fPlotType << " (non muSR plot)" << std::endl;
break;
@@ -2321,6 +2337,11 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, std::map<UInt_t, TString
fout << "logy" << std::endl;
}
// lifetimecorrection
if (fPlots[i].fLifeTimeCorrection) {
fout << "lifetimecorrection" << std::endl;
}
// rrf_packing
if (fPlots[i].fRRFPacking) {
fout << "rrf_packing " << fPlots[i].fRRFPacking << std::endl;
@@ -2983,6 +3004,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 {
@@ -3212,7 +3234,7 @@ Bool_t PMsrHandler::HandleRunEntry(PMsrLines &lines)
UInt_t addT0Counter = 0;
Int_t ival;
Double_t dval;
Double_t dval;
iter = lines.begin();
while ((iter != lines.end()) && !error) {
@@ -3324,6 +3346,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 {
@@ -4591,6 +4614,7 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
case MSR_PLOT_SINGLE_HISTO: // like: runs 1 5 13
case MSR_PLOT_SINGLE_HISTO_RRF:
case MSR_PLOT_ASYM:
case MSR_PLOT_BNMR:
case MSR_PLOT_ASYM_RRF:
case MSR_PLOT_NON_MUSR:
case MSR_PLOT_MU_MINUS:
@@ -4778,6 +4802,8 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
param.fLogX = true;
} else if (iter1->fLine.Contains("logy", TString::kIgnoreCase)) {
param.fLogY = true;
} else if (iter1->fLine.Contains("lifetimecorrection", TString::kIgnoreCase)) {
param.fLifeTimeCorrection = true;
} else if (iter1->fLine.Contains("view_packing", TString::kIgnoreCase)) {
tokens = iter1->fLine.Tokenize(" \t");
if (!tokens) {
@@ -5008,6 +5034,7 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
std::cerr << std::endl << ">> 2=forward-backward asym,";
std::cerr << std::endl << ">> 3=forward-backward RRF asym,";
std::cerr << std::endl << ">> 4=mu minus single histo,";
std::cerr << std::endl << ">> 5=forward-backward beta-NMR asym,";
std::cerr << std::endl << ">> 8=non muSR.";
std::cerr << std::endl << ">> <run_list> is the list of runs, e.g. runs 1 3";
std::cerr << std::endl << ">> range is optional";
@@ -5865,6 +5892,58 @@ Bool_t PMsrHandler::CheckRunBlockIntegrity()
fRuns[i].SetPacking(1);
}
break;
case PRUN_ASYMMETRY_BNMR:
// check alpha
// if ((fRuns[i].GetAlphaParamNo() == -1) && !fFourierOnly) {
// std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
// std::cerr << std::endl << ">> alpha parameter number missing which is needed for an asymmetry fit.";
// std::cerr << std::endl << ">> Consider to check the manual ;-)" << std::endl;
// return false;
// }
// check that there is a forward parameter number
if (fRuns[i].GetForwardHistoNo() == -1) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
std::cerr << std::endl << ">> forward histogram number not defined. Necessary for asymmetry fits." << std::endl;
return false;
}
// check that there is a backward parameter number
if (fRuns[i].GetBackwardHistoNo() == -1) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
std::cerr << std::endl << ">> backward histogram number not defined. Necessary for asymmetry fits." << std::endl;
return false;
}
// check fit range
if (!fRuns[i].IsFitRangeInBin()) { // fit range given as times in usec
if ((fRuns[i].GetFitRange(0) == PMUSR_UNDEFINED) || (fRuns[i].GetFitRange(1) == PMUSR_UNDEFINED)) {
if ((fGlobal.GetFitRange(0) == PMUSR_UNDEFINED) || (fGlobal.GetFitRange(1) == PMUSR_UNDEFINED)) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
std::cerr << std::endl << ">> Fit range is not defined, also NOT present in the GLOBAL block. Necessary for asymmetry fits." << std::endl;
return false;
}
}
}
// check number of T0's provided
if ((fRuns[i].GetT0BinSize() > 2*fRuns[i].GetForwardHistoNoSize()) &&
(fGlobal.GetT0BinSize() > 2*fRuns[i].GetForwardHistoNoSize())) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
std::cerr << std::endl << ">> Found " << fRuns[i].GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetForwardHistoNoSize() << " in forward. Needs to be fixed." << std::endl;
std::cerr << std::endl << ">> In GLOBAL block: " << fGlobal.GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetForwardHistoNoSize() << ". Needs to be fixed." << std::endl;
return false;
}
if ((fRuns[i].GetT0BinSize() > 2*fRuns[i].GetBackwardHistoNoSize()) &&
(fGlobal.GetT0BinSize() > 2*fRuns[i].GetBackwardHistoNoSize())) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
std::cerr << std::endl << ">> Found " << fRuns[i].GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetBackwardHistoNoSize() << " in backward. Needs to be fixed." << std::endl;
std::cerr << std::endl << ">> In GLOBAL block: " << fGlobal.GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetBackwardHistoNoSize() << ". Needs to be fixed." << std::endl;
return false;
}
// check packing
if ((fRuns[i].GetPacking() == -1) && (fGlobal.GetPacking() == -1)) {
std::cerr << std::endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **WARNING** in RUN block number " << i+1;
std::cerr << std::endl << ">> Packing is neither defined here, nor in the GLOBAL block, will set it to 1." << std::endl;
fRuns[i].SetPacking(1);
}
break;
case PRUN_ASYMMETRY_RRF:
// check alpha
if ((fRuns[i].GetAlphaParamNo() == -1) && !fFourierOnly) {
@@ -6221,7 +6300,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()) {
std::cerr << std::endl << ">> PMsrHandler::CheckHistoGrouping: **ERROR** # of forward histos != # of backward histos.";
std::cerr << std::endl << ">> Run #" << i+1;

13
src/classes/PMusr.cpp
View File

@@ -1878,6 +1878,19 @@ void PMsrRunBlock::SetMapGlobal(UInt_t idx, Int_t ival)
return;
}
//--------------------------------------------------------------------------
// SetEstimatedAlpha (public)
//--------------------------------------------------------------------------
/**
* <p> set the value of estimated alpha at position idx
*
* \param alpha is the estimated value
*/
void PMsrRunBlock::SetEstimatedAlpha(Double_t dval)
{
fAlpha = dval;
}
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// implementation PStringNumberList
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

62
src/classes/PMusrCanvas.cpp
View File

@@ -712,6 +712,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
std::cerr << std::endl << ">> PMusrCanvas::UpdateDataTheoryPad(): **ERROR** couldn't obtain run no " << runNo << " for a beta-NMR asymmetry plot";
std::cerr << std::endl;
return;
}
// handle data
HandleDataSet(i, runNo, data);
break;
case MSR_FITTYPE_ASYM_RRF:
data = fRunList->GetAsymmetryRRF(runNo, PRunListCollection::kRunNo);
if (!data) { // something wrong
@@ -904,7 +916,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);
@@ -1599,6 +1612,7 @@ void PMusrCanvas::ExportData(const Char_t *fileName)
case MSR_PLOT_SINGLE_HISTO:
case MSR_PLOT_SINGLE_HISTO_RRF:
case MSR_PLOT_ASYM:
case MSR_PLOT_BNMR:
case MSR_PLOT_ASYM_RRF:
case MSR_PLOT_MU_MINUS:
if (fDifferenceView) { // difference view plot
@@ -1899,7 +1913,7 @@ void PMusrCanvas::ExportData(const Char_t *fileName)
xmax = fMultiGraphData->GetXaxis()->GetBinCenter(xmaxBin);
// fill ascii dump data
for (UInt_t i=0; i<fNonMusrData.size(); i++) { // go through all the histogramms
for (UInt_t i=0; i<fNonMusrData.size(); i++) { // go through all the graphs
// clean up dump
dump.dataX.clear();
dump.data.clear();
@@ -1917,6 +1931,15 @@ void PMusrCanvas::ExportData(const Char_t *fileName)
}
}
// if anything found keep it
if (dump.dataX.size() > 0)
dumpVector.push_back(dump);
// clean up dump
dump.dataX.clear();
dump.data.clear();
dump.dataErr.clear();
// go through all theory bins
for (Int_t j=0; j<fNonMusrData[i].theory->GetN(); j++) {
// get x and y value
@@ -2758,7 +2781,7 @@ void PMusrCanvas::HandleDataSet(UInt_t plotNo, UInt_t runNo, PRunData *data)
size = data->GetValue()->size();
dataSet.dataRange->SetXRange(start, end); // full possible range
// make sure that for asymmetry the y-range is initialized reasonably
if (fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_ASYM)
if ((fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_ASYM) || (fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_BNMR))
dataSet.dataRange->SetYRange(-0.4, 0.4);
// extract necessary range information
if ((fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fTmin.size() == 0) &&
@@ -2774,6 +2797,7 @@ void PMusrCanvas::HandleDataSet(UInt_t plotNo, UInt_t runNo, PRunData *data)
fXmax = end;
}
if ((fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_ASYM) ||
(fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_BNMR) ||
(fMsrHandler->GetMsrRunList()->at(runNo).IsLifetimeCorrected())) {
fYRangePresent = true;
fYmin = -0.4;
@@ -2822,6 +2846,7 @@ void PMusrCanvas::HandleDataSet(UInt_t plotNo, UInt_t runNo, PRunData *data)
// make sure that for asymmetry the y-range is initialized reasonably
if ((fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_ASYM) ||
(fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fPlotType == MSR_PLOT_BNMR) ||
(fMsrHandler->GetMsrRunList()->at(runNo).IsLifetimeCorrected())) {
dataSet.dataRange->SetYRange(-0.4, 0.4);
}
@@ -4686,7 +4711,13 @@ void PMusrCanvas::PlotData(Bool_t unzoom)
fDataTheoryPad->SetLogy(1);
// set x-axis label
fHistoFrame->GetXaxis()->SetTitle("Time (#mus)");
if (fPlotType == MSR_PLOT_BNMR ) {
// For BNMR/BNQR runs use seconds
fHistoFrame->GetXaxis()->SetTitle("time (s)");
} else {
fHistoFrame->GetXaxis()->SetTitle("time (#mus)");
}
// set y-axis label
TString yAxisTitle;
PMsrRunList *runList = fMsrHandler->GetMsrRunList();
@@ -4708,6 +4739,9 @@ void PMusrCanvas::PlotData(Bool_t unzoom)
case MSR_PLOT_ASYM:
yAxisTitle = "Asymmetry";
break;
case MSR_PLOT_BNMR:
yAxisTitle = "Asymmetry";
break;
case MSR_PLOT_MU_MINUS:
yAxisTitle = "N(t) per bin";
break;
@@ -4804,7 +4838,8 @@ void PMusrCanvas::PlotData(Bool_t unzoom)
}
}
if (fMsrHandler->GetMsrPlotList()->at(fPlotNumber).fLogY) {
dataYmin = 1.0e-4 * dataYmax;
if (dataYmin <= 0.0)
dataYmin = 1.0e-4 * dataYmax;
}
// create fMultiGraphData, and add all data and theory
@@ -4969,7 +5004,12 @@ void PMusrCanvas::PlotDifference(Bool_t unzoom)
fHistoFrame->SetBins(noOfPoints, dataXmin, dataXmax);
// set x-axis label
if (fPlotType == MSR_PLOT_BNMR) {
// For BNMR/BNQR runs use seconds
fHistoFrame->GetXaxis()->SetTitle("time (s)");
} else {
fHistoFrame->GetXaxis()->SetTitle("time (#mus)");
}
// set y-axis label
fHistoFrame->GetYaxis()->SetTitleOffset(1.3);
fHistoFrame->GetYaxis()->SetTitle("data-theory");
@@ -6087,7 +6127,12 @@ void PMusrCanvas::PlotAverage(Bool_t unzoom)
// define x-axis title
TString xAxisTitle("");
if (fCurrentPlotView == PV_DATA) {
if (fPlotType == MSR_PLOT_BNMR) {
// For BNMR/BNQR runs use seconds
xAxisTitle = TString("time (s)");
} else {
xAxisTitle = TString("time (#mus)");
}
} else { // all the Fourier
if (fFourier.fUnits == FOURIER_UNIT_GAUSS) {
xAxisTitle = TString("Field (G)");
@@ -6120,6 +6165,9 @@ void PMusrCanvas::PlotAverage(Bool_t unzoom)
case MSR_PLOT_ASYM:
yAxisTitle = "<asymmetry>";
break;
case MSR_PLOT_BNMR:
yAxisTitle = "<asymmetry>";
break;
case MSR_PLOT_MU_MINUS:
yAxisTitle = "<N(t)> per bin";
break;
@@ -6476,8 +6524,8 @@ UInt_t PMusrCanvas::GetNeededAccuracy(PMsrParamStructure param)
// find last significant digit
for (Int_t i=strlen(str)-1; i>=0; i--) {
if (str[i] != '0') {
if ((i-decimalPoint) < precLimit)
accuracy = i-decimalPoint;
if (((UInt_t)i-decimalPoint) < precLimit)
accuracy = (UInt_t)i-decimalPoint;
else
accuracy = precLimit;
break;

2
src/classes/PRunAsymmetry.cpp
View File

@@ -1776,7 +1776,7 @@ Bool_t PRunAsymmetry::GetProperDataRange(PRawRunData* runData, UInt_t histoNo[2]
}
// 3rd check if end is within proper bounds
if (end[i] < 0) {
std::cerr << std::endl << ">> PRunAsymmetry::GetProperDataRange(): **ERROR** end data bin ("<< end[i] << ") doesn't make any sense!";
std::cerr << std::endl << ">> PRunAsymmetry::GetProperDataRange(): **ERROR** end data bin (" << end[i] << ") doesn't make any sense!";
std::cerr << std::endl;
return false;
}

1775
src/classes/PRunAsymmetryBNMR.cpp Normal file
View File

File diff suppressed because it is too large Load Diff

4
src/classes/PRunAsymmetryRRF.cpp
View File

@@ -596,7 +596,7 @@ Bool_t PRunAsymmetryRRF::PrepareData()
for (UInt_t i=1; i<forwardHistoNo.size(); i++) { // loop over the groupings
for (UInt_t j=0; j<runData->GetDataBin(forwardHistoNo[i])->size(); j++) { // loop over the bin indices
// make sure that the index stays within proper range
if ((j+fT0s[2*i]-fT0s[0] >= 0) && (j+fT0s[2*i]-fT0s[0] < runData->GetDataBin(forwardHistoNo[i])->size())) {
if (((Int_t)j+fT0s[2*i]-fT0s[0] >= 0) && (j+fT0s[2*i]-fT0s[0] < runData->GetDataBin(forwardHistoNo[i])->size())) {
fForward[j] += forward[i][j+(Int_t)fT0s[2*i]-(Int_t)fT0s[0]];
}
}
@@ -770,7 +770,7 @@ Bool_t PRunAsymmetryRRF::SubtractEstimatedBkg()
}
// check if start is within histogram bounds
if ((start[0] >= fForward.size()) || (start[1] >= fBackward.size())) {
if ((start[0] >= fForward.size()) || (start[1] >= fBackward.size())) {
std::cerr << std::endl << ">> PRunAsymmetryRRF::SubtractEstimatedBkg(): **ERROR** background bin values out of bound!";
std::cerr << std::endl << ">> histo lengths (f/b) = (" << fForward.size() << "/" << fBackward.size() << ").";
std::cerr << std::endl << ">> background start (f/b) = (" << start[0] << "/" << start[1] << ").";

3
src/classes/PRunBase.cpp
View File

@@ -34,11 +34,8 @@
#include <TString.h>
#include <TObjArray.h>
#include <TObjString.h>
//#include <TFile.h>
#include <TFolder.h>
#include "TLemRunHeader.h"
#include "PRunBase.h"
//--------------------------------------------------------------------------

7
src/classes/PRunDataHandler.cpp
View File

@@ -3045,8 +3045,9 @@ Bool_t PRunDataHandler::ReadMudFile()
setup += TString(str) + TString("/");
if (TString(str) == "BNQR" || TString(str) == "BNMR") {
std::cerr << "PRunDataHandler::ReadMudFile: **INFORMATION** this run was performed on " << str << std::endl;
std::cerr << "PRunDataHandler::ReadMudFile: **INFORMATION** apply correction to time resolution" << std::endl;
// identified BNMR/BNQR, correct time resolution.
timeResMultiplier = 1.0e15;
timeResMultiplier = 1.0e9;
}
}
success = MUD_getApparatus( fh, str, sizeof(str) );
@@ -3128,8 +3129,8 @@ Bool_t PRunDataHandler::ReadMudFile()
}
}
runData.SetTimeResolution(static_cast<Double_t>(timeResolution) * timeResMultiplier); // s -> ns
// An additional factor of 1e6 needed for bNMR
runData.SetTimeResolution(static_cast<Double_t>(timeResolution) * timeResMultiplier); // s -> ns or s -> ms for bNMR
// Other possibility:
// Check if it is a bNMR run and fix it or check if "timeres" line
// was introduced in the msr file

122
src/classes/PRunListCollection.cpp
View File

@@ -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)
//--------------------------------------------------------------------------
@@ -372,6 +408,18 @@ Double_t PRunListCollection::GetSingleRunChisqExpected(const std::vector<Double_
case PRUN_SINGLE_HISTO:
expectedChisq = fRunSingleHistoList[subIdx]->CalcChiSquareExpected(par);
break;
case PRUN_SINGLE_HISTO_RRF:
expectedChisq = fRunSingleHistoRRFList[subIdx]->CalcChiSquareExpected(par);
break;
case PRUN_ASYMMETRY:
expectedChisq = fRunAsymmetryList[subIdx]->CalcChiSquareExpected(par);
break;
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;
@@ -430,6 +478,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;
@@ -529,6 +580,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)
//--------------------------------------------------------------------------
@@ -710,6 +783,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;
@@ -749,6 +825,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();
@@ -887,6 +966,49 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
return data;
}
//--------------------------------------------------------------------------
// GetAsymmetryBNMR (public)
//--------------------------------------------------------------------------
/**
* <p>Get a processed asymmetry from beta-NMR 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()) {
std::cerr << std::endl << ">> PRunListCollection::GetAsymmetryBNMR(): **ERROR** index = " << index << " out of bounds";
std::cerr << std::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;
}
//--------------------------------------------------------------------------
// GetAsymmetryRRF (public)
//--------------------------------------------------------------------------