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

Merge remote-tracking branch 'remotes/muonspin/musrfit/master'

Browse Source
This commit is contained in:
salman
2016-01-22 17:31:27 +01:00
parent e5f98215e9 d0da0dd063
commit 3e9bfbc5be
98 changed files with 6894 additions and 4764 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
src/classes/Makefile.am
View File

@@ -15,12 +15,14 @@ h_sources = \
../include/PMusrT0.h \
../include/PPrepFourier.h \
../include/PRunAsymmetry.h \
../include/PRunAsymmetryRRF.h \
../include/PRunBase.h \
../include/PRunDataHandler.h \
../include/PRunListCollection.h \
../include/PRunNonMusr.h \
../include/PRunMuMinus.h \
../include/PRunSingleHisto.h \
../include/PRunSingleHistoRRF.h \
../include/PStartupHandler.h \
../include/PTheory.h
@@ -59,12 +61,14 @@ cpp_sources = \
PMusrT0.cpp \
PPrepFourier.cpp \
PRunAsymmetry.cpp \
PRunAsymmetryRRF.cpp \
PRunBase.cpp \
PRunDataHandler.cpp \
PRunListCollection.cpp \
PRunNonMusr.cpp \
PRunMuMinus.cpp \
PRunSingleHisto.cpp \
PRunSingleHistoRRF.cpp \
PStartupHandler.cpp \
PTheory.cpp

2
src/classes/PFitter.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

6
src/classes/PFitterFcn.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -77,12 +77,16 @@ Double_t PFitterFcn::operator()(const std::vector<Double_t>& par) const
if (fUseChi2) { // chi square
value += fRunListCollection->GetSingleHistoChisq(par);
value += fRunListCollection->GetSingleHistoRRFChisq(par);
value += fRunListCollection->GetAsymmetryChisq(par);
value += fRunListCollection->GetAsymmetryRRFChisq(par);
value += fRunListCollection->GetMuMinusChisq(par);
value += fRunListCollection->GetNonMusrChisq(par);
} else { // max likelihood
value += fRunListCollection->GetSingleHistoMaximumLikelihood(par);
value += fRunListCollection->GetSingleHistoRRFMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryRRFMaximumLikelihood(par);
value += fRunListCollection->GetMuMinusMaximumLikelihood(par);
value += fRunListCollection->GetNonMusrMaximumLikelihood(par);
}

2
src/classes/PFourier.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2015 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PFourierCanvas.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2015 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PFunction.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PFunctionHandler.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PMsr2Data.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2009-2015 by Bastian M. Wojek / Andreas Suter *
* Copyright (C) 2009-2016 by Bastian M. Wojek / Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

