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

some more work towards musrt0. Furthermore improved data structure

Browse Source
This commit is contained in:
nemu 2009-03-31 07:15:43 +00:00
parent f632aa5ba9
commit 41264e9eee
6 changed files with 199 additions and 71 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

113
src/classes/PMsrHandler.cpp
View File

@ -80,6 +80,17 @@ PMsrHandler::~PMsrHandler()
fParam.clear();
fTheory.clear();
fFunctions.clear();
for (unsigned int i=0; i<fRuns.size(); i++) {
fRuns[i].fRunName.clear();
fRuns[i].fBeamline.clear();
fRuns[i].fInstitute.clear();
fRuns[i].fFileFormat.clear();
fRuns[i].fMap.clear();
fRuns[i].fBkgFix.clear();
fRuns[i].fBkgRange.clear();
fRuns[i].fDataRange.clear();
fRuns[i].fT0.clear();
}
fRuns.clear();
fCommands.clear();
fPlots.clear();
@ -586,38 +597,31 @@ int PMsrHandler::WriteMsrLogFile()
CheckAndWriteComment(f, ++lineNo);
}
// backgr.fix
if (fRuns[i].fBkgFixPresent[0]) {
if (fRuns[i].fBkgFix.size() > 0) {
f.width(15);
f << endl << left << "backgr.fix";
f.precision(prec);
f.width(12);
f << left << fRuns[i].fBkgFix[0];
if (fRuns[i].fBkgFixPresent[1]) {
for (unsigned int j=0; j<fRuns[i].fBkgFix.size(); j++) {
f.precision(prec);
f.width(12);
f << left << fRuns[i].fBkgFix[1];
f << left << fRuns[i].fBkgFix[j];
}
CheckAndWriteComment(f, ++lineNo);
}
// background
if (fRuns[i].fBkgRange[0] != -1) {
if (fRuns[i].fBkgRange.size() > 0) {
f.width(16);
f << endl << left << "background";
for (unsigned int j=0; j<4; j++) {
if (fRuns[i].fBkgRange[j] == -1)
break;
for (unsigned int j=0; j<fRuns[i].fBkgRange.size(); j++) {
f.width(8);
f << left << fRuns[i].fBkgRange[j];
}
CheckAndWriteComment(f, ++lineNo);
}
// data
if (fRuns[i].fDataRange[0] != -1) {
if (fRuns[i].fDataRange.size() > 0) {
f.width(16);
f << endl << left << "data";
for (unsigned int j=0; j<4; j++) {
if (fRuns[i].fDataRange[j] == -1)
break;
for (unsigned int j=0; j<fRuns[i].fDataRange.size(); j++) {
f.width(8);
f << left << fRuns[i].fDataRange[j];
}
@ -943,6 +947,60 @@ bool PMsrHandler::SetMsrParamPosError(unsigned int i, double value)
return true;
}
//--------------------------------------------------------------------------
// SetMsrDataRangeEntry (public)
//--------------------------------------------------------------------------
/**
* <p>
*
* \param runNo
* \param idx
* \param bin
*/
void PMsrHandler::SetMsrDataRangeEntry(unsigned int runNo, unsigned int idx, int bin)
{
if ((runNo < 0) || (runNo > fRuns.size())) { // error
cout << endl << "**ERROR** in PMsrHandler::SetMsrDataRangeEntry: runNo = " << runNo << ", is out of valid range 0.." << fRuns.size();
cout << endl;
return;
}
if ((idx < 0) || (idx > fRuns[runNo].fDataRange.size())) { // error
cout << endl << "**ERROR** in PMsrHandler::SetMsrDataRangeEntry: idx = " << idx << ", is out of valid range 0.." << fRuns[runNo].fDataRange.size();
cout << endl;
return;
}
fRuns[runNo].fDataRange[idx] = bin;
}
//--------------------------------------------------------------------------
// SetMsrBkgRangeEntry (public)
//--------------------------------------------------------------------------
/**
* <p>
*
* \param runNo
* \param idx
* \param bin
*/
void PMsrHandler::SetMsrBkgRangeEntry(unsigned int runNo, unsigned int idx, int bin)
{
if ((runNo < 0) || (runNo > fRuns.size())) { // error
cout << endl << "**ERROR** in PMsrHandler::SetMsrBkgRangeEntry: runNo = " << runNo << ", is out of valid range 0.." << fRuns.size();
cout << endl;
return;
}
if ((idx < 0) || (idx > fRuns[runNo].fBkgRange.size())) { // error
cout << endl << "**ERROR** in PMsrHandler::SetMsrBkgRangeEntry: idx = " << idx << ", is out of valid range 0.." << fRuns[runNo].fBkgRange.size();
cout << endl;
return;
}
fRuns[runNo].fBkgRange[idx] = bin;
}
//--------------------------------------------------------------------------
// ParameterInUse (public)
//--------------------------------------------------------------------------
@ -1464,14 +1522,14 @@ bool PMsrHandler::HandleRunEntry(PMsrLines &lines)
// backgr.fix ----------------------------------------------
if (iter->fLine.BeginsWith("backgr.fix", TString::kIgnoreCase)) {
if ((tokens->GetEntries() != 2) && (tokens->GetEntries() != 3)) {
if (tokens->GetEntries() < 2) {
error = true;
} else {
for (int i=1; i<tokens->GetEntries(); i++) {
ostr = dynamic_cast<TObjString*>(tokens->At(i));
str = ostr->GetString();
if (str.IsFloat())
param.fBkgFix[i-1] = str.Atof();
param.fBkgFix.push_back(str.Atof());
else
error = true;
}
@ -1480,14 +1538,14 @@ bool PMsrHandler::HandleRunEntry(PMsrLines &lines)
// background ---------------------------------------------
if (iter->fLine.BeginsWith("background", TString::kIgnoreCase)) {
if ((tokens->GetEntries() != 3) && (tokens->GetEntries() != 5)) {
if ((tokens->GetEntries() < 3) || (tokens->GetEntries() % 2 != 1)) { // odd number (>=3) of entries needed
error = true;
} else {
for (int i=1; i<tokens->GetEntries(); i++) {
ostr = dynamic_cast<TObjString*>(tokens->At(i));
str = ostr->GetString();
if (str.IsDigit())
param.fBkgRange[i-1] = str.Atoi();
param.fBkgRange.push_back(str.Atoi());
else
error = true;
}
@ -1496,14 +1554,14 @@ bool PMsrHandler::HandleRunEntry(PMsrLines &lines)
// data --------------------------------------------------
if (iter->fLine.BeginsWith("data", TString::kIgnoreCase)) {
if ((tokens->GetEntries() != 3) && (tokens->GetEntries() != 5)) {
if ((tokens->GetEntries() < 3) || (tokens->GetEntries() % 2 != 1)) { // odd number (>=3) of entries needed
error = true;
} else {
for (int i=1; i<tokens->GetEntries(); i++) {
ostr = dynamic_cast<TObjString*>(tokens->At(i));
str = ostr->GetString();
if (str.IsDigit())
param.fDataRange[i-1] = str.Atoi();
param.fDataRange.push_back(str.Atoi());
else
error = true;
}
@ -1690,9 +1748,9 @@ bool PMsrHandler::HandleRunEntry(PMsrLines &lines)
bool found;
if (fRuns[i].fBkgFitParamNo >= 0) { // check if backgr.fit is given
found = true;
} else if (fRuns[i].fBkgFixPresent[0]) { // check if backgr.fix is given
} else if (fRuns[i].fBkgFix.size() > 0) { // check if backgr.fix is given
found = true;
} else if (fRuns[i].fBkgRange[0] >= 0) { // check if background window is given
} else if (fRuns[i].fBkgRange.size() > 0) { // check if background window is given
found = true;
} else {
found = false;
@ -1736,14 +1794,9 @@ void PMsrHandler::InitRunParameterStructure(PMsrRunStructure &param)
param.fMap.clear(); // empty list
param.fForwardHistoNo = -1;
param.fBackwardHistoNo = -1;
for (int i=0; i<2; i++) {
param.fBkgFixPresent[i] = false;
param.fBkgFix[i] = 0.0;
}
for (int i=0; i<4; i++)
param.fBkgRange[i] = -1;
for (int i=0; i<4; i++)
param.fDataRange[i] = -1;
param.fBkgFix.clear();
param.fBkgRange.clear();
param.fDataRange.clear();
param.fT0.clear();
for (int i=0; i<4; i++)
param.fFitRange[i] = -1;

76
src/classes/PMusrT0.cpp
View File

@ -420,6 +420,24 @@ void PMusrT0::SetDataFirstChannel()
fDataRange[0] = fHisto->GetXaxis()->FindFixBin(x);
// set the data first bin in msr-Handler
unsigned int idx = 0;
switch(fDetectorTag) {
case DETECTOR_TAG_FORWARD:
idx = fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_BACKWARD:
idx = 2 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_LEFT:
idx = 4 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_RIGHT:
idx = 6 + fAddRunNo * fAddRunOffset;
break;
default:
break;
}
fMsrHandler->SetMsrDataRangeEntry(fRunNo, idx, fDataRange[0]);
// shift line to the proper position
fFirstDataLine->SetX1(x);
@ -461,7 +479,25 @@ void PMusrT0::SetDataLastChannel()
// get binx to set the data last channel corresponding to fPx
fDataRange[1] = fHisto->GetXaxis()->FindFixBin(x);
// set the data last bin in msr-Handler
// set the data first bin in msr-Handler
unsigned int idx = 0;
switch(fDetectorTag) {
case DETECTOR_TAG_FORWARD:
idx = 1 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_BACKWARD:
idx = 3 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_LEFT:
idx = 5 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_RIGHT:
idx = 7 + fAddRunNo * fAddRunOffset;
break;
default:
break;
}
fMsrHandler->SetMsrDataRangeEntry(fRunNo, idx, fDataRange[1]);
// shift line to the proper position
fLastDataLine->SetX1(x);
@ -504,6 +540,24 @@ void PMusrT0::SetBkgFirstChannel()
fBkgRange[0] = fHisto->GetXaxis()->FindFixBin(x);
// set the background first bin in msr-Handler
unsigned int idx = 0;
switch(fDetectorTag) {
case DETECTOR_TAG_FORWARD:
idx = fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_BACKWARD:
idx = 2 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_LEFT:
idx = 4 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_RIGHT:
idx = 6 + fAddRunNo * fAddRunOffset;
break;
default:
break;
}
fMsrHandler->SetMsrBkgRangeEntry(fRunNo, idx, fBkgRange[0]);
// shift line to the proper position
fFirstBkgLine->SetX1(x);
@ -545,7 +599,25 @@ void PMusrT0::SetBkgLastChannel()
// get binx to set the background last channel corresponding to fPx
fBkgRange[1] = fHisto->GetXaxis()->FindFixBin(x);
// set the background last bin in msr-Handler
// set the background first bin in msr-Handler
unsigned int idx = 0;
switch(fDetectorTag) {
case DETECTOR_TAG_FORWARD:
idx = 1 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_BACKWARD:
idx = 3 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_LEFT:
idx = 5 + fAddRunNo * fAddRunOffset;
break;
case DETECTOR_TAG_RIGHT:
idx = 7 + fAddRunNo * fAddRunOffset;
break;
default:
break;
}
fMsrHandler->SetMsrBkgRangeEntry(fRunNo, idx, fBkgRange[1]);
// shift line to the proper position
fLastBkgLine->SetX1(x);

5
src/classes/PRunAsymmetry.cpp
View File

@ -436,8 +436,9 @@ bool PRunAsymmetry::PrepareData()
}
// subtract background from histogramms ------------------------------------------
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given
if (fRunInfo->fBkgRange[0] != 0) {
cout << endl << ">> fBkgFix.size() = " << fRunInfo->fBkgFix.size() << ", fBkgRange.size() = " << fRunInfo->fBkgRange.size() << endl;
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given
if (fRunInfo->fBkgRange.size() != 0) {
if (!SubtractEstimatedBkg())
return false;
} else { // no background given to do the job

14
src/classes/PRunSingleHisto.cpp
View File

@ -119,7 +119,7 @@ double PRunSingleHisto::CalcChiSquare(const std::vector<double>& par)
// get background
double bkg;
if (fRunInfo->fBkgFitParamNo == -1) { // bkg not fitted
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given (background interval)
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given (background interval)
bkg = fBackground;
} else { // fixed bkg given
bkg = fRunInfo->fBkgFix[0];
@ -193,7 +193,7 @@ double PRunSingleHisto::CalcMaxLikelihood(const std::vector<double>& par)
// get background
double bkg;
if (fRunInfo->fBkgFitParamNo == -1) { // bkg not fitted
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given (background interval)
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given (background interval)
bkg = fBackground;
} else { // fixed bkg given
bkg = fRunInfo->fBkgFix[0];
@ -267,7 +267,7 @@ void PRunSingleHisto::CalcTheory()
// get background
double bkg;
if (fRunInfo->fBkgFitParamNo == -1) { // bkg not fitted
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given (background interval)
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given (background interval)
bkg = fBackground;
} else { // fixed bkg given
bkg = fRunInfo->fBkgFix[0];
@ -605,8 +605,8 @@ bool PRunSingleHisto::PrepareFitData(PRawRunData* runData, const unsigned int hi
// check how the background shall be handled
if (fRunInfo->fBkgFitParamNo == -1) { // bkg shall **NOT** be fitted
// subtract background from histogramms ------------------------------------------
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given
if (fRunInfo->fBkgRange[0] != 0) {
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given
if (fRunInfo->fBkgRange.size() != 0) {
if (!EstimateBkg(histoNo))
return false;
} else { // no background given to do the job
@ -771,7 +771,7 @@ cout << endl << ">> data start time = " << fData.fDataTimeStart;
// get background
double bkg;
if (fRunInfo->fBkgFitParamNo == -1) { // bkg not fitted
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given (background interval)
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given (background interval)
if (!EstimateBkg(histoNo))
return false;
bkg = fBackground;
@ -885,7 +885,7 @@ bool PRunSingleHisto::PrepareViewData(PRawRunData* runData, const unsigned int h
// get background
double bkg;
if (fRunInfo->fBkgFitParamNo == -1) { // bkg not fitted
if (!fRunInfo->fBkgFixPresent[0]) { // no fixed background given (background interval)
if (fRunInfo->fBkgFix.size() == 0) { // no fixed background given (background interval)
if (!EstimateBkg(histoNo))
return false;
bkg = fBackground;

3
src/include/PMsrHandler.h
View File

@ -71,6 +71,9 @@ class PMsrHandler
virtual bool SetMsrParamPosErrorPresent(unsigned int i, bool value);
virtual bool SetMsrParamPosError(unsigned int i, double value);
virtual void SetMsrDataRangeEntry(unsigned int runNo, unsigned int idx, int bin);
virtual void SetMsrBkgRangeEntry(unsigned int runNo, unsigned int idx, int bin);
virtual void SetMsrStatisticMin(double min) { fStatistic.fMin = min; }
virtual void SetMsrStatisticNdf(unsigned int ndf) { fStatistic.fNdf = ndf; }

59
src/include/PMusr.h
View File

@ -242,36 +242,35 @@ typedef vector<PMsrParamStructure> PMsrParamList;
*
*/
typedef struct {
PStringVector fRunName; ///< name of the run file
PStringVector fBeamline; ///< e.g. mue4, mue1, pim3, emu, m15, ... (former: run type)
PStringVector fInstitute; ///< e.g. psi, ral, triumf (former: run format)
PStringVector fFileFormat; ///< e.g. root, nexus, psi-bin, mud, ascii, db
int fFitType; ///< fit type: 0=single histo fit, 2=asymmetry fit, 4=asymmetry in RRF, 8=non muSR
int fAlphaParamNo; ///< alpha parameter number (fit type 2, 4)
int fBetaParamNo; ///< beta parameter number (fit type 2, 4)
int fNormParamNo; ///< N0 parameter number (fit type 0)
int fBkgFitParamNo; ///< background fit parameter number (fit type 0)
int fPhaseParamNo; ///< ??? NEEDED ??? NEEDS TO BE CHECKED !!!
int fLifetimeParamNo; ///< muon lifetime parameter number (fit type 0)
bool fLifetimeCorrection; ///< lifetime correction flag for viewing (fit type 0)
PIntVector fMap; ///< map vector needed to switch parameters for different runs within a single theory
int fForwardHistoNo; ///< forward histogram number (fit type 0, 2, 4)
int fBackwardHistoNo; ///< backward histogram number (fit type 2, 4)
bool fBkgFixPresent[2]; ///< flag showing if a fixed background is present (fit type 0, 2, 4).
double fBkgFix[2]; ///< fixed background in (1/ns) (fit type 0, 2, 4)
int fBkgRange[4]; ///< background bin range (fit type 0, 2, 4)
int fDataRange[4]; ///< data bin range (fit type 0, 2, 4)
PIntVector fT0; ///< t0 bins (fit type 0, 2, 4). if fit type 0 -> f0, f1, f2, ...; if fit type 2, 4 -> f0, b0, f1, b1, ...
double fFitRange[2]; ///< fit range in (us)
int fPacking; ///< packing/rebinning
double fRRFFreq; ///< rotating reference frequency (fit type 4)
int fRRFPacking; ///< rotating reference packing (fit type 4)
int fAlpha2ParamNo; ///< rotating reference alpha2 (fit type 4)
int fBeta2ParamNo; ///< rotating reference beta2 (fit type 4)
int fRightHistoNo; ///< rotating reference right histogram number (fit type 4)
int fLeftHistoNo; ///< rotating reference left histogram number (fit type 4)
int fXYDataIndex[2]; ///< used to get the data indices when using db-files (fit type 8)
TString fXYDataLabel[2]; ///< used to get the indices via labels when using db-files (fit type 8)
PStringVector fRunName; ///< name of the run file
PStringVector fBeamline; ///< e.g. mue4, mue1, pim3, emu, m15, ... (former: run type)
PStringVector fInstitute; ///< e.g. psi, ral, triumf (former: run format)
PStringVector fFileFormat; ///< e.g. root, nexus, psi-bin, mud, ascii, db
int fFitType; ///< fit type: 0=single histo fit, 2=asymmetry fit, 4=asymmetry in RRF, 8=non muSR
int fAlphaParamNo; ///< alpha parameter number (fit type 2, 4)
int fBetaParamNo; ///< beta parameter number (fit type 2, 4)
int fNormParamNo; ///< N0 parameter number (fit type 0)
int fBkgFitParamNo; ///< background fit parameter number (fit type 0)
int fPhaseParamNo; ///< ??? NEEDED ??? NEEDS TO BE CHECKED !!!
int fLifetimeParamNo; ///< muon lifetime parameter number (fit type 0)
bool fLifetimeCorrection; ///< lifetime correction flag for viewing (fit type 0)
PIntVector fMap; ///< map vector needed to switch parameters for different runs within a single theory
int fForwardHistoNo; ///< forward histogram number (fit type 0, 2, 4)
int fBackwardHistoNo; ///< backward histogram number (fit type 2, 4)
PDoubleVector fBkgFix; ///< fixed background in (1/ns) (fit type 0, 2, 4)
PIntVector fBkgRange; ///< background bin range (fit type 0, 2, 4)
PIntVector fDataRange; ///< data bin range (fit type 0, 2, 4)
PIntVector fT0; ///< t0 bins (fit type 0, 2, 4). if fit type 0 -> f0, f1, f2, ...; if fit type 2, 4 -> f0, b0, f1, b1, ...
double fFitRange[2]; ///< fit range in (us)
int fPacking; ///< packing/rebinning
double fRRFFreq; ///< rotating reference frequency (fit type 4)
int fRRFPacking; ///< rotating reference packing (fit type 4)
int fAlpha2ParamNo; ///< rotating reference alpha2 (fit type 4)
int fBeta2ParamNo; ///< rotating reference beta2 (fit type 4)
int fRightHistoNo; ///< rotating reference right histogram number (fit type 4)
int fLeftHistoNo; ///< rotating reference left histogram number (fit type 4)
int fXYDataIndex[2]; ///< used to get the data indices when using db-files (fit type 8)
TString fXYDataLabel[2]; ///< used to get the indices via labels when using db-files (fit type 8)
} PMsrRunStructure;
//-------------------------------------------------------------