Index: /trunk/MagicSoft/Mars/Changelog
===================================================================
--- /trunk/MagicSoft/Mars/Changelog	(revision 4835)
+++ /trunk/MagicSoft/Mars/Changelog	(revision 4836)
@@ -40,5 +40,7 @@
 
    * mhist/MHAlpha.cc:
-     - paint significance
+     - paint significance and othe rinformations
+     - unified fit in Finalize and Paint
+     - replaced significance calculation by Li/Ma
 
    * mhvstime/MHVsTime.[h,cc]:
@@ -53,4 +55,15 @@
    * msignal/MArrivalTime.[h,cc]:
      - added Print()
+
+   * manalysis/MEventRateCalc.[h,cc]:
+     - added the difference in time between two events into the output
+     - made setup more flexible
+
+   * mbase/MContinue.cc:
+     - fixed a bug which caused a problem if MContinue was not in the
+       main tasklist
+
+   * mimage/MHImagePar.[h,cc], mimage/MHNewImagePar.[h,cc]:
+     - added Paint function to support logarithmic y-axis scles
 
 
Index: /trunk/MagicSoft/Mars/mhist/MHAlpha.cc
===================================================================
--- /trunk/MagicSoft/Mars/mhist/MHAlpha.cc	(revision 4835)
+++ /trunk/MagicSoft/Mars/mhist/MHAlpha.cc	(revision 4836)
@@ -101,17 +101,4 @@
 }
 
