some tweaking in order to supress ghost lines.

This commit is contained in:
suter_a 2016-01-07 16:10:47 +01:00
parent 097a6db6b0
commit 52e2c20527
2 changed files with 18 additions and 6 deletions

View File

@ -630,12 +630,24 @@ Bool_t PRunSingleHistoRRF::PrepareFitData(PRawRunData* runData, const UInt_t his
Double_t n0 = EstimateN0(errN0); Double_t n0 = EstimateN0(errN0);
// 4a) A(t) = exp(+t/tau) [N(t)-Nbkg] / N0 - 1.0 // 4a) A(t) = exp(+t/tau) [N(t)-Nbkg] / N0 - 1.0
for (Int_t i=fGoodBins[0]; i<fGoodBins[1]; i++) { for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
fForward[i] = fForward[i] / n0 - 1.0; fForward[i] = fForward[i] / n0 - 1.0;
} }
// make a DC correction of A(t) --> this is introduced to reduce ghost lines.
// This is a dirty trick and should be put on a better quantitative base!
Double_t dc_offset=0.0;
for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
dc_offset += fForward[i];
}
dc_offset /= (Double_t)(fGoodBins[1]-fGoodBins[0]+1);
cout << "info> dc_offset = " << dc_offset << endl;
for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
fForward[i] -= dc_offset;
}
// 4b) error estimate of A(t): errA(t) = exp(+t/tau)/N0 sqrt( N(t) + ([N(t)-N_bkg]/N0)^2 errN0^2 ) // 4b) error estimate of A(t): errA(t) = exp(+t/tau)/N0 sqrt( N(t) + ([N(t)-N_bkg]/N0)^2 errN0^2 )
for (Int_t i=fGoodBins[0]; i<fGoodBins[1]; i++) { for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
time_tau = (startTime + fTimeResolution * (i - fGoodBins[0])) / PMUON_LIFETIME; time_tau = (startTime + fTimeResolution * (i - fGoodBins[0])) / PMUON_LIFETIME;
exp_t_tau = exp(time_tau); exp_t_tau = exp(time_tau);
fAerr.push_back(exp_t_tau/n0*sqrt(rawNt[i]+pow(((rawNt[i]-fBackground)/n0)*errN0,2.0))); fAerr.push_back(exp_t_tau/n0*sqrt(rawNt[i]+pow(((rawNt[i]-fBackground)/n0)*errN0,2.0)));
@ -646,14 +658,14 @@ Bool_t PRunSingleHistoRRF::PrepareFitData(PRawRunData* runData, const UInt_t his
Double_t wRRF = globalBlock->GetRRFFreq("Mc"); Double_t wRRF = globalBlock->GetRRFFreq("Mc");
Double_t phaseRRF = globalBlock->GetRRFPhase()*TMath::TwoPi()/180.0; Double_t phaseRRF = globalBlock->GetRRFPhase()*TMath::TwoPi()/180.0;
Double_t time = 0.0; Double_t time = 0.0;
for (Int_t i=fGoodBins[0]; i<fGoodBins[1]; i++) { for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
time = startTime + fTimeResolution * ((Double_t)i - (Double_t)fGoodBins[0]); time = startTime + fTimeResolution * ((Double_t)i - (Double_t)fGoodBins[0]);
fForward[i] *= 2.0*cos(wRRF * time + phaseRRF); fForward[i] *= 2.0*cos(wRRF * time + phaseRRF);
} }
// 6) RRF packing // 6) RRF packing
Double_t dval=0.0; Double_t dval=0.0;
for (Int_t i=fGoodBins[0]; i<fGoodBins[1]; i++) { for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
if (fRRFPacking == 1) { if (fRRFPacking == 1) {
fData.AppendValue(fForward[i]); fData.AppendValue(fForward[i]);
} else { // RRF packing > 1 } else { // RRF packing > 1
@ -671,7 +683,7 @@ Bool_t PRunSingleHistoRRF::PrepareFitData(PRawRunData* runData, const UInt_t his
// the error estimate of the unpacked RRF asymmetry is: errA_RRF(t) \simeq exp(t/tau)/N0 sqrt( [N(t) + ((N(t)-N_bkg)/N0)^2 errN0^2] ) // the error estimate of the unpacked RRF asymmetry is: errA_RRF(t) \simeq exp(t/tau)/N0 sqrt( [N(t) + ((N(t)-N_bkg)/N0)^2 errN0^2] )
dval = 0.0; dval = 0.0;
// the packed RRF asymmetry error // the packed RRF asymmetry error
for (Int_t i=fGoodBins[0]; i<fGoodBins[1]; i++) { for (Int_t i=fGoodBins[0]; i<=fGoodBins[1]; i++) {
if (((i-fGoodBins[0]) % fRRFPacking == 0) && (i != fGoodBins[0])) { // fill data if (((i-fGoodBins[0]) % fRRFPacking == 0) && (i != fGoodBins[0])) { // fill data
fData.AppendErrorValue(sqrt(2.0*dval)/fRRFPacking); // the factor 2.0 is needed since the high frequency part is suppressed. fData.AppendErrorValue(sqrt(2.0*dval)/fRRFPacking); // the factor 2.0 is needed since the high frequency part is suppressed.
dval = 0.0; dval = 0.0;

View File

@ -58,7 +58,7 @@ class PRunSingleHistoRRF : public PRunBase
virtual Bool_t PrepareViewData(PRawRunData* runData, const UInt_t histoNo); virtual Bool_t PrepareViewData(PRawRunData* runData, const UInt_t histoNo);
private: private:
static const Double_t N0EstimateEndTime = 2.0; ///< end time in (us) over which N0 is estimated. Should evetually go into the musrfit_startup.xml static const Double_t N0EstimateEndTime = 1.0; ///< end time in (us) over which N0 is estimated. Should eventually be estimated automatically ...
UInt_t fNoOfFitBins; ///< number of bins to be fitted UInt_t fNoOfFitBins; ///< number of bins to be fitted
Double_t fBackground; ///< needed if background range is given (units: 1/bin) Double_t fBackground; ///< needed if background range is given (units: 1/bin)