353
src/classes/PMsrHandler.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2015 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -311,6 +311,11 @@ Int_t PMsrHandler::ReadMsrFile()
if (!CheckAddRunParameters())
result = PMUSR_MSR_SYNTAX_ERROR;
// check that if RRF settings are present, the RUN block settings do correspond
if (result == PMUSR_SUCCESS)
if (!CheckRRFSettings())
result = PMUSR_MSR_SYNTAX_ERROR;
if (result == PMUSR_SUCCESS) {
CheckMaxLikelihood(); // check if the user wants to use max likelihood with asymmetry/non-muSR fit (which is not implemented)
CheckLegacyLifetimecorrection(); // check if lifetimecorrection is found in RUN blocks, if yes transfer it to PLOT blocks
@@ -348,7 +353,7 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
{
const UInt_t prec = 6; // default output precision for float/doubles
UInt_t neededPrec = 0;
UInt_t neededWidth = 9;
UInt_t neededWidth = 9;
Int_t tag, lineNo = 0, number;
Int_t runNo = -1, addRunNo = 0;
@@ -612,9 +617,15 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_FITTYPE_SINGLE_HISTO:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO << " (single histogram fit)" << endl;
break;
case MSR_FITTYPE_SINGLE_HISTO_RRF:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO_RRF << " (single histogram RRF fit)" << endl;
break;
case MSR_FITTYPE_ASYM:
fout << left << "fittype" << MSR_FITTYPE_ASYM << " (asymmetry fit)" << endl ;
break;
case MSR_FITTYPE_ASYM_RRF:
fout << left << "fittype" << MSR_FITTYPE_ASYM_RRF << " (asymmetry RRF fit)" << endl ;
break;
case MSR_FITTYPE_MU_MINUS:
fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << endl ;
break;
@@ -624,6 +635,27 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
default:
break;
}
} else if (sstr.BeginsWith("rrf_freq", TString::kIgnoreCase) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << left << "rrf_freq ";
fout.width(8);
neededPrec = LastSignificant(fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data()),10);
fout.precision(neededPrec);
fout << left << fixed << fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data());
fout << " " << fGlobal.GetRRFUnit();
fout << endl;
} else if (sstr.BeginsWith("rrf_phase", TString::kIgnoreCase) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << "rrf_phase ";
fout.width(8);
fout << left << fGlobal.GetRRFPhase();
fout << endl;
} else if (sstr.BeginsWith("rrf_packing", TString::kIgnoreCase) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << "rrf_packing ";
fout.width(8);
fout << left << fGlobal.GetRRFPacking();
fout << endl;
} else if (sstr.BeginsWith("data")) {
fout.width(16);
fout << left << "data";
@@ -749,9 +781,15 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_FITTYPE_SINGLE_HISTO:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO << " (single histogram fit)" << endl;
break;
case MSR_FITTYPE_SINGLE_HISTO_RRF:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO_RRF << " (single histogram RRF fit)" << endl;
break;
case MSR_FITTYPE_ASYM:
fout << left << "fittype" << MSR_FITTYPE_ASYM << " (asymmetry fit)" << endl ;
break;
case MSR_FITTYPE_ASYM_RRF:
fout << left << "fittype" << MSR_FITTYPE_ASYM_RRF << " (asymmetry RRF fit)" << endl ;
break;
case MSR_FITTYPE_MU_MINUS:
fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << endl ;
break;
@@ -1107,9 +1145,15 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_PLOT_SINGLE_HISTO:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (single histo plot)" << endl;
break;
case MSR_PLOT_SINGLE_HISTO_RRF:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (single histo RRF plot)" << endl;
break;
case MSR_PLOT_ASYM:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (asymmetry plot)" << endl;
break;
case MSR_PLOT_ASYM_RRF:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (asymmetry RRF plot)" << endl;
break;
case MSR_PLOT_MU_MINUS:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (mu minus plot)" << endl;
break;
@@ -1604,9 +1648,15 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
case MSR_FITTYPE_SINGLE_HISTO:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO << " (single histogram fit)" << endl;
break;
case MSR_FITTYPE_SINGLE_HISTO_RRF:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO_RRF << " (single histogram RRF fit)" << endl;
break;
case MSR_FITTYPE_ASYM:
fout << left << "fittype" << MSR_FITTYPE_ASYM << " (asymmetry fit)" << endl ;
break;
case MSR_FITTYPE_ASYM_RRF:
fout << left << "fittype" << MSR_FITTYPE_ASYM_RRF << " (asymmetry RRF fit)" << endl ;
break;
case MSR_FITTYPE_MU_MINUS:
fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << endl ;
break;
@@ -1618,6 +1668,30 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
}
}
// RRF related stuff
if ((fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data()) > 0.0) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << left << "rrf_freq ";
fout.width(8);
fout << left << fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data());
fout << " " << fGlobal.GetRRFUnit();
fout << endl;
}
if ((fGlobal.GetRRFPhase() != 0.0) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << "rrf_phase ";
fout.width(8);
fout << left << fGlobal.GetRRFPhase();
fout << endl;
}
if ((fGlobal.GetRRFPacking() != -1) && (fGlobal.GetFitType() == MSR_FITTYPE_SINGLE_HISTO_RRF)) {
fout.width(16);
fout << "rrf_packing ";
fout.width(8);
fout << left << fGlobal.GetRRFPacking();
fout << endl;
}
// data range
if ((fGlobal.GetDataRange(0) != -1) || (fGlobal.GetDataRange(1) != -1) || (fGlobal.GetDataRange(2) != -1) || (fGlobal.GetDataRange(3) != -1)) {
fout.width(16);
@@ -1765,9 +1839,15 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
case MSR_FITTYPE_SINGLE_HISTO:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO << " (single histogram fit)" << endl;
break;
case MSR_FITTYPE_SINGLE_HISTO_RRF:
fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO_RRF << " (single histogram RRF fit)" << endl;
break;
case MSR_FITTYPE_ASYM:
fout << left << "fittype" << MSR_FITTYPE_ASYM << " (asymmetry fit)" << endl ;
break;
case MSR_FITTYPE_ASYM_RRF:
fout << left << "fittype" << MSR_FITTYPE_ASYM_RRF << " (asymmetry RRF fit)" << endl ;
break;
case MSR_FITTYPE_MU_MINUS:
fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << endl ;
break;
@@ -2108,9 +2188,15 @@ Int_t PMsrHandler::WriteMsrFile(const Char_t *filename, map<UInt_t, TString> *co
case MSR_PLOT_SINGLE_HISTO:
fout << "PLOT " << fPlots[i].fPlotType << " (single histo plot)" << endl;
break;
case MSR_PLOT_SINGLE_HISTO_RRF:
fout << "PLOT " << fPlots[i].fPlotType << " (single histo RRF plot)" << endl;
break;
case MSR_PLOT_ASYM:
fout << "PLOT " << fPlots[i].fPlotType << " (asymmetry plot)" << endl;
break;
case MSR_PLOT_ASYM_RRF:
fout << "PLOT " << fPlots[i].fPlotType << " (asymmetry RRF plot)" << endl;
break;
case MSR_PLOT_MU_MINUS:
fout << "PLOT " << fPlots[i].fPlotType << " (mu minus plot)" << endl;
break;
@@ -2796,8 +2882,8 @@ Bool_t PMsrHandler::HandleGlobalEntry(PMsrLines &lines)
TString str;
TObjArray *tokens = 0;
TObjString *ostr = 0;
Int_t ival;
Double_t dval;
Int_t ival = 0;
Double_t dval = 0.0;
UInt_t addT0Counter = 0;
// since this routine is called, a GLOBAL block is present
@@ -2828,7 +2914,9 @@ Bool_t PMsrHandler::HandleGlobalEntry(PMsrLines &lines)
if (str.IsDigit()) {
Int_t fittype = str.Atoi();
if ((fittype == MSR_FITTYPE_SINGLE_HISTO) ||
(fittype == MSR_FITTYPE_SINGLE_HISTO_RRF) ||
(fittype == MSR_FITTYPE_ASYM) ||
(fittype == MSR_FITTYPE_ASYM_RRF) ||
(fittype == MSR_FITTYPE_MU_MINUS) ||
(fittype == MSR_FITTYPE_NON_MUSR)) {
global.SetFitType(fittype);
@@ -2839,6 +2927,55 @@ Bool_t PMsrHandler::HandleGlobalEntry(PMsrLines &lines)
error = true;
}
}
} else if (iter->fLine.BeginsWith("rrf_freq", TString::kIgnoreCase)) {
if (tokens->GetEntries() < 3) {
error = true;
} else {
ostr = dynamic_cast<TObjString*>(tokens->At(1));
str = ostr->GetString();
if (str.IsFloat()) {
dval = str.Atof();
if (dval <= 0.0)
error = true;
}
if (!error) {
ostr = dynamic_cast<TObjString*>(tokens->At(2));
str = ostr->GetString();
global.SetRRFFreq(dval, str.Data());
if (global.GetRRFFreq(str.Data()) == RRF_FREQ_UNDEF)
error = true;
}
}
} else if (iter->fLine.BeginsWith("rrf_packing", TString::kIgnoreCase)) {
if (tokens->GetEntries() < 2) {
error = true;
} else {
ostr = dynamic_cast<TObjString*>(tokens->At(1));
str = ostr->GetString();
if (str.IsDigit()) {
ival = str.Atoi();
if (ival > 0) {
global.SetRRFPacking(ival);
} else {
error = true;
}
} else {
error = true;
}
}
} else if (iter->fLine.BeginsWith("rrf_phase", TString::kIgnoreCase)) {
if (tokens->GetEntries() < 2) {
error = true;
} else {
ostr = dynamic_cast<TObjString*>(tokens->At(1));
str = ostr->GetString();
if (str.IsFloat()) {
dval = str.Atof();
global.SetRRFPhase(dval);
} else {
error = true;
}
}
} else if (iter->fLine.BeginsWith("data", TString::kIgnoreCase)) { // data
if (tokens->GetEntries() < 3) {
error = true;
@@ -3113,7 +3250,9 @@ Bool_t PMsrHandler::HandleRunEntry(PMsrLines &lines)
if (str.IsDigit()) {
Int_t fittype = str.Atoi();
if ((fittype == MSR_FITTYPE_SINGLE_HISTO) ||
(fittype == MSR_FITTYPE_SINGLE_HISTO_RRF) ||
(fittype == MSR_FITTYPE_ASYM) ||
(fittype == MSR_FITTYPE_ASYM_RRF) ||
(fittype == MSR_FITTYPE_MU_MINUS) ||
(fittype == MSR_FITTYPE_NON_MUSR)) {
param.SetFitType(fittype);
@@ -4027,7 +4166,9 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
error = true;
break;
case MSR_PLOT_SINGLE_HISTO: // like: runs 1 5 13
case MSR_PLOT_SINGLE_HISTO_RRF:
case MSR_PLOT_ASYM:
case MSR_PLOT_ASYM_RRF:
case MSR_PLOT_NON_MUSR:
case MSR_PLOT_MU_MINUS:
rl = new PStringNumberList(iter1->fLine.Data());
@@ -4439,8 +4580,10 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
cerr << endl << ">> [use_fit_ranges [ymin ymax]]";
cerr << endl << ">> [view_packing n]";
cerr << endl;
cerr << endl << ">> where <plot_type> is: 0=single histo asym,";
cerr << endl << ">> where <plot_type> is: 0=single histo,";
cerr << endl << ">> 1=RRF single histo,";
cerr << endl << ">> 2=forward-backward asym,";
cerr << endl << ">> 3=forward-backward RRF asym,";
cerr << endl << ">> 4=mu minus single histo,";
cerr << endl << ">> 8=non muSR.";
cerr << endl << ">> <run_list> is the list of runs, e.g. runs 1 3";
@@ -5193,6 +5336,53 @@ Bool_t PMsrHandler::CheckRunBlockIntegrity()
fRuns[i].SetPacking(1);
}
break;
case PRUN_SINGLE_HISTO_RRF:
// check that there is a forward parameter number
if (fRuns[i].GetForwardHistoNo() == -1) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> forward parameter number not defined. Necessary for single histogram RRF fits." << endl;
return false;
}
if ((fRuns[i].GetNormParamNo() > static_cast<Int_t>(fParam.size())) && !fFourierOnly) {
// check if forward histogram number is a function
if (fRuns[i].GetNormParamNo() - MSR_PARAM_FUN_OFFSET > static_cast<Int_t>(fParam.size())) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> forward histogram number " << fRuns[i].GetNormParamNo() << " is larger than the number of fit parameters (" << fParam.size() << ").";
cerr << endl << ">> Consider to check the manual ;-)" << endl;
return false;
}
}
// check fit range
if (!fRuns[i].IsFitRangeInBin() && !fFourierOnly) { // fit range given as times in usec (RUN block)
if ((fRuns[i].GetFitRange(0) == PMUSR_UNDEFINED) || (fRuns[i].GetFitRange(1) == PMUSR_UNDEFINED)) { // check fit range in RUN block
if (!fGlobal.IsFitRangeInBin()) { // fit range given as times in usec (GLOBAL block)
if ((fGlobal.GetFitRange(0) == PMUSR_UNDEFINED) || (fGlobal.GetFitRange(1) == PMUSR_UNDEFINED)) { // check fit range in GLOBAL block
cerr << endl << "PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << " Fit range is not defined. Necessary for single histogram fits." << endl;
return false;
}
}
}
}
// check number of T0's provided
if ((fRuns[i].GetT0BinSize() > fRuns[i].GetForwardHistoNoSize()) &&
(fGlobal.GetT0BinSize() > fRuns[i].GetForwardHistoNoSize())) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> Found " << fRuns[i].GetT0BinSize() << " T0 entries. Expecting only " << fRuns[i].GetForwardHistoNoSize() << ". Needs to be fixed." << endl;
cerr << endl << ">> In GLOBAL block: " << fGlobal.GetT0BinSize() << " T0 entries. Expecting only " << fRuns[i].GetForwardHistoNoSize() << ". Needs to be fixed." << endl;
return false;
}
// check that RRF frequency is given
if (fGlobal.GetRRFUnitTag() == RRF_UNIT_UNDEF) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** no RRF frequency found in the GLOBAL block." << endl;
return false;
}
// check that RRF packing is given
if (fGlobal.GetRRFPacking() == -1) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** no RRF packing found in the GLOBAL block." << endl;
return false;
}
break;
case PRUN_ASYMMETRY:
// check alpha
if ((fRuns[i].GetAlphaParamNo() == -1) && !fFourierOnly) {
@@ -5245,6 +5435,62 @@ Bool_t PMsrHandler::CheckRunBlockIntegrity()
fRuns[i].SetPacking(1);
}
break;
case PRUN_ASYMMETRY_RRF:
// check alpha
if ((fRuns[i].GetAlphaParamNo() == -1) && !fFourierOnly) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> alpha parameter number missing which is needed for an asymmetry RRF fit.";
cerr << endl << ">> Consider to check the manual ;-)" << endl;
return false;
}
// check that there is a forward parameter number
if (fRuns[i].GetForwardHistoNo() == -1) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> forward histogram number not defined. Necessary for asymmetry RRF fits." << endl;
return false;
}
// check that there is a backward parameter number
if (fRuns[i].GetBackwardHistoNo() == -1) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> backward histogram number not defined. Necessary for asymmetry RRF fits." << 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)) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> Fit range is not defined, also NOT present in the GLOBAL block. Necessary for asymmetry RRF fits." << endl;
return false;
}
}
}
// check number of T0's provided
if ((fRuns[i].GetT0BinSize() > 2*fRuns[i].GetForwardHistoNoSize()) &&
(fGlobal.GetT0BinSize() > 2*fRuns[i].GetForwardHistoNoSize())) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> Found " << fRuns[i].GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetForwardHistoNoSize() << " in forward. Needs to be fixed." << endl;
cerr << endl << ">> In GLOBAL block: " << fGlobal.GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetForwardHistoNoSize() << ". Needs to be fixed." << endl;
return false;
}
if ((fRuns[i].GetT0BinSize() > 2*fRuns[i].GetBackwardHistoNoSize()) &&
(fGlobal.GetT0BinSize() > 2*fRuns[i].GetBackwardHistoNoSize())) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** in RUN block number " << i+1;
cerr << endl << ">> Found " << fRuns[i].GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetBackwardHistoNoSize() << " in backward. Needs to be fixed." << endl;
cerr << endl << ">> In GLOBAL block: " << fGlobal.GetT0BinSize() << " T0 entries. Expecting only " << 2*fRuns[i].GetBackwardHistoNoSize() << ". Needs to be fixed." << endl;
return false;
}
// check that RRF frequency is given
if (fGlobal.GetRRFUnitTag() == RRF_UNIT_UNDEF) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** no RRF frequency found in the GLOBAL block." << endl;
return false;
}
// check that RRF packing is given
if (fGlobal.GetRRFPacking() == -1) {
cerr << endl << ">> PMsrHandler::CheckRunBlockIntegrity(): **ERROR** no RRF packing found in the GLOBAL block." << endl;
return false;
}
break;
case PRUN_MU_MINUS:
// needs eventually to be implemented
break;
@@ -5631,6 +5877,103 @@ void PMsrHandler::CheckMaxLikelihood()
}
}
//--------------------------------------------------------------------------
// CheckRRFSettings (public)
//--------------------------------------------------------------------------
/**
* <p>Make sure that if RRF settings are found in the GLOBAL section, the fit types
* in the RUN blocks correspond.
*/
Bool_t PMsrHandler::CheckRRFSettings()
{
Bool_t result = true;
Int_t fittype = fGlobal.GetFitType();
// first set of tests: if RRF parameters are set, check if RRF fit is chosen.
if (fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data()) != RRF_FREQ_UNDEF) {
if (fittype != -1) { // check if GLOBAL fittype is set
if ((fittype != MSR_FITTYPE_SINGLE_HISTO_RRF) &&
(fittype != MSR_FITTYPE_ASYM_RRF)) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** found GLOBAL fittype " << fittype << " and";
cerr << endl << ">> RRF settings in the GLOBAL section. This is NOT compatible. Fix it first.";
result = false;
}
} else { // GLOBAL fittype is NOT set
for (UInt_t i=0; i<fRuns.size(); i++) {
fittype = fRuns[i].GetFitType();
if ((fittype != MSR_FITTYPE_SINGLE_HISTO_RRF) &&
(fittype != MSR_FITTYPE_ASYM_RRF)) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** found RUN with fittype " << fittype << " and";
cerr << endl << ">> RRF settings in the GLOBAL section. This is NOT compatible. Fix it first.";
result = false;
break;
}
}
}
} else {
if (fGlobal.GetRRFPacking() != -1) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **WARNING** found in the GLOBAL section rrf_packing, without";
cerr << endl << ">> rrf_freq. Doesn't make any sense. Will drop rrf_packing";
cerr << endl << endl;
fGlobal.SetRRFPacking(-1);
}
if (fGlobal.GetRRFPhase() != 0.0) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **WARNING** found in the GLOBAL section rrf_phase, without";
cerr << endl << ">> rrf_freq. Doesn't make any sense. Will drop rrf_phase";
cerr << endl << endl;
fGlobal.SetRRFPhase(0.0);
}
}
// if not a RRF fit, done at this point
if ((fittype != MSR_FITTYPE_SINGLE_HISTO_RRF) &&
(fittype != MSR_FITTYPE_ASYM_RRF)) {
return true;
}
// second set of tests: if RRF fit is chosen, do I find the necessary RRF parameters?
fittype = fGlobal.GetFitType();
if ((fittype == MSR_FITTYPE_SINGLE_HISTO_RRF) ||
(fittype == MSR_FITTYPE_ASYM_RRF)) { // make sure RRF freq and RRF packing are set
if (fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data()) == RRF_FREQ_UNDEF) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** RRF fit chosen, but";
cerr << endl << ">> no RRF frequency found in the GLOBAL section! Fix it.";
return false;
}
if (fGlobal.GetRRFPacking() == -1) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** RRF fit chosen, but";
cerr << endl << ">> no RRF packing found in the GLOBAL section! Fix it.";
return false;
}
} else { // check single runs for RRF
UInt_t rrfFitCounter = 0;
for (UInt_t i=0; i<fRuns.size(); i++) {
fittype = fRuns[i].GetFitType();
if ((fittype == MSR_FITTYPE_SINGLE_HISTO_RRF) ||
(fittype == MSR_FITTYPE_ASYM_RRF)) { // make sure RRF freq and RRF packing are set
rrfFitCounter++;
}
}
if (rrfFitCounter != fRuns.size()) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** #Runs (" << fRuns.size() << ") != # RRF fits found (" << rrfFitCounter << ")";
cerr << endl << ">> This is currently not supported.";
return false;
}
if (fGlobal.GetRRFFreq(fGlobal.GetRRFUnit().Data()) == RRF_FREQ_UNDEF) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** RRF fit chosen, but";
cerr << endl << ">> no RRF frequency found in the GLOBAL section! Fix it.";
return false;
}
if (fGlobal.GetRRFPacking() == -1) {
cerr << endl << ">> PMsrHandler::CheckRRFSettings: **ERROR** RRF fit chosen, but";
cerr << endl << ">> no RRF packing found in the GLOBAL section! Fix it.";
return false;
}
}
return result;
}
//--------------------------------------------------------------------------
// GetGroupingString (public)
//--------------------------------------------------------------------------