-// --------------------------------------------------------------------------
-//
-// Calculate Significance as
-// significance = (s-b)/sqrt(s+k*k*b) mit k=s/b
-// 
-Double_t MHAlpha::Significance(Double_t s, Double_t b)
-{
-    const Double_t k = b==0 ? 0 : s/b;
-    const Double_t f = s+k*k*b;
-
-    return f==0 ? 0 : (s-b)/TMath::Sqrt(f);
-}
-
 Bool_t MHAlpha::SetupFill(const MParList *pl)
 {
@@ -158,8 +145,57 @@
 void MHAlpha::Paint(Option_t *opt)
 {
+    DoFit();
+}
+
+// --------------------------------------------------------------------------
+//
+// Draw the histogram
+//
+void MHAlpha::Draw(Option_t *opt)
+{
+    TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
+    pad->SetBorderMode(0);
+
+    fHist.Draw();
+
+    AppendPad("");
+}
+
+// --------------------------------------------------------------------------
+//
+// This is a preliminary implementation of a alpha-fit procedure for
+// all possible source positions. It will be moved into its own
+// more powerfull class soon.
+//
+// The fit function is "gaus(0)+pol2(3)" which is equivalent to:
+//   [0]*exp(-0.5*((x-[1])/[2])^2) + [3] + [4]*x + [5]*x^2
+// or
+//   A*exp(-0.5*((x-mu)/sigma)^2) + a + b*x + c*x^2
+//
+// Parameter [1] is fixed to 0 while the alpha peak should be
+// symmetric around alpha=0.
+//
+// Parameter [4] is fixed to 0 because the first derivative at
+// alpha=0 should be 0, too.
+//
+// In a first step the background is fitted between bgmin and bgmax,
+// while the parameters [0]=0 and [2]=1 are fixed.
+//
+// In a second step the signal region (alpha<sigmax) is fittet using
+// the whole function with parameters [1], [3], [4] and [5] fixed.
+//
+// The number of excess and background events are calculated as
+//   s = int(0, sigint, gaus(0)+pol2(3))
+//   b = int(0, sigint,         pol2(3))
+//
+// The Significance is calculated using the Significance() member
+// function.
+//
+Bool_t MHAlpha::DoFit(Double_t &sig, Bool_t paint)
+{
     Float_t sigint=fSigInt;
     Float_t sigmax=fSigMax;
-    Float_t bgmin =fBgMin;
-    Float_t bgmax =fBgMax;
+    Float_t bgmin=fBgMin;
+    Float_t bgmax=fBgMax;
     Byte_t polynom=fPolynom;
 
@@ -175,11 +211,10 @@
 
     TH1 *h=&fHist;
-    //const Double_t alpha0 = h->GetBinContent(1);
-    //const Double_t alphaw = h->GetXaxis()->GetBinWidth(1);
+    const Double_t alpha0 = h->GetBinContent(1);
+    const Double_t alphaw = h->GetXaxis()->GetBinWidth(1);
 
     // Check for the regios which is not filled...
-    const Double_t alpha0 = h->GetBinContent(1);
     if (alpha0==0)
-        return;
+        return kFALSE; //*fLog << warn << "Histogram empty." << endl;
 
     // First fit a polynom in the off region
@@ -191,12 +226,16 @@
         func.ReleaseParameter(i);
 
-    //h->Fit(&func, "N0Q", "", bgmin, bgmax);
     h->Fit(&func, "N0Q", "", bgmin, bgmax);
-    func.SetRange(0, 90);
-    func.SetLineColor(kRed);
-    func.Paint("same");
-
-    //h4a.Fill(func.GetChisquare());
-    //h5a.Fill(func.GetProb());
+
+    const Double_t chisq1 = func.GetChisquare()/func.GetNDF();
+
+    // ------------------------------------
+    if (paint)
+    {
+        func.SetRange(0, 90);
+        func.SetLineColor(kRed);
+        func.Paint("same");
+    }
+    // ------------------------------------
 
     func.ReleaseParameter(0);
@@ -217,339 +256,51 @@
     func.SetParameter(2, sigmax*0.75);
 
-    //h->Fit(&func, "N0Q", "", 0, sigmax);
     h->Fit(&func, "N0Q", "", 0, sigmax);
-    func.SetLineColor(kGreen);
-    func.Paint("same");
-
-    const Double_t w=fHist.GetXaxis()->GetBinWidth(1);
-
-    const Double_t s   = func.Integral(0, sigint)/w;
+
+    const Double_t chisq2    = func.GetChisquare()/func.GetNDF();
+    const Double_t sigmagaus = func.GetParameter(2);
+
+    // ------------------------------------
+    if (paint)
+    {
+        func.SetLineColor(kGreen);
+        func.Paint("same");
+    }
+    // ------------------------------------
+
+    const Double_t s   = func.Integral(0, sigint)/alphaw;
     func.SetParameter(0, 0);
     func.SetParameter(2, 1);
-    const Double_t b   = func.Integral(0, sigint)/w;
-    const Double_t sig = MMath::SignificanceLiMaSigned(s, b);
-
-    TLatex text;
-    text.PaintLatex(35, fHist.GetMaximum()*1.1, 0, 0.05, Form("\\sigma=%.1f E=%d", sig, (int)(s-b)));
-}
-
-// --------------------------------------------------------------------------
-//
-// Draw the histogram
-//
-void MHAlpha::Draw(Option_t *opt)
-{
-    TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
-    pad->SetBorderMode(0);
-
-    fHist.Draw();
-
-    AppendPad("");
-}
-
-// --------------------------------------------------------------------------
-//
-// This is a preliminary implementation of a alpha-fit procedure for
-// all possible source positions. It will be moved into its own
-// more powerfull class soon.
-//
-// The fit function is "gaus(0)+pol2(3)" which is equivalent to:
-//   [0]*exp(-0.5*((x-[1])/[2])^2) + [3] + [4]*x + [5]*x^2
-// or
-//   A*exp(-0.5*((x-mu)/sigma)^2) + a + b*x + c*x^2
-//
-// Parameter [1] is fixed to 0 while the alpha peak should be
-// symmetric around alpha=0.
-//
-// Parameter [4] is fixed to 0 because the first derivative at
-// alpha=0 should be 0, too.
-//
-// In a first step the background is fitted between bgmin and bgmax,
-// while the parameters [0]=0 and [2]=1 are fixed.
-//
-// In a second step the signal region (alpha<sigmax) is fittet using
-// the whole function with parameters [1], [3], [4] and [5] fixed.
-//
-// The number of excess and background events are calculated as
-//   s = int(0, sigint, gaus(0)+pol2(3))
-//   b = int(0, sigint,         pol2(3))
-//
-// The Significance is calculated using the Significance() member
-// function.
-//
+    const Double_t b   = func.Integral(0, sigint)/alphaw;
+
+    sig = MMath::SignificanceLiMaSigned(s, b);
+
+    // ------------------------------------
+    if (paint)
+    {
+        TLatex text;
+        text.PaintLatex(9, fHist.GetMaximum()*1.11, 0, 0.04,
+                        Form("\\sigma_{Li/Ma}=%.1f  \\omega=%.1f\\circ  E=%d  (\\chi_{b}^{2}/N=%.1f  \\chi_{s}^{2}/N=%.1f)",
+                             sig, sigmagaus, (int)(s-b), chisq1, chisq2));
+    }
+    // ------------------------------------
+
+    return kTRUE;
+}
+
 Bool_t MHAlpha::Finalize()
 {
-    //*fLog << dbg << "Fitting..." << endl;
-
-    Float_t sigint=fSigInt;
-    Float_t sigmax=fSigMax;
-    Float_t bgmin=fBgMin;
-    Float_t bgmax=fBgMax;
-    Byte_t polynom=fPolynom;
-
-    //                      xmin, xmax, npar
-    //TF1 func("MyFunc", fcn, 0, 90, 6);
-    // Implementing the function yourself is only about 5% faster
-    TF1 func("", Form("gaus(0) + pol%d(3)", polynom), 0, 90);
-    TArrayD maxpar(func.GetNpar());
-
-    //*fLog << "Setup pars..." << endl;
-
-    func.FixParameter(1, 0);
-    func.FixParameter(4, 0);
-    func.SetParLimits(2, 0, 90);
-    func.SetParLimits(3, -1, 1);
-
-    //    TStopwatch clk;
-    //    clk.Start();
-   /*
-    *fLog << inf;
-    *fLog << "Signal fit:     alpha < " << sigmax << endl;
-    *fLog << "Integration:    alpha < " << sigint << endl;
-    *fLog << "Background fit: " << bgmin << " < alpha < " << bgmax << endl;
-    *fLog << "Polynom order:  " << (int)polynom << endl;
-    *fLog << "Fitting False Source Plot..." << flush;
-    */
-    TH1 *h=&fHist;
-    const Double_t alpha0 = h->GetBinContent(1);
-    const Double_t alphaw = h->GetXaxis()->GetBinWidth(1);
-
-    // Check for the regios which is not filled...
-    if (alpha0==0)
+    Double_t sig = 0;
+
+    if (!DoFit(sig))
     {
         *fLog << warn << "Histogram empty." << endl;
-        return kTRUE;
-    }
-
-    // First fit a polynom in the off region
-    func.FixParameter(0, 0);
-    func.FixParameter(2, 1);
-    func.ReleaseParameter(3);
-
-    //*fLog << "Release pars..." << endl;
-
-    for (int i=5; i<func.GetNpar(); i++)
-        func.ReleaseParameter(i);
-
-    //*fLog << dbg << "Fit 1" << endl;
-    //h->Fit(&func, "N0Q", "", bgmin, bgmax);
-    h->Fit(&func, "N0Q", "", bgmin, bgmax);
-    //*fLog << dbg << ix << "/" << iy << ":  " << func.GetParameter(3) << "    " << func.GetParameter(4) << endl;
-
-    //h4a.Fill(func.GetChisquare());
-    //h5a.Fill(func.GetProb());
-
-    //func.SetParLimits(0, 0.5*h->GetBinContent(1), 1.5*h->GetBinContent(1));
-    //func.SetParLimits(2, 5, 90);
-
-    //*fLog << dbg << "Change params..." << endl;
-
-    func.ReleaseParameter(0);
-    //func.ReleaseParameter(1);
-    func.ReleaseParameter(2);
-    func.FixParameter(3, func.GetParameter(3));
-    for (int i=5; i<func.GetNpar(); i++)
-        func.FixParameter(i, func.GetParameter(i));
-
-    // Do not allow signals smaller than the background
-    const Double_t A  = alpha0-func.GetParameter(3);
-    const Double_t dA = TMath::Abs(A);
-    func.SetParLimits(0, -dA*4, dA*4);
-    func.SetParLimits(2, 0, 90);
-
-    // Now fit a gaus in the on region on top of the polynom
-    func.SetParameter(0, A);
-    func.SetParameter(2, sigmax*0.75);
-
-    //*fLog << dbg << "Fit 2..." << endl;
-    //h->Fit(&func, "N0Q", "", 0, sigmax);
-    h->Fit(&func, "N0Q", "", 0, sigmax);
-    //*fLog << dbg << "     " << func.GetParameter(0) << "    " << func.GetParameter(1) << "    " << func.GetParameter(2) << endl;
-
-    //*fLog << dbg << "Copy par..." << endl;
-    TArrayD p(func.GetNpar(), func.GetParameters());
-
-    // Fill results into some histograms
-    /*
-    h0a.Fill(p[0]);
-    h0b.Fill(p[3]);
-    h1.Fill(p[1]);
-    h2.Fill(p[2]);
-    if (polynom>1)
-        h3.Fill(p[5]);
-    h4b.Fill(func.GetChisquare());
-    //h5b.Fill(func.GetProb());
-    */
-    // Implementing the integral as analytical function
-    // gives the same result in the order of 10e-5
-    // and it is not faster at all...
-
-    //*fLog << dbg << "Int1 par..." << endl;
-    const Double_t s = func.Integral(0, sigint)/alphaw;
-
-    func.SetParameter(0, 0);
-    //func.SetParameter(2, 1);
-
-    //*fLog << dbg << "Int2 par..." << endl;
-    const Double_t b   = func.Integral(0, sigint)/alphaw;
-    const Double_t sig = Significance(s, b);
-
-    //*fLog << dbg << "Set Val "<<sig<<" ..." << fResult << endl;
+        return kFALSE;
+    }
+
     if (fResult)
         fResult->SetVal(sig);
 
