diff --git a/src/classes/PMsrHandler.cpp b/src/classes/PMsrHandler.cpp index 2f903fe4..a8a617cc 100644 --- a/src/classes/PMsrHandler.cpp +++ b/src/classes/PMsrHandler.cpp @@ -604,7 +604,37 @@ Int_t PMsrHandler::WriteMsrLogFile(const Bool_t messages) } break; case MSR_TAG_GLOBAL: - // STILL MISSING + // STILL QUITE A FEW ITEMS MISSING !!!!!! + sstr = str; + if (sstr.BeginsWith("fittype")) { + fout.width(16); + switch (fGlobal.GetFitType()) { + case MSR_FITTYPE_SINGLE_HISTO: + fout << left << "fittype" << MSR_FITTYPE_SINGLE_HISTO << " (single histogram fit)" << endl; + break; + case MSR_FITTYPE_ASYM: + fout << left << "fittype" << MSR_FITTYPE_ASYM << " (asymmetry fit)" << endl ; + break; + 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 ; + break; + default: + break; + } + } else if (sstr.BeginsWith("data")) { + // STILL MISSING + } else if (sstr.BeginsWith("t0")) { + // STILL MISSING + } else if (sstr.BeginsWith("fit")) { + // STILL MISSING + } else if (sstr.BeginsWith("packing")) { + // STILL MISSING + } else { + fout << str.Data() << endl; + } break; case MSR_TAG_RUN: sstr = str; @@ -2646,13 +2676,49 @@ Bool_t PMsrHandler::HandleGlobalEntry(PMsrLines &lines) if (tokens->GetEntries() < 3) { error = true; } else { // fit given in time, i.e. fit , where , are given as doubles - for (Int_t i=1; i<3; i++) { - ostr = dynamic_cast(tokens->At(i)); + if (iter->fLine.Contains("fgb", TString::kIgnoreCase)) { // fit given in bins, i.e. fit fgb+n0 lgb-n1 + // check 1st entry, i.e. fgb[+n0] + ostr = dynamic_cast(tokens->At(1)); str = ostr->GetString(); - if (str.IsFloat()) - global.SetFitRange(str.Atof(), i-1); - else - error = true; + Ssiz_t idx = str.First("+"); + TString numStr = str; + if (idx > -1) { // '+' present hence extract n0 + numStr.Remove(0,idx+1); + if (numStr.IsFloat()) { + global.SetFitRangeOffset(numStr.Atoi(), 0); + } else { + error = true; + } + } else { // n0 == 0 + global.SetFitRangeOffset(0, 0); + } + // check 2nd entry, i.e. lgb[-n1] + ostr = dynamic_cast(tokens->At(2)); + str = ostr->GetString(); + idx = str.First("-"); + numStr = str; + if (idx > -1) { // '-' present hence extract n1 + numStr.Remove(0,idx+1); + if (numStr.IsFloat()) { + global.SetFitRangeOffset(numStr.Atoi(), 1); + } else { + error = true; + } + } else { // n0 == 0 + global.SetFitRangeOffset(0, 0); + } + + if (!error) + global.SetFitRangeInBins(true); + } else { // fit given in time, i.e. fit , where , are given as doubles + for (Int_t i=1; i<3; i++) { + ostr = dynamic_cast(tokens->At(i)); + str = ostr->GetString(); + if (str.IsFloat()) + global.SetFitRange(str.Atof(), i-1); + else + error = true; + } } } } else if (iter->fLine.BeginsWith("packing", TString::kIgnoreCase)) { // packing diff --git a/src/classes/PMusr.cpp b/src/classes/PMusr.cpp index e2ca724b..0756c25f 100644 --- a/src/classes/PMusr.cpp +++ b/src/classes/PMusr.cpp @@ -709,8 +709,11 @@ PMsrGlobalBlock::PMsrGlobalBlock() for (UInt_t i=0; i<4; i++) { fDataRange[i] = -1; // undefined data bin range } + fFitRangeInBins = false; // default is that fit range is given in time NOT bins fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range + fFitRangeOffset[0] = -1; // undefined start fit range offset + fFitRangeOffset[1] = -1; // undefined end fit range offset fPacking = -1; // undefined packing/rebinning } @@ -833,6 +836,43 @@ void PMsrGlobalBlock::SetFitRange(Double_t dval, UInt_t idx) fFitRange[idx] = dval; } +//-------------------------------------------------------------------------- +// GetFitRangeOffset (public) +//-------------------------------------------------------------------------- +/** + *

