Index: /trunk/MagicSoft/Mars/Changelog
===================================================================
--- /trunk/MagicSoft/Mars/Changelog	(revision 2989)
+++ /trunk/MagicSoft/Mars/Changelog	(revision 2990)
@@ -46,4 +46,13 @@
    * mraw/MRawRunHeader.[h,cc]:
      - initialize fRunType to new enum kRTNone
+
+   * mfilter/MFGeomag.cc:
+     - set fResult to kFALSE at the beginning of Process
+     - according to this removed setting fResult=kFALSE before return
+     - replaced some float by Float_t
+     - added some const-qualifiers
+     - replaced (rig<0.5/pr*R*(r2-r1)+r1) by (rig-r1)*pr<rnd
+       with rnd = R * (r2-r1)/2 to make sure that we cannot
+       devide by 0
 
 
Index: /trunk/MagicSoft/Mars/mfilter/MFGeomag.cc
===================================================================
--- /trunk/MagicSoft/Mars/mfilter/MFGeomag.cc	(revision 2989)
+++ /trunk/MagicSoft/Mars/mfilter/MFGeomag.cc	(revision 2990)
@@ -145,19 +145,18 @@
 //
 Int_t MFGeomag::Process()
-   {
+{
+    fResult = kFALSE;
+
     const Float_t en =  fMcEvt->GetEnergy();       // for rigidity (set P = E)
-    float rig = en;
+    Float_t rig = en;
     const Float_t az =  fMcEvt->GetTelescopePhi(); // charge theta phi are entries in table
     const Float_t th =  fMcEvt->GetTelescopeTheta();
 
-    Int_t indadd=0;         //first part of table (positive particles)
+    Int_t indadd=0;              //first part of table (positive particles)
     switch (fMcEvt->GetPartId())
     {
     case kGAMMA:
-        if (!fGammaElectron)         //accept gammas if not set to electrons
-	  {
-	    fResult = kFALSE;
-	    return kTRUE;
-	  }
+        if (!fGammaElectron)     //accept gammas if not set to electrons
+            return kTRUE;
         indadd = 1152;           //second part of table (negative particles)
         break;
@@ -168,11 +167,9 @@
 
     case kPROTON:                //protons
+    case kPOSITRON:              //positrons
         break;
 
     case kELECTRON:              //electrons
         indadd = 1152;           //second part of table (negative particles)
-        break;          
-
-    case kPOSITRON:              //positrons
         break;
 
@@ -181,28 +178,31 @@
         return kFALSE;
     }
-             // here is the cut for charged particles using the table
-      int it=(int)(th*11.459156);    // table steps are in 5 deg = 1/11.459 rads
-      int ia=(int)(az*11.459156);
-      ia = (ia+36) % 72;             // azimuth definitions differ by 180 deg
-
-      float r1=fRigMin[72*it+ia+indadd];
-      if (rig<=r1) {
-          fResult=kTRUE;    // reject
-          return kTRUE;
-      }
-
-      float r2=fRigMax[72*it+ia+indadd];
-      if (rig>=r2) {
-          fResult=kFALSE;   // accept
-          return kTRUE;
-      }
-
-      float R = gRandom->Rndm(0);        // accept if above intermediate threshold
-      float pr=fProb[72*it+ia+indadd];
-      fResult = kFALSE;
-
-      if (rig < 0.5/pr*R*(r2-r1) + r1)
-          fResult = kTRUE;               // pretty good approximation
-
-      return kTRUE;
-   }
+
+    // here is the cut for charged particles using the table
+
+    int it=(int)(th*11.459156);    // table steps are in 5 deg = 1/11.459 rads
+    int ia=(int)(az*11.459156);
+
+    ia = (ia+36) % 72;             // azimuth definitions differ by 180 deg
+
+    const Float_t r1=fRigMin[72*it+ia+indadd];
+    if (rig<=r1)
+    {
+        fResult=kTRUE;   // reject
+        return kTRUE;
+    }
+
+    const Float_t r2=fRigMax[72*it+ia+indadd];
+    if (rig>=r2)
+        return kTRUE;   // accept
+
+    const Float_t pr=fProb[72*it+ia+indadd];
+
+    // accept if above intermediate threshold
+    const Float_t rnd = (r2-r1)/2 * gRandom->Rndm(0);
+
+    if ((rig-r1)*pr < rnd)
+        fResult = kTRUE;                // pretty good approximation
+
+    return kTRUE;
+}