-    *fLog << inf << "Found Significance: " << sig << "  (width=";
-    *fLog << func.GetParameter(2) << "°, excess=" << s-b << ")" << endl;
-
-    //*fLog << dbg << "Done." << endl;
     return kTRUE;
-    /*
-    if (sig!=0)
-        h6.Fill(sig);
-
-    if (sig<maxs)
-    {
-        maxs = sig;
-        maxx = ix;
-        maxy = iy;
-        maxpar = p;
-        }
-       */
-/*
-
-    *fLog << inf << "done." << endl;
-
-    h0a.GetXaxis()->SetRangeUser(0, maxalpha0*1.5);
-    h0b.GetXaxis()->SetRangeUser(0, maxalpha0*1.5);
-
-    //hists->SetMinimum(GetMinimumGT(*hists));
-    histb->SetMinimum(GetMinimumGT(*histb));
-
-//    clk.Stop();
-//    clk.Print();
-    TCanvas *c=new TCanvas;
-
-    c->Divide(3,2, 0, 0);
-    c->cd(1);
-    gPad->SetBorderMode(0);
-    hists->Draw("colz");
-    hists->SetBit(kCanDelete);
-    c->cd(2);
-    gPad->SetBorderMode(0);
-    hist->Draw("colz");
-    hist->SetBit(kCanDelete);
-    c->cd(3);
-    gPad->SetBorderMode(0);
-    histb->Draw("colz");
-    histb->SetBit(kCanDelete);
-    c->cd(4);
-    gPad->Divide(1,3, 0, 0);
-    TVirtualPad *p=gPad;
-    p->SetBorderMode(0);
-    p->cd(1);
-    gPad->SetBorderMode(0);
-    h0b.DrawCopy();
-    h0a.DrawCopy("same");
-    p->cd(2);
-    gPad->SetBorderMode(0);
-    h3.DrawCopy();
-    p->cd(3);
-    gPad->SetBorderMode(0);
-    h2.DrawCopy();
-    c->cd(6);
-    gPad->Divide(1,2, 0, 0);
-    TVirtualPad *q=gPad;
-    q->SetBorderMode(0);
-    q->cd(1);
-    gPad->SetBorderMode(0);
-    gPad->SetBorderMode(0);
-    h4b.DrawCopy();
-    h4a.DrawCopy("same");
-    h6.DrawCopy("same");
-    q->cd(2);
-    gPad->SetBorderMode(0);
-    //h5b.DrawCopy();
-    //h5a.DrawCopy("same");
-    c->cd(5);
-    gPad->SetBorderMode(0);
-    if (maxx>0 && maxy>0)
-    {
-        const char *title = Form(" \\alpha for x=%.2f y=%.2f (\\sigma_{max}=%.1f) ",
-                                 fHist.GetXaxis()->GetBinCenter(maxx),
-                                 fHist.GetYaxis()->GetBinCenter(maxy), maxs);
-
-        TH1 *result = fHist.ProjectionZ("AlphaFit", maxx, maxx, maxy, maxy);
-        result->SetDirectory(NULL);
-        result->SetNameTitle("Result \\alpha", title);
-        result->SetBit(kCanDelete);
-        result->SetXTitle("\\alpha [\\circ]");
-        result->SetYTitle("Counts");
-        result->Draw();
-
-        TF1 f1("", func.GetExpFormula(), 0, 90);
-        TF1 f2("", func.GetExpFormula(), 0, 90);
-        f1.SetParameters(maxpar.GetArray());
-        f2.SetParameters(maxpar.GetArray());
-        f2.FixParameter(0, 0);
-        f2.FixParameter(1, 0);
-        f2.FixParameter(2, 1);
-        f1.SetLineColor(kGreen);
-        f2.SetLineColor(kRed);
-
-        f2.DrawCopy("same");
-        f1.DrawCopy("same");
-
-        TPaveText *leg = new TPaveText(0.35, 0.10, 0.90, 0.35, "brNDC");
-        leg->SetBorderSize(1);
-        leg->SetTextSize(0.04);
-        leg->AddText(0.5, 0.82, Form("A * exp(-(\\frac{x-\\mu}{\\sigma})^{2}/2) + pol%d", polynom))->SetTextAlign(22);
-        //leg->AddText(0.5, 0.82, "A * exp(-(\\frac{x-\\mu}{\\sigma})^{2}/2) + b*x^{2} + a")->SetTextAlign(22);
-        leg->AddLine(0, 0.65, 0, 0.65);
-        leg->AddText(0.06, 0.54, Form("A=%.2f", maxpar[0]))->SetTextAlign(11);
-        leg->AddText(0.06, 0.34, Form("\\sigma=%.2f", maxpar[2]))->SetTextAlign(11);
-        leg->AddText(0.06, 0.14, Form("\\mu=%.2f (fix)", maxpar[1]))->SetTextAlign(11);
-        leg->AddText(0.60, 0.54, Form("a=%.2f", maxpar[3]))->SetTextAlign(11);
-        leg->AddText(0.60, 0.34, Form("b=%.2f (fix)", maxpar[4]))->SetTextAlign(11);
-        if (polynom>1)
-            leg->AddText(0.60, 0.14, Form("c=%.2f", maxpar[5]))->SetTextAlign(11);
-        leg->SetBit(kCanDelete);
-        leg->Draw();
-
-        q->cd(2);
-
-        TGraph *g = new TGraph;
-        g->SetPoint(0, 0, 0);
-
-        Int_t max=0;
-        Float_t maxsig=0;
-        for (int i=1; i<89; i++)
-        {
-            const Double_t s = f1.Integral(0, (float)i);
-            const Double_t b = f2.Integral(0, (float)i);
-
-            const Double_t sig = Significance(s, b);
-
-            g->SetPoint(g->GetN(), i, sig);
-
-            if (sig>maxsig)
-            {
-                max = i;
-                maxsig = sig;
-            }
-        }
-
-        g->SetNameTitle("SigVs\\alpha", "Significance vs \\alpha");
-        g->GetHistogram()->SetXTitle("\\alpha_{0} [\\circ]");
-        g->GetHistogram()->SetYTitle("Significance");
-        g->SetBit(kCanDelete);
-        g->Draw("AP");
-
-        leg = new TPaveText(0.35, 0.10, 0.90, 0.25, "brNDC");
-        leg->SetBorderSize(1);
-        leg->SetTextSize(0.1);
-        leg->AddText(Form("\\sigma_{max}=%.1f at \\alpha_{max}=%d\\circ", maxsig, max));
-        leg->SetBit(kCanDelete);
-        leg->Draw();
-    }*/
 }
 
Index: /trunk/MagicSoft/Mars/mhist/MHAlpha.h
===================================================================
--- /trunk/MagicSoft/Mars/mhist/MHAlpha.h	(revision 4835)
+++ /trunk/MagicSoft/Mars/mhist/MHAlpha.h	(revision 4836)
@@ -33,4 +33,8 @@
     Int_t fMap;               //!
 
+    Bool_t DoFit(Double_t &, Bool_t);
+    Bool_t DoFit(Double_t &d) { return DoFit(d, kFALSE); }
+    void   DoFit() { Double_t d; DoFit(d, kTRUE); }
+
 public:
     MHAlpha(const char *name=NULL, const char *title=NULL);
@@ -46,6 +50,4 @@
     void StopMapping();
 
-    static Double_t Significance(Double_t s, Double_t b);
-
     ClassDef(MHAlpha, 1) // Alpha-Plot which is fitted online
 };