get fit range offset value at position idx. idx: 0=fit range offset start, 1=fit range offset end. + * + * return: + * - fit range offset value, if idx is within proper boundaries + * - -1, otherwise + * + * \param idx index of the fit range value to be returned + */ +Int_t PMsrGlobalBlock::GetFitRangeOffset(UInt_t idx) +{ + if (idx >= 2) + return -1; + + return fFitRangeOffset[idx]; +} + +//-------------------------------------------------------------------------- +// SetFitRangeOffset (public) +//-------------------------------------------------------------------------- +/** + *

set fit range offset value at position idx. Illegale values will be ignored. + * + * \param ival value to be set + * \param idx index of the fit range value to be set + */ +void PMsrGlobalBlock::SetFitRangeOffset(Int_t ival, UInt_t idx) +{ + if (idx >= 2) + return; + + fFitRangeOffset[idx] = ival; +} + //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ // implementation PMsrRunBlock //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ @@ -863,8 +903,8 @@ PMsrRunBlock::PMsrRunBlock() fFitRangeInBins = false; // default is that fit range is given in time NOT bins fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range - fFitRangeOffset[0] = 0; // undefined start fit range offset - fFitRangeOffset[1] = 0; // undefined end fit range offset + fFitRangeOffset[0] = -1; // undefined start fit range offset + fFitRangeOffset[1] = -1; // undefined end fit range offset fPacking = -1; // undefined packing fXYDataIndex[0] = -1; // undefined x data index (NonMusr) fXYDataIndex[1] = -1; // undefined y data index (NonMusr) @@ -916,8 +956,8 @@ void PMsrRunBlock::CleanUp() fFitRangeInBins = false; // default is that fit range is given in time NOT bins fFitRange[0] = PMUSR_UNDEFINED; // undefined start fit range fFitRange[1] = PMUSR_UNDEFINED; // undefined end fit range - fFitRangeOffset[0] = 0; // undefined start fit range offset - fFitRangeOffset[1] = 0; // undefined end fit range offset + fFitRangeOffset[0] = -1; // undefined start fit range offset + fFitRangeOffset[1] = -1; // undefined end fit range offset fPacking = -1; // undefined packing fXYDataIndex[0] = -1; // undefined x data index (NonMusr) fXYDataIndex[1] = -1; // undefined y data index (NonMusr) diff --git a/src/classes/PRunBase.cpp b/src/classes/PRunBase.cpp index 265fc86f..22706345 100644 --- a/src/classes/PRunBase.cpp +++ b/src/classes/PRunBase.cpp @@ -54,8 +54,8 @@ PRunBase::PRunBase() fRawData = 0; fTimeResolution = -1.0; - fFitStartTime = 0.0; - fFitEndTime = 0.0; + fFitStartTime = PMUSR_UNDEFINED; + fFitEndTime = PMUSR_UNDEFINED; fValid = true; fHandleTag = kEmpty; diff --git a/src/classes/PRunSingleHisto.cpp b/src/classes/PRunSingleHisto.cpp index 3023b182..8f91a8dc 100644 --- a/src/classes/PRunSingleHisto.cpp +++ b/src/classes/PRunSingleHisto.cpp @@ -58,6 +58,7 @@ PRunSingleHisto::PRunSingleHisto() : PRunBase() fScaleN0AndBkg = true; fNoOfFitBins = 0; fBackground = 0; + fPacking = -1; // the 2 following variables are need in case fit range is given in bins, and since // the fit range can be changed in the command block, these variables need to be accessible @@ -81,6 +82,18 @@ PRunSingleHisto::PRunSingleHisto(PMsrHandler *msrInfo, PRunDataHandler *rawData, fScaleN0AndBkg = IsScaleN0AndBkg(); fNoOfFitBins = 0; + fPacking = fRunInfo->GetPacking(); + if (fPacking == -1) { // i.e. packing is NOT given in the RUN-block, it must be given in the GLOBAL-block + fPacking = fMsrInfo->GetMsrGlobal()->GetPacking(); + } + if (fPacking == -1) { // this should NOT happen, somethin is severely wrong + cerr << endl << ">> PRunSingleHisto::PRunSingleHisto: **SEVERE ERROR**: Couldn't find any packing information!"; + cerr << endl << ">> This is very bad :-(, will quit ..."; + cerr << endl; + fValid = false; + return; + } + // the 2 following variables are need in case fit range is given in bins, and since // the fit range can be changed in the command block, these variables need to be accessible fGoodBins[0] = -1; @@ -192,7 +205,7 @@ Double_t PRunSingleHisto::CalcChiSquare(const std::vector& par) // the correction factor is need since the data scales like pack*t_res, // whereas the error scales like sqrt(pack*t_res) if (fScaleN0AndBkg) - chisq *= fRunInfo->GetPacking() * (fTimeResolution * 1.0e3); + chisq *= fPacking * (fTimeResolution * 1.0e3); return chisq; } @@ -285,7 +298,7 @@ Double_t PRunSingleHisto::CalcChiSquareExpected(const std::vector& par // the correction factor is need since the data scales like pack*t_res, // whereas the error scales like sqrt(pack*t_res) if (fScaleN0AndBkg) - chisq *= fRunInfo->GetPacking() * (fTimeResolution * 1.0e3); + chisq *= fPacking * (fTimeResolution * 1.0e3); return chisq; } @@ -352,7 +365,7 @@ Double_t PRunSingleHisto::CalcMaxLikelihood(const std::vector& par) Double_t normalizer = 1.0; if (fScaleN0AndBkg) - normalizer = fRunInfo->GetPacking() * (fTimeResolution * 1.0e3); + normalizer = fPacking * (fTimeResolution * 1.0e3); // In order not to have an IF in the next loop, determine the start and end bins for the fit range now Int_t startTimeBin = static_cast(ceil((fFitStartTime - fData.GetDataTimeStart())/fData.GetDataTimeStep())); @@ -608,6 +621,9 @@ Bool_t PRunSingleHisto::PrepareData() { Bool_t success = true; + // keep the Global block info + PMsrGlobalBlock *globalBlock = fMsrInfo->GetMsrGlobal(); + // get the proper run PRawRunData* runData = fRawData->GetRunData(*fRunInfo->GetRunName()); if (!runData) { // couldn't get run @@ -644,13 +660,21 @@ Bool_t PRunSingleHisto::PrepareData() fT0s[i] = fRunInfo->GetT0Bin(i); } + // fill in the T0's from the GLOBAL block section (if present) + for (UInt_t i=0; iGetT0BinSize(); i++) { + if (fT0s[i] == -1) { // i.e. not given in the RUN block section + fT0s[i] = globalBlock->GetT0Bin(i); + } + } + // fill in the T0's from the data file, if not already present in the msr-file for (UInt_t i=0; iGetT0Bin(histoNo[i]) > 0.0) { fT0s[i] = runData->GetT0Bin(histoNo[i]); fRunInfo->SetT0Bin(fT0s[i], i); // keep value for the msr-file } + } } // fill in the T0's gaps, i.e. in case the T0's are NOT in the msr-file and NOT in the data file @@ -827,6 +851,15 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN Int_t end; start = fRunInfo->GetDataRange(0); end = fRunInfo->GetDataRange(1); + + // check if data range has been provided, and if not try to get it from the GLOBAL block section + if (start < 0) { + start = fMsrInfo->GetMsrGlobal()->GetDataRange(0); + } + if (end < 0) { + end = fMsrInfo->GetMsrGlobal()->GetDataRange(1); + } + // check if data range has been provided, and if not try to estimate them if (start < 0) { Int_t offset = (Int_t)(10.0e-3/fTimeResolution); @@ -864,7 +897,7 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN return false; } - // keep good bins for potential latter use + // keep good bins for potential later use fGoodBins[0] = start; fGoodBins[1] = end; @@ -877,6 +910,11 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN fRunInfo->SetFitRange(fFitEndTime, 1); } + // check if fit range is given in the run block. If not, i.e. it has to be found in the global section! + if (fRunInfo->GetFitRange(0) == PMUSR_UNDEFINED) { +// ANYTHING NEEDED AT THIS POINT??? as35 + } + // check how the background shall be handled if (fRunInfo->GetBkgFitParamNo() == -1) { // bkg shall **NOT** be fitted // subtract background from histogramms ------------------------------------------ @@ -908,13 +946,13 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN // in order that after rebinning the fit does not need to be redone (important for plots) // the value is normalize to per 1 nsec if scaling is whished if (fScaleN0AndBkg) - normalizer = fRunInfo->GetPacking() * (fTimeResolution * 1.0e3); // fTimeResolution us->ns + normalizer = fPacking * (fTimeResolution * 1.0e3); // fTimeResolution us->ns // data start at data_start-t0 // time shifted so that packing is included correctly, i.e. t0 == t0 after packing - fData.SetDataTimeStart(fTimeResolution*((Double_t)start-(Double_t)t0+(Double_t)(fRunInfo->GetPacking()-1)/2.0)); - fData.SetDataTimeStep(fTimeResolution*fRunInfo->GetPacking()); + fData.SetDataTimeStart(fTimeResolution*((Double_t)start-(Double_t)t0+(Double_t)(fPacking-1)/2.0)); + fData.SetDataTimeStep(fTimeResolution*fPacking); for (Int_t i=start; iGetPacking() == 1) { + if (fPacking == 1) { value = fForward[i]; value /= normalizer; fData.AppendValue(value); @@ -922,8 +960,8 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN fData.AppendErrorValue(1.0/normalizer); else fData.AppendErrorValue(TMath::Sqrt(value)); - } else { // packed data, i.e. fRunInfo->GetPacking() > 1 - if (((i-start) % fRunInfo->GetPacking() == 0) && (i != start)) { // fill data + } else { // packed data, i.e. fPacking > 1 + if (((i-start) % fPacking == 0) && (i != start)) { // fill data value /= normalizer; fData.AppendValue(value); if (value == 0.0) @@ -964,7 +1002,7 @@ Bool_t PRunSingleHisto::PrepareFitData(PRawRunData* runData, const UInt_t histoN Bool_t PRunSingleHisto::PrepareRawViewData(PRawRunData* runData, const UInt_t histoNo) { // check if view_packing is wished - Int_t packing = fRunInfo->GetPacking(); + Int_t packing = fPacking; if (fMsrInfo->GetMsrPlotList()->at(0).fViewPacking > 0) { packing = fMsrInfo->GetMsrPlotList()->at(0).fViewPacking; } @@ -974,7 +1012,7 @@ Bool_t PRunSingleHisto::PrepareRawViewData(PRawRunData* runData, const UInt_t hi if (fScaleN0AndBkg) { dataNorm = 1.0/ (packing * (fTimeResolution * 1.0e3)); // fTimeResolution us->ns } else if (!fScaleN0AndBkg && (fMsrInfo->GetMsrPlotList()->at(0).fViewPacking > 0)) { - theoryNorm = (Double_t)fMsrInfo->GetMsrPlotList()->at(0).fViewPacking/(Double_t)fRunInfo->GetPacking(); + theoryNorm = (Double_t)fMsrInfo->GetMsrPlotList()->at(0).fViewPacking/(Double_t)packing; } // raw data, since PMusrCanvas is doing ranging etc. @@ -1157,7 +1195,7 @@ Bool_t PRunSingleHisto::PrepareRawViewData(PRawRunData* runData, const UInt_t hi Bool_t PRunSingleHisto::PrepareViewData(PRawRunData* runData, const UInt_t histoNo) { // check if view_packing is wished. This is a global option for all PLOT blocks! - Int_t packing = fRunInfo->GetPacking(); + Int_t packing = fPacking; if (fMsrInfo->GetMsrPlotList()->at(0).fViewPacking > 0) { packing = fMsrInfo->GetMsrPlotList()->at(0).fViewPacking; } @@ -1171,7 +1209,7 @@ Bool_t PRunSingleHisto::PrepareViewData(PRawRunData* runData, const UInt_t histo if (fScaleN0AndBkg) { dataNorm = 1.0/ (packing * (fTimeResolution * 1.0e3)); // fTimeResolution us->ns } else if (!fScaleN0AndBkg && (fMsrInfo->GetMsrPlotList()->at(0).fViewPacking > 0)) { - theoryNorm = (Double_t)fMsrInfo->GetMsrPlotList()->at(0).fViewPacking/(Double_t)fRunInfo->GetPacking(); + theoryNorm = (Double_t)fMsrInfo->GetMsrPlotList()->at(0).fViewPacking/(Double_t)packing; } // transform raw histo data. This is done the following way (for details see the manual): @@ -1461,7 +1499,7 @@ void PRunSingleHisto::EstimateN0() if (fScaleN0AndBkg) { N0 /= fTimeResolution*1.0e3; } else { - N0 *= fRunInfo->GetPacking(); + N0 *= fPacking; } cout << ">> PRunSingleHisto::EstimateN0: found N0=" << param->at(paramNo-1).fValue << ", will set it to N0=" << N0 << endl; @@ -1507,10 +1545,10 @@ Bool_t PRunSingleHisto::EstimateBkg(UInt_t histoNo) // calculate proper background range if (beamPeriod != 0.0) { - Double_t timeBkg = (Double_t)(end-start)*(fTimeResolution*fRunInfo->GetPacking()); // length of the background intervall in time + Double_t timeBkg = (Double_t)(end-start)*(fTimeResolution*fPacking); // length of the background intervall in time UInt_t fullCycles = (UInt_t)(timeBkg/beamPeriod); // how many proton beam cylces can be placed within the proposed background intervall // correct the end of the background intervall such that the background is as close as possible to a multiple of the proton cylce - end = start + (UInt_t) ((fullCycles*beamPeriod)/(fTimeResolution*fRunInfo->GetPacking())); + end = start + (UInt_t) ((fullCycles*beamPeriod)/(fTimeResolution*fPacking)); cout << endl << "PRunSingleHisto::EstimatBkg(): Background " << start << ", " << end; if (end == start) end = fRunInfo->GetBkgRange(1); @@ -1545,7 +1583,7 @@ Bool_t PRunSingleHisto::EstimateBkg(UInt_t histoNo) if (fScaleN0AndBkg) fBackground = bkg / (fTimeResolution * 1e3); // keep background (per 1 nsec) for chisq, max.log.likelihood, fTimeResolution us->ns else - fBackground = bkg * fRunInfo->GetPacking(); // keep background (per bin) + fBackground = bkg * fPacking; // keep background (per bin) fRunInfo->SetBkgEstimated(fBackground, 0); diff --git a/src/include/PMusr.h b/src/include/PMusr.h index 3b5a56cb..a2b7e634 100644 --- a/src/include/PMusr.h +++ b/src/include/PMusr.h @@ -537,21 +537,27 @@ class PMsrGlobalBlock { virtual Int_t GetDataRange(UInt_t idx); virtual UInt_t GetT0BinSize() { return fT0.size(); } virtual Double_t GetT0Bin(UInt_t idx=0); + virtual Bool_t IsFitRangeInBin() { return fFitRangeInBins; } virtual Double_t GetFitRange(UInt_t idx); + virtual Int_t GetFitRangeOffset(UInt_t idx); virtual Int_t GetPacking() { return fPacking; } virtual void SetFitType(Int_t ival) { fFitType = ival; } virtual void SetDataRange(Int_t ival, Int_t idx); virtual void SetT0Bin(Double_t dval, Int_t idx=-1); + virtual void SetFitRangeInBins(Bool_t bval) { fFitRangeInBins = bval; } virtual void SetFitRange(Double_t dval, UInt_t idx); + virtual void SetFitRangeOffset(Int_t ival, UInt_t idx); virtual void SetPacking(Int_t ival) { fPacking = ival; } private: - Int_t fFitType; ///< fit type: 0=single histo fit, 2=asymmetry fit, 4=mu^- single histo fit, 8=non muSR fit - Int_t fDataRange[4]; ///< data bin range (fit type 0, 2, 4) - PDoubleVector fT0; ///< t0 bins (fit type 0, 2, 4). if fit type 0 -> f0, f1, f2, ...; if fit type - Double_t fFitRange[2]; ///< fit range in (us) - Int_t fPacking; ///< packing/rebinning + Int_t fFitType; ///< fit type: 0=single histo fit, 2=asymmetry fit, 4=mu^- single histo fit, 8=non muSR fit + Int_t fDataRange[4]; ///< data bin range (fit type 0, 2, 4) + PDoubleVector fT0; ///< t0 bins (fit type 0, 2, 4). if fit type 0 -> f0, f1, f2, ...; if fit type + Bool_t fFitRangeInBins; ///< flag telling if fit range is given in time or in bins + Double_t fFitRange[2]; ///< fit range in (us) + Int_t fFitRangeOffset[2]; ///< if fit range is given in bins it can have the form fit fgb+n0 lgb-n1. This variable holds the n0 and n1. + Int_t fPacking; ///< packing/rebinning }; //------------------------------------------------------------- diff --git a/src/include/PRunSingleHisto.h b/src/include/PRunSingleHisto.h index 2640f0be..118df87e 100644 --- a/src/include/PRunSingleHisto.h +++ b/src/include/PRunSingleHisto.h @@ -62,6 +62,7 @@ class PRunSingleHisto : public PRunBase Bool_t fScaleN0AndBkg; ///< true=scale N0 and background to 1/ns, otherwise 1/bin UInt_t fNoOfFitBins; ///< number of bins to be fitted Double_t fBackground; ///< needed if background range is given (units: 1/bin) + Int_t fPacking; ///< packing for this particular run. Either given in the RUN- or GLOBAL-block. Int_t fGoodBins[2]; ///< keep first/last good bins. 0=fgb, 1=lgb