Index: trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.cc
===================================================================
--- trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.cc	(revision 853)
+++ trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.cc	(revision 853)
@@ -0,0 +1,150 @@
+/* ======================================================================== *\
+!
+! *
+! * This file is part of MARS, the MAGIC Analysis and Reconstruction
+! * Software. It is distributed to you in the hope that it can be a useful
+! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
+! * It is distributed WITHOUT ANY WARRANTY.
+! *
+! * Permission to use, copy, modify and distribute this software and its
+! * documentation for any purpose is hereby granted without fee,
+! * provided that the above copyright notice appear in all copies and
+! * that both that copyright notice and this permission notice appear
+! * in supporting documentation. It is provided "as is" without express
+! * or implied warranty.
+! *
+!
+!
+!   Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de)
+!              Thomas Bretz  12/2000 (tbretz@uni-sw.gwdg.de)
+!
+!   Copyright: MAGIC Software Development, 2000-2001
+!
+!
+\* ======================================================================== */
+
+#include "MHMcCollectionArea.h" 
+
+#include <MLog.h>
+#include <TH2.h> 
+
+ClassImp(MHMcCollectionArea)
+
+MHMcCollectionArea::MHMcCollectionArea(const char *name, const char *title)
+{ 
+    //
+    //   default constructor
+    //
+
+    //   initialize the histogram for the distribution r vs E
+    //
+    //   we set the energy range from 1 Gev to 10000 GeV (in log 5 orders
+    //   of magnitude) and for each order we take 10 subdivision --> 50 xbins
+    //
+    //   we set the radius range from 0 m to 500 m with 10 m bin --> 50 ybins
+
+  
+    *fName  = name  ? name  : "MHMcCollectionArea";
+    *fTitle = title ? title : "Data to Calculate Coll-Area";
+
+
+    fHistAll = new TH2D("collAreaAll", "all showers - Radius vs log(E) distribution",
+                        50, 0., 5.,
+                        50, 0., 500.);
+
+    fHistSel = new TH2D("collAreaSel", "selected showers - Radius vs log(E) distribution",
+                        50, 0., 5.,
+                        50, 0., 500.);
+
+    fHistCol = new TH1D("collArea", "Collection Area",
+                        50, 0., 5.);
+}
+
+MHMcCollectionArea::~MHMcCollectionArea()
+{
+    delete fHistAll;
+    delete fHistSel;
+    delete fHistCol;
+}
+
+void MHMcCollectionArea::FillAll(Float_t log10E, Float_t radius)
+{
+    fHistAll->Fill(log10E, radius);
+}
+
+void MHMcCollectionArea::FillSel(Float_t log10E, Float_t radius)
+{
+    fHistSel->Fill(log10E, radius);
+}
+
+void MHMcCollectionArea::DrawAll()
+{
+    fHistAll->Draw();
+}
+
+void MHMcCollectionArea::DrawSel()
+{
+    fHistSel->Draw();
+}
+
+void MHMcCollectionArea::Draw(Option_t* option)
+{
+    fHistCol->Draw(option);
+}
+
+void MHMcCollectionArea::CalcEfficiency()
+{
+    // Description!
+
+    //
+    //  first of all calculate the efficency
+    //  do it here by hand to get the right error of efficency
+    //
+    const Int_t nbinx = fHistSel->GetXaxis()->GetNbins();
+    const Int_t nbiny = fHistSel->GetYaxis()->GetNbins();
+
+    for (Int_t ix=1; ix<=nbinx; ix++)
+    {
+        for (Int_t iy=1; iy<=nbiny; iy++)
+        {
+            const Float_t Ns = fHistSel->GetCellContent(ix, iy);
+            const Float_t Na = fHistAll->GetCellContent(ix, iy);
+
+            if (Na <= 0)
+                continue;
+
+            const Double_t eff = Ns/Na;
+            const Double_t err = sqrt(Na + Na*Ns - Ns*Ns - Ns) / (Na*Na);
+
+            fHistSel->SetCellContent(ix, iy, eff);
+            fHistSel->SetCellError(ix, iy, err);
+        }
+    }
+
+    //
+    //  now calculate the Collection area for different
+    //  energies
+    //
+    for (Int_t ix=1; ix<=nbinx; ix++)
+    {
+        Double_t errA = 0;
+        Double_t colA = 0;
+
+        for (Int_t iy=1; iy<=nbiny; iy++)
+        {
+            const Double_t r1  = fHistSel->GetYaxis()->GetBinLowEdge(iy);
+            const Double_t r2  = fHistSel->GetYaxis()->GetBinLowEdge(iy+1);
+
+            const Double_t A   = kPI * (r2*r2 - r1*r1);
+
+            const Double_t eff = fHistSel->GetCellContent(ix, iy);
+            const Double_t err = fHistSel->GetCellError(ix, iy);
+
+            colA += eff*A;
+            errA += A*A * err*err;
+        }
+
+        fHistCol->SetBinContent(ix, colA);
+        fHistCol->SetBinError(ix, sqrt(errA));
+    }
+}
Index: trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.h
===================================================================
--- trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.h	(revision 853)
+++ trunk/MagicSoft/Mars/mhist/MHMcCollectionArea.h	(revision 853)
@@ -0,0 +1,42 @@
+#ifndef MHMCCOLLECTIONAREA_H
+#define MHMCCOLLECTIONAREA_H
+
+#ifndef MAGIC_H
+#include "MAGIC.h"
+#endif
+#ifndef MPARCONTAINER_H
+#include "MParContainer.h"
+#endif
+
+//
+// because of some strange reason this cannot be put into MonteCarloIncl
+//
+#ifndef ROOT_TH1
+#include <TH1.h>
+#endif
+
+class TH2D;
+
+class MHMcCollectionArea : public MParContainer
+{
+private:
+    TH2D  *fHistAll ; //! all simulated showers
+    TH2D  *fHistSel ; //! the selected showers
+    TH1D  *fHistCol ; //  the collection area
+
+public:
+
+    MHMcCollectionArea(const char *name=NULL, const char *title=NULL) ;
+    ~MHMcCollectionArea() ;
+
+    void FillAll(Float_t log10E, Float_t radius) ;
+    void FillSel(Float_t log10E, Float_t radius) ;
+    void DrawAll() ;
+    void DrawSel() ;
+    void Draw(Option_t* option = "") ;
+    void CalcEfficiency() ;
+
+    ClassDef(MHMcCollectionArea, 1)  //  Data Container to calculate Collection Area
+};
+
+#endif
