Changeset 689 for trunk/MagicSoft/Mars

Timestamp:

03/07/01 11:20:28 (24 years ago)

Author:

harald

Message:

calculation of errors for collection area included!

Location:

trunk/MagicSoft/Mars

Files:

• macros/getCollArea.C (modified) (3 diffs)
• mmontecarlo/MCollArea.cc (modified) (2 diffs)

trunk/MagicSoft/Mars/macros/getCollArea.C

@@ -1 @@
-{{ 
+void getCollArea(char *filename = "/big0/Maggi/CamData/Gamma/gamma_15_on.root" ) 
+{ 
   
   MParList  *parlist  = new MParList() ; 
@@ -14 +15 @@
 
 
-  MReadTree reader("/big0/Maggi/CamData/Gamma/gamma_15_on.root","Events" ) ; 
+  MReadTree reader( filename ,"Events" ) ; 
   tasklist->AddToList( &reader ) ; 
 
@@ -34 +35 @@
   collArea->Draw() ; 
   
-}} 
+}

trunk/MagicSoft/Mars/mmontecarlo/MCollArea.cc

@@ -74 +74 @@
 { 
   //  first of all calculate the efficency
-  fHistSel->Divide( fHistAll) ; 
+  //  do it here by hand to get the right error of efficency
+ 
+  Int_t iBinx = ( (TAxis *) fHistSel->GetXaxis())->GetNbins() ; 
+  Int_t iBiny = ( (TAxis *) fHistSel->GetYaxis())->GetNbins() ; 
+ 
+  Float_t N, Nall ; 
+  Double_t effi, error ; 
+  for (Int_t ix=1; ix<=iBiny; ix++ ) 
+    {
+      for (Int_t iy=1; iy<=iBiny; iy++ ) 
+        { 
+          effi = error = 0. ; 
+
+          N    = fHistSel->GetCellContent(ix, iy) ; 
+          Nall = fHistAll->GetCellContent(ix, iy) ; 
+          
+          if ( Nall > 0 ) { 
+            effi   = N / Nall ; 
+            error = sqrt ( Nall + Nall * N - N * N - N ) / (Nall * Nall ) ; 
+
+            cout << ix << " " << iy 
+                 << " N " << N 
+                 << " Nall " << Nall
+                 << " effi  " << effi 
+                 << " error " << error 
+                 << endl ; 
+
+            fHistSel->SetCellContent(ix, iy, effi) ; 
+            fHistSel->SetCellError(ix, iy, error) ; 
+          } 
+          else 
+            cout << ix << " " << iy << endl ; 
+            
+          
+
+          
+        }
+    } 
   
   // 
@@ -81 +118 @@
   //   
   
-  Int_t iBinx = ( (TAxis *) fHistSel->GetXaxis())->GetNbins() ; 
-  Int_t iBiny = ( (TAxis *) fHistSel->GetYaxis())->GetNbins() ; 
   
-  Double_t  r1, r2, eff, A ; 
+  Double_t  r1, r2, eff, errEff, A, errA, collA ; 
 
   for (Int_t ix=1; ix<=iBiny; ix++ ) 
     {
-      A = 0. ; 
+      A = errA = collA = errEff = 0. ; 
+      
       for (Int_t iy=1; iy<=iBiny; iy++ ) 
         { 
           r1 = ( (TAxis *) fHistSel->GetYaxis())->GetBinLowEdge(iy) ; 
           r2 = ( (TAxis *) fHistSel->GetYaxis())->GetBinLowEdge(iy+1) ; 
-          eff= fHistSel->GetCellContent(ix, iy) ; 
-          A += eff * 3.141592654 * ( r2*r2 - r1*r1 ) ; 
+          A  = 3.141592654 * ( r2*r2 - r1*r1 ) ; 
+          eff    = fHistSel->GetCellContent(ix, iy) ;
+          errEff = fHistSel->GetCellError(ix, iy) ; 
+          collA += eff * A ; 
+          
+          errA  += ( A * A ) * errEff * errEff ; 
+        }
 
-        }
+      errA = sqrt( errA ) ; 
+
       cout << ix << " --> " << A
            << endl ; 
       
-      fHistColl->SetBinContent(ix, A ) ; 
+      fHistColl->SetBinContent(ix, collA ) ; 
+      fHistColl->SetBinError(ix, errA ) ;