137
src/classes/PMusr.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -34,6 +34,8 @@ using namespace std;
#include <boost/algorithm/string.hpp>
using namespace boost;
#include "TMath.h"
#include "PMusr.h"
//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -708,6 +710,10 @@ void PRawRunData::SetTempError(const UInt_t idx, const Double_t errTemp)
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
@@ -720,6 +726,135 @@ PMsrGlobalBlock::PMsrGlobalBlock()
fPacking = -1; // undefined packing/rebinning
}
//--------------------------------------------------------------------------
// GetRRFFreq (public)
//--------------------------------------------------------------------------
/**
* <p> get RRF frequency value in specific units. If units is unknown, RRF_UNDEF_FREQ will be returned.
*
* \param unit unit string in which the units shall be given
*/
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)
//--------------------------------------------------------------------------
/**
* <p> set RRF frequency value in specific units. If units is unknown, 0.0 will be set.
*
* \param RRF frequency value
* \param unit unit string in which the units shall be given
*/
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 {
cerr << endl << ">> PMsrGlobalBlock::SetRRFFreq: **ERROR** found undefined RRF unit '" << unit << "'!";
cerr << endl << ">> Will set RRF frequency to 0.0." << endl;
fRRFFreq = 0.0;
fRRFUnitTag = RRF_UNIT_UNDEF;
}
fRRFFreq = freq;
fRRFUnitTag = unitTag;
}
//--------------------------------------------------------------------------
// GetRRFUnit (public)
//--------------------------------------------------------------------------
/**
* <p> returns RRF frequency unit.
*/
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)
//--------------------------------------------------------------------------
/**
* <p> set RRF packing.
*
* \param RRF packing
*/
void PMsrGlobalBlock::SetRRFPacking(Int_t pack)
{
if (pack <= 0) {
cerr << endl << "PMsrGlobalBlock::SetRRFPacking: **WARNING** found RRF packing <= 0. Likely doesn't make any sense." << endl;
fRRFPacking = -1; // set to undefined
}
fRRFPacking = pack;
}
//--------------------------------------------------------------------------
// GetDataRange (public)
//--------------------------------------------------------------------------

