Changeset 7720 for trunk


Ignore:
Timestamp:
05/19/06 18:29:35 (18 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/Changelog

    r7719 r7720  
    4343   * mreport/MReportFileReadCC.cc:
    4444     - always output the file format version
     45
     46   * mhflux/MAlphaFitter.[h,cc]:
     47     - implemented fitting of the off-data for background determination
     48     - the result values are not yet accessible
     49
     50   * mhflux/MHThetaSqN.[h,cc]:
     51     - some improvements regarding the treatment of the signal-regions
     52     - allow to set a different cut-level for off-cuts
    4553
    4654

  • trunk/MagicSoft/Mars/NEWS

    r7719 r7720  
    44
    55   - merpp: better handling of problems with the TH and TD part of the
    6        CC-REPORT for files older than 200507190 and 200412210
    7        respectively
     6     CC-REPORT for files older than 200507190 and 200412210
     7     respectively
     8
     9   - ganymed: implemented a new class (MHThetaSqN) which allows to use more
     10     than one off-source region in wobble-mode. To use it add the following to
     11     your ganymed_wobble.rc:
     12       + MJCut.NameHist: MHThetaSqN             (switch the new feature on)
     13       + MHThetaSqN.NumOffSourcePos: 3          (define the number of off-regions)
     14       + MHThetaSqN.DoOffCut: Yes,No            (switch on/off the off-cut)
     15       + Cut1.ThetaCut: None
     16       + MHThetaSqN.SignificanceCutLevel: 2.0   (increase off-cut by 2.0/1.7)
    817
    918

  • trunk/MagicSoft/Mars/mhflux/MAlphaFitter.cc

    r7694 r7720  
    3838//  + Double_t fSignificanceExc;  // significance of a known excess
    3939//
     40// Version 3:
     41// ----------
     42//  + TArrayD fErrors;  // errors of coefficients
     43//
    4044//
    4145//////////////////////////////////////////////////////////////////////////////
    … …  
    7579
    7680    fCoefficients.Reset();
     81    fErrors.Reset();
    7782}
    7883
    7984// --------------------------------------------------------------------------
    8085//
    81 // This is a preliminary implementation of a alpha-fit procedure for
    82 // all possible source positions. It will be moved into its own
    83 // more powerfull class soon.
    84 //
    85 // The fit function is "gaus(0)+pol2(3)" which is equivalent to:
    86 //   [0]*exp(-0.5*((x-[1])/[2])^2) + [3] + [4]*x + [5]*x^2
    87 // or
    88 //   A*exp(-0.5*((x-mu)/sigma)^2) + a + b*x + c*x^2
    89 //
    90 // Parameter [1] is fixed to 0 while the alpha peak should be
    91 // symmetric around alpha=0.
    92 //
    93 // Parameter [4] is fixed to 0 because the first derivative at
    94 // alpha=0 should be 0, too.
    95 //
    96 // In a first step the background is fitted between bgmin and bgmax,
    97 // while the parameters [0]=0 and [2]=1 are fixed.
    98 //
    99 // In a second step the signal region (alpha<sigmax) is fittet using
    100 // the whole function with parameters [1], [3], [4] and [5] fixed.
    101 //
    102 // The number of excess and background events are calculated as
    103 //   s = int(hist,    0, 1.25*sigint)
    104 //   b = int(pol2(3), 0, 1.25*sigint)
    105 //
    106 // The Significance is calculated using the Significance() member
    107 // function.
    108 //
    109 Bool_t MAlphaFitter::Fit(TH1D &h, Bool_t paint)
    110 {
    111     Clear();
     86// This function implementes the fit to the off-data as used in Fit()
     87//
     88Bool_t MAlphaFitter::FitOff(TH1D &h, Int_t paint)
     89{
    11290    if (h.GetEntries()==0)
    11391        return kFALSE;
    114 
    115     Double_t sigmax=fSigMax;
    116     Double_t bgmin =fBgMin;
    117     Double_t bgmax =fBgMax;
    118 
    119     const Double_t alpha0 = h.GetBinContent(1);
    120     const Double_t alphaw = h.GetXaxis()->GetBinWidth(1);
    121 
    122     // Check for the regios which is not filled...
    123     if (alpha0==0)
    124         return kFALSE; //*fLog << warn << "Histogram empty." << endl;
    12592
    12693    // First fit a polynom in the off region
    … …  
    138105            fFunc->SetParameter(i, 0);
    139106
     107    if (!fFitBackground)
     108        return kTRUE;
     109
    140110    // options : N  do not store the function, do not draw
    141111    //           I  use integral of function in bin rather than value at bin center
    … …  
    143113    //           Q  quiet mode
    144114    //           E  Perform better Errors estimation using Minos technique
    145     if (fFitBackground)
    146     {
    147         h.Fit(fFunc, "NQI", "", bgmin, bgmax);
    148         fChiSqBg = fFunc->GetChisquare()/fFunc->GetNDF();
    149     }
    150 
     115    h.Fit(fFunc, "NQIE", "", fBgMin, fBgMax);
     116    fChiSqBg = fFunc->GetChisquare()/fFunc->GetNDF();
     117
     118    fCoefficients.Set(fFunc->GetNpar(), fFunc->GetParameters());
     119    fErrors.Set(fFunc->GetNpar());
     120    for (int i=3; i<fFunc->GetNpar(); i++)
     121        fErrors[i] = fFunc->GetParError(i);
    151122
    152123    // ------------------------------------
    153     if (paint && fFitBackground)
    154     {
    155         fFunc->SetRange(0, 90);
    156         fFunc->SetLineColor(kRed);
    157         fFunc->SetLineWidth(2);
     124
     125    if (paint)
     126    {
     127        if (paint==2)
     128        {
     129            fFunc->SetLineColor(kBlack);
     130            fFunc->SetLineWidth(1);
     131        }
     132        else
     133        {
     134            fFunc->SetRange(0, 90);
     135            fFunc->SetLineColor(kRed);
     136            fFunc->SetLineWidth(2);
     137        }
    158138        fFunc->Paint("same");
    159139    }
    160     // ------------------------------------
    161 
    162     fFunc->ReleaseParameter(0);  // It is also released by SetParLimits later on
    163     //func.ReleaseParameter(1);  // It is also released by SetParLimits later on
    164     fFunc->ReleaseParameter(2);
    165     for (int i=3; i<fFunc->GetNpar(); i++)
    166         fFunc->FixParameter(i, fFunc->GetParameter(i));
    167 
    168     // Do not allow signals smaller than the background
    169     const Double_t s  = fSignalFunc==kGauss ? fFunc->GetParameter(3) : TMath::Exp(fFunc->GetParameter(3));
    170     const Double_t A  = alpha0-s;
    171     const Double_t dA = TMath::Abs(A);
    172     fFunc->SetParLimits(0, -dA*4, dA*4);  // SetParLimits also releases the parameter
    173     fFunc->SetParLimits(2, 0, 90);        // SetParLimits also releases the parameter
    174 
    175     // Now fit a gaus in the on region on top of the polynom
    176     fFunc->SetParameter(0, A);
    177     fFunc->SetParameter(2, sigmax*0.75);
    178 
    179     // options : N  do not store the function, do not draw
    180     //           I  use integral of function in bin rather than value at bin center
    181     //           R  use the range specified in the function range
    182     //           Q  quiet mode
    183     //           E  Perform better Errors estimation using Minos technique
    184     h.Fit(fFunc, "NQI", "", 0, sigmax);
    185 
    186     fChiSqSignal = fFunc->GetChisquare()/fFunc->GetNDF();
    187     fCoefficients.Set(fFunc->GetNpar(), fFunc->GetParameters());
    188 
    189     //const Bool_t ok = NDF>0 && chi2<2.5*NDF;
    190 
    191     // ------------------------------------
    192     if (paint)
    193     {
    194         fFunc->SetLineColor(kGreen);
    195         fFunc->SetLineWidth(2);
    196         fFunc->Paint("same");
    197     }
    198     // ------------------------------------
    199 
    200      //const Double_t s = fFunc->Integral(0, fSigInt)/alphaw;
    201      fFunc->SetParameter(0, 0);
    202      fFunc->SetParameter(2, 1);
    203      //const Double_t b = fFunc->Integral(0, fSigInt)/alphaw;
    204      //fSignificance = MMath::SignificanceLiMaSigned(s, b);
     140
     141    return kTRUE;
     142}
     143
     144// --------------------------------------------------------------------------
     145//
     146// Calculate the result of the fit and set the corresponding data members
     147//
     148void MAlphaFitter::FitResult(const TH1D &h)
     149{
     150    const Double_t alphaw = h.GetXaxis()->GetBinWidth(1);
    205151
    206152    const Int_t bin = h.GetXaxis()->FindFixBin(fSigInt*0.999);
    … …  
    218164    if (fEventsExcess<0)
    219165        fEventsExcess=0;
     166}
     167
     168// --------------------------------------------------------------------------
     169//
     170// This is a preliminary implementation of a alpha-fit procedure for
     171// all possible source positions. It will be moved into its own
     172// more powerfull class soon.
     173//
     174// The fit function is "gaus(0)+pol2(3)" which is equivalent to:
     175//   [0]*exp(-0.5*((x-[1])/[2])^2) + [3] + [4]*x + [5]*x^2
     176// or
     177//   A*exp(-0.5*((x-mu)/sigma)^2) + a + b*x + c*x^2
     178//
     179// Parameter [1] is fixed to 0 while the alpha peak should be
     180// symmetric around alpha=0.
     181//
     182// Parameter [4] is fixed to 0 because the first derivative at
     183// alpha=0 should be 0, too.
     184//
     185// In a first step the background is fitted between bgmin and bgmax,
     186// while the parameters [0]=0 and [2]=1 are fixed.
     187//
     188// In a second step the signal region (alpha<sigmax) is fittet using
     189// the whole function with parameters [1], [3], [4] and [5] fixed.
     190//
     191// The number of excess and background events are calculated as
     192//   s = int(hist,    0, 1.25*sigint)
     193//   b = int(pol2(3), 0, 1.25*sigint)
     194//
     195// The Significance is calculated using the Significance() member
     196// function.
     197//
     198Bool_t MAlphaFitter::Fit(TH1D &h, Bool_t paint)
     199{
     200    Clear();
     201
     202    // Check for the region which is not filled...
     203    // if (alpha0==0)
     204    //     return kFALSE;
     205
     206    // Perform fit to the off-data
     207    if (!FitOff(h, paint))
     208        return kFALSE;
     209
     210    fFunc->ReleaseParameter(0);  // It is also released by SetParLimits later on
     211    //func.ReleaseParameter(1);  // It is also released by SetParLimits later on
     212    fFunc->ReleaseParameter(2);
     213    for (int i=3; i<fFunc->GetNpar(); i++)
     214        fFunc->FixParameter(i, fFunc->GetParameter(i));
     215
     216
     217    // Do not allow signals smaller than the background
     218    const Double_t alpha0 = h.GetBinContent(1);
     219    const Double_t s      = fSignalFunc==kGauss ? fFunc->GetParameter(3) : TMath::Exp(fFunc->GetParameter(3));
     220    const Double_t A      = alpha0-s;
     221    //const Double_t dA     = TMath::Abs(A);
     222    //fFunc->SetParLimits(0, -dA*4, dA*4);  // SetParLimits also releases the parameter
     223    fFunc->SetParLimits(2, 0, 90);        // SetParLimits also releases the parameter
     224
     225    // Now fit a gaus in the on region on top of the polynom
     226    fFunc->SetParameter(0, A);
     227    fFunc->SetParameter(2, fSigMax*0.75);
     228
     229    // options : N  do not store the function, do not draw
     230    //           I  use integral of function in bin rather than value at bin center
     231    //           R  use the range specified in the function range
     232    //           Q  quiet mode
     233    //           E  Perform better Errors estimation using Minos technique
     234    h.Fit(fFunc, "NQI", "", 0, fSigMax);
     235
     236    fChiSqSignal = fFunc->GetChisquare()/fFunc->GetNDF();
     237    fCoefficients.Set(fFunc->GetNpar(), fFunc->GetParameters());
     238    for (int i=0; i<3; i++)
     239        fErrors[i] = fFunc->GetParError(i);
     240    //const Bool_t ok = NDF>0 && chi2<2.5*NDF;
     241
     242    // ------------------------------------
     243    if (paint)
     244    {
     245        fFunc->SetLineColor(kGreen);
     246        fFunc->SetLineWidth(2);
     247        fFunc->Paint("same");
     248    }
     249    // ------------------------------------
     250
     251     //const Double_t s = fFunc->Integral(0, fSigInt)/alphaw;
     252     fFunc->SetParameter(0, 0);
     253     fFunc->SetParameter(2, 1);
     254     //const Double_t b = fFunc->Integral(0, fSigInt)/alphaw;
     255     //fSignificance = MMath::SignificanceLiMaSigned(s, b);
     256
     257     // Calculate the fit result and set the corresponding data members
     258     FitResult(h);
     259
     260     return kTRUE;
     261}
     262
     263Double_t MAlphaFitter::DoOffFit(const TH1D &hon, const TH1D &hof, Bool_t paint)
     264{
     265    if (fSignalFunc!=kThetaSq)
     266        return 0;
     267
     268    // ----------------------------------------------------------------------------
     269
     270    const Int_t bin = hon.GetXaxis()->FindFixBin(fSigInt*0.999);
     271
     272
     273    MAlphaFitter fit(*this);
     274    fit.EnableBackgroundFit();
     275    fit.SetBackgroundFitMin(0);
     276
     277    // produce a histogram containing the off-samples from on-source and
     278    // off-source in the off-source region and the on-data in the source-region
     279    TH1D h(hof);
     280    h.Add(&hon);
     281    h.Scale(0.5);
     282    for (int i=1; i<=bin+3; i++)
     283    {
     284        h.SetBinContent(i, hof.GetBinContent(i));
     285        h.SetBinError(  i, hof.GetBinError(i));
     286    }
     287
     288    // Now fit the off-data
     289    if (!fit.FitOff(h, paint?2:0)) // FIXME: Show fit!
     290        return -1;
     291
     292    // Calculate fit-result
     293    fit.FitResult(h);
     294
     295    // Do a gaussian error propagation to calculated the error of
     296    // the background estimated from the fit
     297    const Double_t ea = fit.GetErrors()[3];
     298    const Double_t eb = fit.GetErrors()[4];
     299    const Double_t a  = fit.GetCoefficients()[3];
     300    const Double_t b  = fit.GetCoefficients()[4];
     301
     302    const Double_t t  = fIntegralMax;
     303
     304    const Double_t ex  = TMath::Exp(t*b);
     305    const Double_t eab = TMath::Exp(a)/b;
     306
     307    const Double_t eA = ex-1;
     308    const Double_t eB = t*ex - eA/b;
     309
     310    const Double_t w  = h.GetXaxis()->GetBinWidth(1);
     311
     312    // Error of estimated background
     313    const Double_t er = TMath::Abs(eab)*TMath::Hypot(eA*ea, eB*eb)/w;
     314
     315    // Calculate arbitrary scale factor from propagated error from the
     316    // condistion: sqrt(alpha*background) = est.background/est.error
     317    // const Double_t bg = hof.Integral(1, bin);
     318    // const Double_t sc = bg * er*er / (fit2.GetEventsBackground()*fit2.GetEventsBackground());
     319    // Assuming that bg and fit2.GetEventsBackground() are rather identical:
     320    const Double_t sc = er*er / fit.GetEventsBackground();
     321
     322    return sc;
     323    /*
     324     cout << MMath::SignificanceLiMaSigned(hon.Integral(1, bin), fit.GetEventsBackground()/sc, sc) << " ";
     325     cout << sc << " ";
     326     cout << fit.fChiSqBg << endl;
     327    */
    220328
    221329    return kTRUE;
    … …  
    230338    MAlphaFitter fit(*this);
    231339    fit.SetPolynomOrder(0);
    232 
    233340    if (alpha<=0 || !fit.Fit(h, paint))
    234341        return kFALSE;
    … …  
    237344    fChiSqBg      = fit.GetChiSqBg();
    238345    fCoefficients = fit.GetCoefficients();
     346    fErrors       = fit.GetErrors();
     347
     348
     349    // ----------------------------------------------------------------------------
     350
     351    const Double_t scale = DoOffFit(hon, hof, paint);
     352    if (scale<0)
     353        return kFALSE;
     354
     355    // ----------------------------------------------------------------------------
    239356
    240357    const Int_t bin = hon.GetXaxis()->FindFixBin(fSigInt*0.999);
    … …  
    321438    f.fCoefficients.Set(fCoefficients.GetSize());
    322439    f.fCoefficients.Reset();
     440    f.fErrors.Set(fCoefficients.GetSize());
     441    f.fErrors.Reset();
    323442
    324443    // Result

  • trunk/MagicSoft/Mars/mhflux/MAlphaFitter.h

    r7710 r7720  
    6666
    6767    TArrayD fCoefficients;      // Fit result
     68    TArrayD fErrors;            // Fit errors
    6869
    6970    // Function
    … …  
    7677    // Minimization strategy
    7778    Strategy_t fStrategy;       // How to calc minimization value
     79
     80    Double_t DoOffFit(const TH1D &hon, const TH1D &hof, Bool_t paint);
     81    Bool_t   FitOff(TH1D &h, Int_t paint);
     82    void     FitResult(const TH1D &h);
    7883
    7984public:
    … …  
    8287        fSigMax(75), fBgMin(45), fBgMax(85), fScaleMin(40), fScaleMax(80),
    8388        fPolynomOrder(2), fFitBackground(kTRUE), fSignalFunc(kGauss),
    84         fCoefficients(3+fPolynomOrder+1),
     89        fCoefficients(3+fPolynomOrder+1), fErrors(3+fPolynomOrder+1),
    8590        fFunc(new TF1("", Form("gaus(0) + pol%d(3)", fPolynomOrder), 0, 90)),
    8691        fScaleMode(kOffRegion), fScaleUser(1), fStrategy(kSignificance)
    … …  
    150155        fFunc->SetName("Dummy");
    151156        gROOT->GetListOfFunctions()->Remove(fFunc);
     157
    152158        fCoefficients.Set(3+fPolynomOrder+1);
    153159        fCoefficients.Reset();
     160
     161        fErrors.Set(3+fPolynomOrder+1);
     162        fErrors.Reset();
    154163    }
    155164    void EnableBackgroundFit(Bool_t b=kTRUE) { fFitBackground=b; }
    … …  
    176185    Double_t GetCoefficient(Int_t i) const { return fCoefficients[i]; }
    177186    const TArrayD &GetCoefficients() const { return fCoefficients; }
     187    const TArrayD &GetErrors() const       { return fErrors; }
    178188    Double_t Eval(Double_t d) const { return fFunc ? fFunc->Eval(d) : 0; }
    179189
    … …  
    232242    Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE);
    233243
    234     ClassDef(MAlphaFitter, 2)
     244    ClassDef(MAlphaFitter, 3)
    235245};
    236246

  • trunk/MagicSoft/Mars/mhflux/MHThetaSqN.cc

    r7718 r7720  
    6767//
    6868MHThetaSqN::MHThetaSqN(const char *name, const char *title)
    69     : MHAlpha(name, title), fDisp(0), fSrcPosCam(0),
     69    : MHAlpha(name, title), fDisp(0), fSrcPosCam(0), fSignificanceCutLevel(1.7),
    7070    fNumBinsSignal(3), fNumBinsTotal(75), fNumOffSourcePos(3), fDoOffCut(kTRUE)
    7171{
    … …  
    175175        fFit.SetScaleUser(1./fNumOffSourcePos);
    176176
     177    fThetaSqCut = fSignificanceCutLevel*fFit.GetSignalIntegralMax()/1.7;
     178
    177179    return kTRUE;
    178180}
    … …  
    215217    const Float_t rad = TMath::TwoPi()/n;
    216218
    217     for (UInt_t i=0; i<n; i++)
    218     {
    219         const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
    220         const Double_t d   = (src-org).Mod2();
    221 
     219/*
     220    for (UInt_t i=0; i<n; i++)
     221    {
    222222        //   off: is in src region   on: is in off regions
    223223        /// if (!fOffData && i==0) || (fOffData && i!=0)
    224224        if ((bool)fOffData ^ (i==0) )
     225            continue;
     226
     227        const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
     228        SetVal((src-org).Mod2());
     229
     230        if (!MHAlpha::Fill(NULL, weight))
     231            return kFALSE;
     232    }
     233    */
     234
     235    // Calculate distance (theta-sq) to all (off-)source regions
     236    TArrayD dist(n);
     237    for (UInt_t i=0; i<n; i++)
     238    {
     239        const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
     240        dist[i] = (src-org).Mod2();
     241    }
     242
     243    // Processing off-data
     244    // Check if event's origin is in the on-regions
     245    if (!fOffData && fDoOffCut && dist[0]<fThetaSqCut)
     246        return kTRUE;
     247
     248    for (UInt_t i=0; i<n; i++)
     249    {
     250        //   off: is in src region   on: is in off regions
     251        /// if (!fOffData && i==0) || (fOffData && i!=0)
     252        if ((bool)fOffData ^ (i==0) )
     253            continue;
     254
     255        Stat_t w = weight;
     256
     257        // Processing on-data
     258        if (fOffData && fDoOffCut)
     259        {
     260            /*
     261             static int cnt=0;
     262             if (dist[1+(cnt++%fNumOffSourcePos)]<fFit.GetSignalIntegralMax())
     263                continue;
     264             */
     265
     266            // Check if event's origin is in one of the off-regions
     267            for (UInt_t j=1; j<n; j++)
     268                if (dist[j]<fThetaSqCut)
     269                {
     270                    w *= (float)(fNumOffSourcePos-1)/fNumOffSourcePos;
     271                    break;
     272                }
     273        }
     274
     275        SetVal(dist[i]);
     276
     277        if (!MHAlpha::Fill(NULL, w))
     278            return kFALSE;
     279    }
     280
     281    /*
     282
     283    // Calculate distance (theta-sq) to all (off-)source regions
     284    TArrayD dist(n);
     285    for (UInt_t i=0; i<n; i++)
     286    {
     287        const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
     288        dist[i] = (src-org).Mod2();
     289    }
     290
     291    for (UInt_t i=0; i<n; i++)
     292    {
     293        //   off: is in src region   on: is in off regions
     294        /// if (!fOffData && i==0) || (fOffData && i!=0)
     295        if ((bool)fOffData ^ (i==0) )
     296            continue;
     297
     298        if (fDoOffCut)
     299        {
     300                UInt_t j;
     301                for (j=0; j<n; j++)
     302                {
     303                    if (i==j)
     304                        continue;
     305
     306                    // *1.5: 99%
     307                    // *1.4: 98%
     308                    // *1.3: 97%
     309                    // *1.2: 96%
     310                    // *1.1: 94%
     311                    // *1.0: 91%
     312
     313                    if (dist[j]<fThetaSqCut)
     314                        break;
     315                }
     316                if (j<n)
     317                    continue;
     318        }
     319
     320        SetVal(dist[i]);
     321
     322        if (!MHAlpha::Fill(NULL, weight))
     323            return kFALSE;
     324    }
     325  */
     326    /*
     327    for (UInt_t i=0; i<n; i++)
     328    {
     329        //   off: is in src region   on: is in off regions
     330        /// if (!fOffData && i==0) || (fOffData && i!=0)
     331        if ((bool)fOffData ^ (i==0) )
     332            continue;
     333
     334        if (fDoOffCut)
     335        {
     336           const TVector2 src1 = const_cast<TVector2&>(src0).Rotate(i*rad+TMath::Pi());
     337           const Double_t d1   = (src1-org).Mod2();
     338
     339           if (d1<fFit.GetSignalIntegralMax())
     340               continue;
     341        }
     342
     343        const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
     344        const Double_t d0  = (src-org).Mod2();
     345
     346        SetVal(d0);
     347
     348        if (!MHAlpha::Fill(NULL, weight))
     349            return kFALSE;
     350    }
     351    */
     352
     353/*
     354    for (UInt_t i=0; i<n; i++)
     355    {
     356        const TVector2 src = const_cast<TVector2&>(src0).Rotate(i*rad);
     357        const Double_t d   = (src-org).Mod2();
     358
     359        //   off: is in src region   on: is in off regions
     360        /// if (!fOffData && i==0) || (fOffData && i!=0)
     361        if ((bool)fOffData ^ (i==0) )
    225362        {
    226363            if (d<fFit.GetSignalIntegralMax() && fDoOffCut)
    … …  
    236373            return kFALSE;
    237374    }
    238 
     375  */
    239376    if (!fOffData)
    240377    {
    … …  
    248385        if (dist<TMath::Sqrt(fFit.GetSignalIntegralMax()))
    249386        {
    250             *fLog << warn << "WARNING - (Off-)source regions start overlapping: ";
     387            *fLog << warn << "WARNING - Source regions overlap: ";
    251388            *fLog << "distance " << dist << " less than theta-sq cut ";
    252389            *fLog << TMath::Sqrt(fFit.GetSignalIntegralMax()) << "!" << endl;
    … …  
    334471        rc = kTRUE;
    335472    }
     473    if (IsEnvDefined(env, prefix, "SignificanceCutLevel", print))
     474    {
     475        SetSignificanceCutLevel(GetEnvValue(env, prefix, "SignificanceCutLevel", fSignificanceCutLevel));
     476        rc = kTRUE;
     477    }
    336478    return rc;
    337479}

  • trunk/MagicSoft/Mars/mhflux/MHThetaSqN.h

    r7718 r7720  
    1717    MSrcPosCam  *fSrcPosCam; //!
    1818
    19     Double_t fMm2Deg;
     19    Double_t fMm2Deg;        //!
     20    Double_t fThetaSqCut;    //!
     21    Double_t fSignificanceCutLevel;
    2022
    2123    UInt_t fNumBinsSignal;
    … …  
    5961    void SetNumOffSourcePos(UInt_t n=3)  { fNumOffSourcePos=TMath::Max(n, 1U); }
    6062
     63    void SetSignificanceCutLevel(Double_t l=1.7)  { fSignificanceCutLevel=l; }
     64
    6165    void SetDoOffCut(Bool_t b=kTRUE)     { fDoOffCut = b; }
    6266
Note: See TracChangeset for help on using the changeset viewer.