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

replace RRF Asymmetry Fit dummy, by Mu Minus Single Histo Fit dummy

Browse Source
This commit is contained in:
nemu
2010-01-12 16:00:12 +00:00
parent fe853c1025
commit 41c5d45d8b
14 changed files with 101 additions and 111 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
src/classes/Makefile.am
View File

@ -17,7 +17,7 @@ h_sources = \
../include/PRunDataHandler.h \
../include/PRunListCollection.h \
../include/PRunNonMusr.h \
../include/PRunRRF.h \
../include/PRunMuMinus.h \
../include/PRunSingleHisto.h \
../include/PStartupHandler.h \
../include/PTheory.h \
@ -51,7 +51,7 @@ cpp_sources = \
PRunDataHandler.cpp \
PRunListCollection.cpp \
PRunNonMusr.cpp \
PRunRRF.cpp \
PRunMuMinus.cpp \
PRunSingleHisto.cpp \
PStartupHandler.cpp \
PTheory.cpp \

7
src/classes/PFitterFcn.cpp
View File

@ -54,12 +54,11 @@ PFitterFcn::PFitterFcn(PRunListCollection *runList, Bool_t useChi2)
fRunListCollection = runList;
// check if max likelihood is used together with asymmetry/RRF/nonMusr data.
// check if max likelihood is used together with asymmetry/nonMusr data.
// if yes place a warning since this option is not implemented and a fall back
// to chi2 will be used.
if (!fUseChi2) {
if ((fRunListCollection->GetNoOfAsymmetry() > 0) ||
(fRunListCollection->GetNoOfRRF() > 0) ||
(fRunListCollection->GetNoOfNonMusr() > 0)) {
cerr << endl << "**WARNING**: Maximum Log Likelihood Fit is only implemented for Single Histogram Fit";
cerr << endl << " Will fall back to Chi Square Fit.";
@ -91,12 +90,12 @@ Double_t PFitterFcn::operator()(const std::vector<Double_t>& par) const
if (fUseChi2) { // chi square
value += fRunListCollection->GetSingleHistoChisq(par);
value += fRunListCollection->GetAsymmetryChisq(par);
value += fRunListCollection->GetRRFChisq(par);
value += fRunListCollection->GetMuMinusChisq(par);
value += fRunListCollection->GetNonMusrChisq(par);
} else { // max likelihood
value += fRunListCollection->GetSingleHistoMaximumLikelihood(par);
value += fRunListCollection->GetAsymmetryMaximumLikelihood(par);
value += fRunListCollection->GetRRFMaximumLikelihood(par);
value += fRunListCollection->GetMuMinusMaximumLikelihood(par);
value += fRunListCollection->GetNonMusrMaximumLikelihood(par);
}

29
src/classes/PMsrHandler.cpp
View File

@ -590,8 +590,8 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
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 << " (RRF asymmetry fit)" << endl ;
case MSR_FITTYPE_MU_MINUS:
fout << left << "fittype" << MSR_FITTYPE_MU_MINUS << " (mu minus fit)" << endl ;
break;
case MSR_FITTYPE_NON_MUSR:
fout << left << "fittype" << MSR_FITTYPE_NON_MUSR << " (non muSR fit)" << endl ;
@ -837,8 +837,8 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
case MSR_PLOT_ASYM:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (asymmetry plot)" << endl;
break;
case MSR_PLOT_ASYM_RRF:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (rotating reference frame plot)" << endl;
case MSR_PLOT_MU_MINUS:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (mu minus plot)" << endl;
break;
case MSR_PLOT_NON_MUSR:
fout << "PLOT " << fPlots[plotNo].fPlotType << " (non muSR plot)" << endl;
@ -850,12 +850,8 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages)
fout << "runs ";
fout.precision(0);
for (UInt_t j=0; j<fPlots[plotNo].fRuns.size(); j++) {
if (fPlots[plotNo].fPlotType != MSR_PLOT_ASYM_RRF) { // all but MSR_PLOT_ASYM_RRF
fout.width(4);
fout << fPlots[plotNo].fRuns[j].Re();
} else { // MSR_PLOT_ASYM_RRF
fout << fPlots[plotNo].fRuns[j].Re() << "," << fPlots[plotNo].fRuns[j].Im() << " ";
}
fout.width(4);
fout << fPlots[plotNo].fRuns[j];
}
fout << endl;
} else if (sstr.BeginsWith("range")) {
@ -1603,7 +1599,7 @@ Bool_t PMsrHandler::HandleRunEntry(PMsrLines &lines)
Int_t fittype = str.Atoi();
if ((fittype == MSR_FITTYPE_SINGLE_HISTO) ||
(fittype == MSR_FITTYPE_ASYM) ||
(fittype == MSR_FITTYPE_ASYM_RRF) ||
(fittype == MSR_FITTYPE_MU_MINUS) ||
(fittype == MSR_FITTYPE_NON_MUSR)) {
param.SetFitType(fittype);
} else {
@ -2434,7 +2430,7 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
tokens = 0;
}
break;
case MSR_PLOT_ASYM_RRF: // like: runs 1,1 1,2
case MSR_PLOT_MU_MINUS: // like: runs 1,1 1,2
tokens = iter1->fLine.Tokenize(" \t");
if (!tokens) {
cerr << endl << ">> PMsrHandler::HandlePlotEntry: **SEVERE ERROR** Couldn't tokenize PLOT in line " << iter1->fLineNo;
@ -2798,14 +2794,9 @@ Bool_t PMsrHandler::HandlePlotEntry(PMsrLines &lines)
cerr << endl;
cerr << endl << "where <plot_type> is: 0=single histo asym,";
cerr << endl << " 2=forward-backward asym,";
cerr << endl << " 4=RRF asym (not implemented yet),";
cerr << endl << " 4=mu minus singhle histo (not implemented yet),";
cerr << endl << " 8=non muSR.";
cerr << endl << "<run_list> is the list of runs";
cerr << endl << " for <plot_type> 0,2,8 it is a list of run numbers, e.g. runs 1 3";
cerr << endl << " for <plot_type> 4 it is a list of 'complex' numbers, where";
cerr << endl << " the real part is the run number, and the";
cerr << endl << " imaginary one is 1=real part or 2=imag part, e.g.";
cerr << endl << " runs 1,1 1,2";
cerr << endl << "<run_list> is the list of runs, e.g. runs 1 3";
cerr << endl << "range is optional";
cerr << endl << "sub_ranges (if present) will plot the N given runs each on its own sub-range";
cerr << endl << "logx, logy (if present) will present the x-, y-axis in log-scale";

28
src/classes/PMusrCanvas.cpp
View File

@ -435,14 +435,14 @@ void PMusrCanvas::UpdateDataTheoryPad()
fPlotType = plotInfo.fPlotType;
for (UInt_t i=0; i<plotInfo.fRuns.size(); i++) {
// first check that plot number is smaller than the maximal number of runs
if ((Int_t)plotInfo.fRuns[i].Re() > (Int_t)runs.size()) {
if ((Int_t)plotInfo.fRuns[i] > (Int_t)runs.size()) {
fValid = false;
cerr << endl << "PMusrCanvas::UpdateDataTheoryPad: **ERROR** run plot number " << (Int_t)plotInfo.fRuns[i].Re() << " is larger than the number of runs " << runs.size();
cerr << endl << "PMusrCanvas::UpdateDataTheoryPad: **ERROR** run plot number " << (Int_t)plotInfo.fRuns[i] << " is larger than the number of runs " << runs.size();
cerr << endl;
return;
}
// check that the plottype and the fittype do correspond
runNo = (UInt_t)plotInfo.fRuns[i].Re()-1;
runNo = (UInt_t)plotInfo.fRuns[i]-1;
//cout << endl << ">> runNo = " << runNo;
//cout << endl;
if (fPlotType != runs[runNo].GetFitType()) {
@ -457,7 +457,7 @@ void PMusrCanvas::UpdateDataTheoryPad()
for (UInt_t i=0; i<plotInfo.fRuns.size(); i++) {
// get run data and create a histogram
data = 0;
runNo = (UInt_t)plotInfo.fRuns[i].Re()-1;
runNo = (UInt_t)plotInfo.fRuns[i]-1;
// get data depending on the fittype
switch (runs[runNo].GetFitType()) {
case MSR_FITTYPE_SINGLE_HISTO:
@ -484,12 +484,12 @@ void PMusrCanvas::UpdateDataTheoryPad()
// handle data
HandleDataSet(i, runNo, data);
break;
case MSR_FITTYPE_ASYM_RRF:
data = fRunList->GetRRF(runNo, PRunListCollection::kRunNo);
case MSR_FITTYPE_MU_MINUS:
data = fRunList->GetMuMinus(runNo, PRunListCollection::kRunNo);
if (!data) { // something wrong
fValid = false;
// error message
cerr << endl << "PMusrCanvas::UpdateDataTheoryPad: **ERROR** couldn't obtain run no " << runNo << " for a RRF plot";
cerr << endl << "PMusrCanvas::UpdateDataTheoryPad: **ERROR** couldn't obtain run no " << runNo << " for a mu minus single histogram plot";
cerr << endl;
return;
}
@ -583,7 +583,7 @@ void PMusrCanvas::UpdateInfoPad()
PMsrRunList runs = *fMsrHandler->GetMsrRunList();
for (UInt_t i=0; i<fData.size(); i++) {
// run label = run_name/histo/T=0K/B=0G/E=0keV/...
runNo = (UInt_t)plotInfo.fRuns[i].Re()-1;
runNo = (UInt_t)plotInfo.fRuns[i]-1;
if (runs[runNo].GetRunNameSize() > 1)
tstr = "++" + *runs[runNo].GetRunName() + TString(","); // run_name
else
@ -2569,9 +2569,11 @@ void PMusrCanvas::PlotData()
}
break;
case MSR_PLOT_ASYM:
case MSR_PLOT_ASYM_RRF:
yAxisTitle = "asymmetry";
break;
case MSR_PLOT_MU_MINUS:
yAxisTitle = "??";
break;
default:
yAxisTitle = "??";
break;
@ -2601,7 +2603,7 @@ void PMusrCanvas::PlotData()
PMsrRunList runs = *fMsrHandler->GetMsrRunList();
PMsrPlotStructure plotInfo = fMsrHandler->GetMsrPlotList()->at(fPlotNumber);
UInt_t runNo = (UInt_t)plotInfo.fRuns[0].Re()-1;
UInt_t runNo = (UInt_t)plotInfo.fRuns[0]-1;
TString xAxisTitle = fRunList->GetXAxisTitle(*runs[runNo].GetRunName(), runNo);
TString yAxisTitle = fRunList->GetYAxisTitle(*runs[runNo].GetRunName(), runNo);
@ -2660,7 +2662,7 @@ void PMusrCanvas::PlotData()
assert(fMultiGraphLegend != 0);
PStringVector legendLabel;
for (UInt_t i=0; i<plotInfo.fRuns.size(); i++) {
runNo = (UInt_t)plotInfo.fRuns[i].Re()-1;
runNo = (UInt_t)plotInfo.fRuns[i]-1;
xAxisTitle = fRunList->GetXAxisTitle(*runs[runNo].GetRunName(), runNo);
yAxisTitle = fRunList->GetYAxisTitle(*runs[runNo].GetRunName(), runNo);
legendLabel.push_back(xAxisTitle + " vs. " + yAxisTitle);
@ -2724,7 +2726,7 @@ void PMusrCanvas::PlotDifference()
PMsrRunList runs = *fMsrHandler->GetMsrRunList();
PMsrPlotStructure plotInfo = fMsrHandler->GetMsrPlotList()->at(fPlotNumber);
UInt_t runNo = (UInt_t)plotInfo.fRuns[0].Re()-1;
UInt_t runNo = (UInt_t)plotInfo.fRuns[0]-1;
TString xAxisTitle = fRunList->GetXAxisTitle(*runs[runNo].GetRunName(), runNo);
// if fMultiGraphDiff is not present create it and add the diff data
@ -3491,7 +3493,7 @@ void PMusrCanvas::SaveDataAscii()
switch (fPlotType) {
case MSR_PLOT_SINGLE_HISTO:
case MSR_PLOT_ASYM:
case MSR_PLOT_ASYM_RRF:
case MSR_PLOT_MU_MINUS:
if (fDifferenceView) { // difference view plot
switch (fCurrentPlotView) {
case PV_DATA:

4
src/classes/PMusrT0.cpp
View File

@ -256,8 +256,8 @@ void PMusrT0::InitDataAndBkg()
fAddRunOffset = 2;
} else if (fitType == MSR_FITTYPE_ASYM) {
fAddRunOffset = 4;
} else if (fitType == MSR_FITTYPE_ASYM_RRF) {
fAddRunOffset = 8;
} else if (fitType == MSR_FITTYPE_MU_MINUS) {
fAddRunOffset = 2;
}
// feed data range histo

53
src/classes/PRunListCollection.cpp
View File

@ -71,12 +71,12 @@ PRunListCollection::~PRunListCollection()
}
fRunAsymmetryList.clear();
//cout << endl << ">> fRunRRFList.size() = " << fRunRRFList.size();
for (UInt_t i=0; i<fRunRRFList.size(); i++) {
fRunRRFList[i]->CleanUp();
fRunRRFList[i]->~PRunRRF();
//cout << endl << ">> fRunMuMinusList.size() = " << fRunMuMinusList.size();
for (UInt_t i=0; i<fRunMuMinusList.size(); i++) {
fRunMuMinusList[i]->CleanUp();
fRunMuMinusList[i]->~PRunMuMinus();
}
fRunRRFList.clear();
fRunMuMinusList.clear();
//cout << endl << ">> fRunNonMusrList.size() = " << fRunNonMusrList.size();
for (UInt_t i=0; i<fRunNonMusrList.size(); i++) {
@ -118,9 +118,9 @@ Bool_t PRunListCollection::Add(Int_t runNo, EPMusrHandleTag tag)
if (!fRunAsymmetryList[fRunAsymmetryList.size()-1]->IsValid())
success = false;
break;
case PRUN_RRF:
fRunRRFList.push_back(new PRunRRF(fMsrInfo, fData, runNo, tag));
if (!fRunRRFList[fRunRRFList.size()-1]->IsValid())
case PRUN_MU_MINUS:
fRunMuMinusList.push_back(new PRunMuMinus(fMsrInfo, fData, runNo, tag));
if (!fRunMuMinusList[fRunMuMinusList.size()-1]->IsValid())
success = false;
break;
case PRUN_NON_MUSR:
@ -169,17 +169,17 @@ Double_t PRunListCollection::GetAsymmetryChisq(const std::vector<Double_t>& par)
}
//--------------------------------------------------------------------------
// GetRRFChisq
// GetMuMinusChisq
//--------------------------------------------------------------------------
/**
* <p>
*/
Double_t PRunListCollection::GetRRFChisq(const std::vector<Double_t>& par) const
Double_t PRunListCollection::GetMuMinusChisq(const std::vector<Double_t>& par) const
{
Double_t chisq = 0.0;
for (UInt_t i=0; i<fRunRRFList.size(); i++)
chisq += fRunRRFList[i]->CalcChiSquare(par);
for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
chisq += fRunMuMinusList[i]->CalcChiSquare(par);
return chisq;
}
@ -234,18 +234,17 @@ Double_t PRunListCollection::GetAsymmetryMaximumLikelihood(const std::vector<Dou
}
//--------------------------------------------------------------------------
// GetRRFMaximumLikelihood
// GetMuMinusMaximumLikelihood
//--------------------------------------------------------------------------
/**
* <p> Since it is not clear yet how to handle RRF fits with max likelihood
* the chi square will be used!
* <p>
*/
Double_t PRunListCollection::GetRRFMaximumLikelihood(const std::vector<Double_t>& par) const
Double_t PRunListCollection::GetMuMinusMaximumLikelihood(const std::vector<Double_t>& par) const
{
Double_t mlh = 0.0;
for (UInt_t i=0; i<fRunRRFList.size(); i++)
mlh += fRunRRFList[i]->CalcChiSquare(par);
for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
mlh += fRunMuMinusList[i]->CalcMaxLikelihood(par);
return mlh;
}
@ -283,8 +282,8 @@ UInt_t PRunListCollection::GetTotalNoOfBinsFitted() const
for (UInt_t i=0; i<fRunAsymmetryList.size(); i++)
counts += fRunAsymmetryList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunRRFList.size(); i++)
counts += fRunRRFList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunMuMinusList.size(); i++)
counts += fRunMuMinusList[i]->GetNoOfFitBins();
for (UInt_t i=0; i<fRunNonMusrList.size(); i++)
counts += fRunNonMusrList[i]->GetNoOfFitBins();
@ -372,7 +371,7 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
}
//--------------------------------------------------------------------------
// GetRRF
// GetMuMinus
//--------------------------------------------------------------------------
/**
* <p>
@ -380,19 +379,25 @@ PRunData* PRunListCollection::GetAsymmetry(UInt_t index, EDataSwitch tag)
* \param index
* \param tag kIndex -> data at index, kRunNo -> data of given run no
*/
PRunData* PRunListCollection::GetRRF(UInt_t index, EDataSwitch tag)
PRunData* PRunListCollection::GetMuMinus(UInt_t index, EDataSwitch tag)
{
PRunData *data = 0;
switch (tag) {
case kIndex:
if ((index < 0) || (index > fRunRRFList.size())) {
cerr << endl << "PRunListCollection::GetRRF: **ERROR** index = " << index << " out of bounds";
if ((index < 0) || (index > fRunMuMinusList.size())) {
cerr << endl << "PRunListCollection::GetMuMinus: **ERROR** index = " << index << " out of bounds";
cerr << endl;
return 0;
}
break;
case kRunNo:
for (UInt_t i=0; i<fRunMuMinusList.size(); i++) {
if (fRunMuMinusList[i]->GetRunNo() == index) {
data = fRunMuMinusList[i]->GetData();
break;
}
}
break;
default: // error
break;

20
src/classes/PRunRRF.cpp → src/classes/PRunMuMinus.cpp
View File

@ -1,6 +1,6 @@
/***************************************************************************
PRunRRF.cpp
PRunMuMinus.cpp
Author: Andreas Suter
e-mail: andreas.suter@psi.ch
@ -31,7 +31,7 @@
#include <iostream>
#include "PRunRRF.h"
#include "PRunMuMinus.h"
//--------------------------------------------------------------------------
// Constructor
@ -40,7 +40,7 @@
* <p>
*
*/
PRunRRF::PRunRRF() : PRunBase()
PRunMuMinus::PRunMuMinus() : PRunBase()
{
fFitStartTime = 0.0;
fFitStopTime = 0.0;
@ -58,7 +58,7 @@ PRunRRF::PRunRRF() : PRunBase()
* \param msrInfo pointer to the msr info structure
* \param runNo number of the run of the msr-file
*/
PRunRRF::PRunRRF(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag) : PRunBase(msrInfo, rawData, runNo, tag)
PRunMuMinus::PRunMuMinus(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, EPMusrHandleTag tag) : PRunBase(msrInfo, rawData, runNo, tag)
{
Bool_t success;
@ -75,7 +75,7 @@ PRunRRF::PRunRRF(PMsrHandler *msrInfo, PRunDataHandler *rawData, UInt_t runNo, E
* <p>
*
*/
PRunRRF::~PRunRRF()
PRunMuMinus::~PRunMuMinus()
{
}
@ -87,7 +87,7 @@ PRunRRF::~PRunRRF()
*
* \param par parameter vector iterated by minuit
*/
Double_t PRunRRF::CalcChiSquare(const std::vector<Double_t>& par)
Double_t PRunMuMinus::CalcChiSquare(const std::vector<Double_t>& par)
{
Double_t chisq = 0.0;
Double_t diff = 0.0;
@ -103,7 +103,7 @@ Double_t PRunRRF::CalcChiSquare(const std::vector<Double_t>& par)
*
* \param par parameter vector iterated by minuit
*/
Double_t PRunRRF::CalcMaxLikelihood(const std::vector<Double_t>& par)
Double_t PRunMuMinus::CalcMaxLikelihood(const std::vector<Double_t>& par)
{
cout << endl << "PRunSingleHisto::CalcMaxLikelihood(): not implemented yet ..." << endl;
@ -117,7 +117,7 @@ Double_t PRunRRF::CalcMaxLikelihood(const std::vector<Double_t>& par)
* <p>
*
*/
void PRunRRF::CalcTheory()
void PRunMuMinus::CalcTheory()
{
}
@ -128,11 +128,11 @@ void PRunRRF::CalcTheory()
* <p>
*
*/
Bool_t PRunRRF::PrepareData()
Bool_t PRunMuMinus::PrepareData()
{
Bool_t success = true;
cout << endl << "in PRunRRF::PrepareData(): will feed fData";
cout << endl << "in PRunMuMinus::PrepareData(): will feed fData";
return success;
}

4
src/classes/PRunNonMusr.cpp
View File

@ -313,9 +313,9 @@ Bool_t PRunNonMusr::PrepareViewData()
for (UInt_t j=0; j<plotBlock.fRuns.size(); j++) {
// cout << endl << ">> j=" << j;
// cout << endl << ">> fRunNo=" << fRunNo;
// cout << endl << ">> plotBlock.fRuns[j].Re()=" << plotBlock.fRuns[j].Re();
// cout << endl << ">> plotBlock.fRuns[j]=" << plotBlock.fRuns[j];
// cout << endl;
if (fRunNo == plotBlock.fRuns[j].Re()-1) { // run found
if (fRunNo == plotBlock.fRuns[j]-1) { // run found
if (first) {
first = false;
xMin = plotBlock.fTmin[0];