2
src/classes/PMusrT0.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

4
src/classes/PPrepFourier.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2015 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -181,7 +181,7 @@ void PPrepFourier::DoBkgCorrection()
if ((fBkgRange[0] != -1) && (fBkgRange[1] != -1)) { // background range is given
// make sure that the bkg range is ok
for (UInt_t i=0; i<fRawData.size(); i++) {
if ((fBkgRange[0] >= fRawData[i].rawData.size()) || (fBkgRange[1] >= fRawData[i].rawData.size())) {
if ((fBkgRange[0] >= (Int_t)fRawData[i].rawData.size()) || (fBkgRange[1] >= (Int_t)fRawData[i].rawData.size())) {
cerr << endl << "PPrepFourier::DoBkgCorrection() **ERROR** bkg-range out of data-range!";
return;
}

8
src/classes/PRunAsymmetry.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -87,7 +87,7 @@ PRunAsymmetry::PRunAsymmetry(PMsrHandler *msrInfo, PRunDataHandler *rawData, UIn
fPacking = fMsrInfo->GetMsrGlobal()->GetPacking();
}
if (fPacking == -1) { // this should NOT happen, somethin is severely wrong
cerr << endl << ">> PRunAsymmetry::PRunAsymmetry: **SEVERE ERROR**: Couldn't find any packing information!";
cerr << endl << ">> PRunAsymmetry::PRunAsymmetry(): **SEVERE ERROR**: Couldn't find any packing information!";
cerr << endl << ">> This is very bad :-(, will quit ...";
cerr << endl;
fValid = false;
@@ -342,7 +342,7 @@ void PRunAsymmetry::SetFitRangeBin(const TString fitRange)
Int_t pos = 2*(fRunNo+1)-1;
if (pos + 1 >= tok->GetEntries()) {
cerr << endl << ">> PRunSingleHisto::SetFitRangeBin(): **ERROR** invalid FIT_RANGE command found: '" << fitRange << "'";
cerr << endl << ">> PRunAsymmetry::SetFitRangeBin(): **ERROR** invalid FIT_RANGE command found: '" << fitRange << "'";
cerr << endl << ">> will ignore it. Sorry ..." << endl;
} else {
// handle fgb+n0 entry
@@ -370,7 +370,7 @@ void PRunAsymmetry::SetFitRangeBin(const TString fitRange)
fFitEndTime = (fGoodBins[1] - offset - fT0s[0]) * fTimeResolution;
}
} else { // error
cerr << endl << ">> PRunSingleHisto::SetFitRangeBin(): **ERROR** invalid FIT_RANGE command found: '" << fitRange << "'";
cerr << endl << ">> PRunAsymmetry::SetFitRangeBin(): **ERROR** invalid FIT_RANGE command found: '" << fitRange << "'";
cerr << endl << ">> will ignore it. Sorry ..." << endl;
}

1511
src/classes/PRunAsymmetryRRF.cpp Normal file
View File

File diff suppressed because it is too large Load Diff

2
src/classes/PRunBase.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PRunDataHandler.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

256
src/classes/PRunListCollection.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *
@@ -58,12 +58,24 @@ PRunListCollection::~PRunListCollection()
}
fRunSingleHistoList.clear();
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++) {
fRunSingleHistoRRFList[i]->CleanUp();
fRunSingleHistoRRFList[i]->~PRunSingleHistoRRF();
}
fRunSingleHistoRRFList.clear();
for (UInt_t i=0; i<fRunAsymmetryList.size(); i++) {
fRunAsymmetryList[i]->CleanUp();
fRunAsymmetryList[i]->~PRunAsymmetry();
}
fRunAsymmetryList.clear();
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++) {
fRunAsymmetryRRFList[i]->CleanUp();
fRunAsymmetryRRFList[i]->~PRunAsymmetryRRF();
}
fRunAsymmetryRRFList.clear();
for (UInt_t i=0; i<fRunMuMinusList.size(); i++) {
fRunMuMinusList[i]->CleanUp();
fRunMuMinusList[i]->~PRunMuMinus();
@@ -106,11 +118,21 @@ Bool_t PRunListCollection::Add(Int_t runNo, EPMusrHandleTag tag)
if (!fRunSingleHistoList[fRunSingleHistoList.size()-1]->IsValid())
success = false;
break;
case PRUN_SINGLE_HISTO_RRF:
fRunSingleHistoRRFList.push_back(new PRunSingleHistoRRF(fMsrInfo, fData, runNo, tag));
if (!fRunSingleHistoRRFList[fRunSingleHistoRRFList.size()-1]->IsValid())
success = false;
break;
case PRUN_ASYMMETRY:
fRunAsymmetryList.push_back(new PRunAsymmetry(fMsrInfo, fData, runNo, tag));
if (!fRunAsymmetryList[fRunAsymmetryList.size()-1]->IsValid())
success = false;
break;
case PRUN_ASYMMETRY_RRF:
fRunAsymmetryRRFList.push_back(new PRunAsymmetryRRF(fMsrInfo, fData, runNo, tag));
if (!fRunAsymmetryRRFList[fRunAsymmetryRRFList.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())
@@ -147,8 +169,12 @@ void PRunListCollection::SetFitRange(const TString fitRange)
{
for (UInt_t i=0; i<fRunSingleHistoList.size(); i++)
fRunSingleHistoList[i]->SetFitRangeBin(fitRange);
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++)
fRunSingleHistoRRFList[i]->SetFitRangeBin(fitRange);
for (UInt_t i=0; i<fRunAsymmetryList.size(); i++)
fRunAsymmetryList[i]->SetFitRangeBin(fitRange);
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
fRunAsymmetryRRFList[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++)
@@ -169,8 +195,12 @@ void PRunListCollection::SetFitRange(const PDoublePairVector fitRange)
{
for (UInt_t i=0; i<fRunSingleHistoList.size(); i++)
fRunSingleHistoList[i]->SetFitRange(fitRange);
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++)
fRunSingleHistoRRFList[i]->SetFitRange(fitRange);
for (UInt_t i=0; i<fRunAsymmetryList.size(); i++)
fRunAsymmetryList[i]->SetFitRange(fitRange);
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
fRunAsymmetryRRFList[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++)
@@ -198,6 +228,27 @@ Double_t PRunListCollection::GetSingleHistoChisq(const std::vector<Double_t>& pa
return chisq;
}
//--------------------------------------------------------------------------
// GetSingleHistoRRFChisq (public)
//--------------------------------------------------------------------------
/**
* <p>Calculates chi-square of <em>all</em> single histogram RRF runs of a msr-file.
*
* <b>return:</b>
* - chi-square of all single histogram RRF runs of the msr-file
*
* \param par fit parameter vector
*/
Double_t PRunListCollection::GetSingleHistoRRFChisq(const std::vector<Double_t>& par) const
{
Double_t chisq = 0.0;
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++)
chisq += fRunSingleHistoRRFList[i]->CalcChiSquare(par);
return chisq;
}
//--------------------------------------------------------------------------
// GetAsymmetryChisq (public)
//--------------------------------------------------------------------------
@@ -219,6 +270,27 @@ Double_t PRunListCollection::GetAsymmetryChisq(const std::vector<Double_t>& par)
return chisq;
}
//--------------------------------------------------------------------------
// GetAsymmetryRRFChisq (public)
//--------------------------------------------------------------------------
/**
* <p>Calculates chi-square of <em>all</em> asymmetry RRF runs of a msr-file.
*
* <b>return:</b>
* - chi-square of all asymmetry RRF runs of the msr-file
*
* \param par fit parameter vector
*/
Double_t PRunListCollection::GetAsymmetryRRFChisq(const std::vector<Double_t>& par) const
{
Double_t chisq = 0.0;
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
chisq += fRunAsymmetryRRFList[i]->CalcChiSquare(par);
return chisq;
}
//--------------------------------------------------------------------------
// GetMuMinusChisq (public)
//--------------------------------------------------------------------------
@@ -299,9 +371,15 @@ Double_t PRunListCollection::GetSingleHistoChisqExpected(const std::vector<Doubl
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_MU_MINUS:
expectedChisq = fRunMuMinusList[subIdx]->CalcChiSquareExpected(par);
break;
@@ -336,16 +414,17 @@ Double_t PRunListCollection::GetSingleRunChisq(const std::vector<Double_t>& par,
return chisq;
}
Int_t subIdx = 0;
Int_t type = fMsrInfo->GetMsrRunList()->at(idx).GetFitType();
if (type == -1) { // i.e. not forun in the RUN block, try the GLOBAL block
if (type == -1) { // i.e. not found in the RUN block, try the GLOBAL block
type = fMsrInfo->GetMsrGlobal()->GetFitType();
}
// count how many entries of this fit-type are present up to idx
UInt_t subIdx = 0;
for (UInt_t i=0; i<idx; i++) {
if (fMsrInfo->GetMsrRunList()->at(i).GetFitType() == type)
subIdx++;
subIdx = idx;
} else { // found in the RUN block
// count how many entries of this fit-type are present up to idx
for (UInt_t i=0; i<idx; i++) {
if (fMsrInfo->GetMsrRunList()->at(i).GetFitType() == type)
subIdx++;
}
}
// return the chisq of the single run
@@ -353,9 +432,15 @@ Double_t PRunListCollection::GetSingleRunChisq(const std::vector<Double_t>& par,
case PRUN_SINGLE_HISTO:
chisq = fRunSingleHistoList[subIdx]->CalcChiSquare(par);
break;
case PRUN_SINGLE_HISTO_RRF:
chisq = fRunSingleHistoRRFList[subIdx]->CalcChiSquare(par);
break;
case PRUN_ASYMMETRY:
chisq = fRunAsymmetryList[subIdx]->CalcChiSquare(par);
break;
case PRUN_ASYMMETRY_RRF:
chisq = fRunAsymmetryRRFList[subIdx]->CalcChiSquare(par);
break;
case PRUN_MU_MINUS:
chisq = fRunMuMinusList[subIdx]->CalcChiSquare(par);
break;
@@ -376,7 +461,7 @@ Double_t PRunListCollection::GetSingleRunChisq(const std::vector<Double_t>& par,
* <p>Calculates log max-likelihood of <em>all</em> single histogram runs of a msr-file.
*
* <b>return:</b>
* - chi-square of all single histogram runs of the msr-file
* - log max-likelihood of all single histogram runs of the msr-file
*
* \param par fit parameter vector
*/
@@ -390,6 +475,27 @@ Double_t PRunListCollection::GetSingleHistoMaximumLikelihood(const std::vector<D
return mlh;
}
//--------------------------------------------------------------------------
// GetSingleHistoRRFMaximumLikelihood (public)
//--------------------------------------------------------------------------
/**
* <p>Calculates log max-likelihood of <em>all</em> single histogram RRF runs of a msr-file.
*
* <b>return:</b>
* - log max-likelihood of all single histogram runs of the msr-file
*
* \param par fit parameter vector
*/
Double_t PRunListCollection::GetSingleHistoRRFMaximumLikelihood(const std::vector<Double_t>& par) const
{
Double_t mlh = 0.0;
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++)
mlh += fRunSingleHistoRRFList[i]->CalcMaxLikelihood(par);
return mlh;
}
//--------------------------------------------------------------------------
// GetAsymmetryMaximumLikelihood (public)
//--------------------------------------------------------------------------
@@ -412,6 +518,28 @@ Double_t PRunListCollection::GetAsymmetryMaximumLikelihood(const std::vector<Dou
return mlh;
}
//--------------------------------------------------------------------------
// GetAsymmetryRRFMaximumLikelihood (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 RRF runs of the msr-file
*
* \param par fit parameter vector
*/
Double_t PRunListCollection::GetAsymmetryRRFMaximumLikelihood(const std::vector<Double_t>& par) const
{
Double_t mlh = 0.0;
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
mlh += fRunAsymmetryRRFList[i]->CalcChiSquare(par);
return mlh;
}
//--------------------------------------------------------------------------
// GetMuMinusMaximumLikelihood (public)
//--------------------------------------------------------------------------
@@ -419,7 +547,7 @@ Double_t PRunListCollection::GetAsymmetryMaximumLikelihood(const std::vector<Dou
* <p>Calculates log max-likelihood of <em>all</em> mu minus runs of a msr-file.
*
* <b>return:</b>
* - chi-square of all mu minus runs of the msr-file
* - log max-likelihood of all mu minus runs of the msr-file
*
* \param par fit parameter vector
*/
@@ -493,9 +621,15 @@ UInt_t PRunListCollection::GetNoOfBinsFitted(const UInt_t idx) const
case PRUN_SINGLE_HISTO:
result = fRunSingleHistoList[subIdx]->GetNoOfFitBins();
break;
case PRUN_SINGLE_HISTO_RRF:
result = fRunSingleHistoRRFList[subIdx]->GetNoOfFitBins();
break;
case PRUN_ASYMMETRY:
result = fRunAsymmetryList[subIdx]->GetNoOfFitBins();
break;
case PRUN_ASYMMETRY_RRF:
result = fRunAsymmetryRRFList[subIdx]->GetNoOfFitBins();
break;
case PRUN_MU_MINUS:
result = fRunMuMinusList[subIdx]->GetNoOfFitBins();
break;
@@ -526,9 +660,15 @@ UInt_t PRunListCollection::GetTotalNoOfBinsFitted() const
for (UInt_t i=0; i<fRunSingleHistoList.size(); i++)
counts += fRunSingleHistoList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++)
counts += fRunSingleHistoRRFList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunAsymmetryList.size(); i++)
counts += fRunAsymmetryList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++)
counts += fRunAsymmetryRRFList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
counts += fRunMuMinusList[i]->GetNoOfFitBins();
@@ -558,7 +698,7 @@ PRunData* PRunListCollection::GetSingleHisto(UInt_t index, EDataSwitch tag)
switch (tag) {
case kIndex:
if ((index < 0) || (index >= fRunSingleHistoList.size())) {
cerr << endl << "PRunListCollection::GetSingleHisto: **ERROR** index = " << index << " out of bounds";
cerr << endl << ">> PRunListCollection::GetSingleHisto(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
@@ -581,6 +721,49 @@ PRunData* PRunListCollection::GetSingleHisto(UInt_t index, EDataSwitch tag)
return data;
}
//--------------------------------------------------------------------------
// GetSingleHistoRRF (public)
//--------------------------------------------------------------------------
/**
* <p>Get a processed single histogram RRF 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::GetSingleHistoRRF(UInt_t index, EDataSwitch tag)
{
PRunData *data = 0;
switch (tag) {
case kIndex:
if ((index < 0) || (index >= fRunSingleHistoRRFList.size())) {
cerr << endl << ">> PRunListCollection::GetSingleHistoRRF(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
fRunSingleHistoRRFList[index]->CalcTheory();
data = fRunSingleHistoRRFList[index]->GetData();
break;
case kRunNo:
for (UInt_t i=0; i<fRunSingleHistoRRFList.size(); i++) {
if (fRunSingleHistoRRFList[i]->GetRunNo() == index) {
data = fRunSingleHistoRRFList[i]->GetData();
break;
}
}
break;
default: // error
break;
}
return data;
}
//--------------------------------------------------------------------------
// GetAsymmetry (public)
//--------------------------------------------------------------------------
@@ -601,7 +784,7 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
switch (tag) {
case kIndex: // called from musrfit when dumping the data
if ((index < 0) || (index > fRunAsymmetryList.size())) {
cerr << endl << "PRunListCollection::GetAsymmetry: **ERROR** index = " << index << " out of bounds";
cerr << endl << ">> PRunListCollection::GetAsymmetry(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
@@ -624,6 +807,49 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
return data;
}
//--------------------------------------------------------------------------
// GetAsymmetryRRF (public)
//--------------------------------------------------------------------------
/**
* <p>Get a processed asymmetry RRF 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::GetAsymmetryRRF(UInt_t index, EDataSwitch tag)
{
PRunData *data = 0;
switch (tag) {
case kIndex: // called from musrfit when dumping the data
if ((index < 0) || (index > fRunAsymmetryRRFList.size())) {
cerr << endl << ">> PRunListCollection::GetAsymmetryRRF(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
fRunAsymmetryRRFList[index]->CalcTheory();
data = fRunAsymmetryRRFList[index]->GetData();
break;
case kRunNo: // called from PMusrCanvas
for (UInt_t i=0; i<fRunAsymmetryRRFList.size(); i++) {
if (fRunAsymmetryRRFList[i]->GetRunNo() == index) {
data = fRunAsymmetryRRFList[i]->GetData();
break;
}
}
break;
default: // error
break;
}
return data;
}
//--------------------------------------------------------------------------
// GetMuMinus (public)
//--------------------------------------------------------------------------
@@ -644,7 +870,7 @@ PRunData* PRunListCollection::GetMuMinus(UInt_t index, EDataSwitch tag)
switch (tag) {
case kIndex:
if ((index < 0) || (index > fRunMuMinusList.size())) {
cerr << endl << "PRunListCollection::GetMuMinus: **ERROR** index = " << index << " out of bounds";
cerr << endl << ">> PRunListCollection::GetMuMinus(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
@@ -686,7 +912,7 @@ PRunData* PRunListCollection::GetNonMusr(UInt_t index, EDataSwitch tag)
switch (tag) {
case kIndex:
if ((index < 0) || (index > fRunNonMusrList.size())) {
cerr << endl << "PRunListCollection::GetNonMusr: **ERROR** index = " << index << " out of bounds";
cerr << endl << ">> PRunListCollection::GetNonMusr(): **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}

2
src/classes/PRunMuMinus.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PRunNonMusr.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PRunSingleHisto.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

1205
src/classes/PRunSingleHistoRRF.cpp Normal file
View File

File diff suppressed because it is too large Load Diff

2
src/classes/PStartupHandler.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PTheory.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PUserFcn.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *

2
src/classes/PUserFcnBase.cpp
View File

@@ -8,7 +8,7 @@
***************************************************************************/
/***************************************************************************
* Copyright (C) 2007-2014 by Andreas Suter *
* Copyright (C) 2007-2016 by Andreas Suter *
* andreas.suter@psi.ch *
* *
* This program is free software; you can redistribute it and/or modify *