Changeset 702 for trunk/MagicSoft/Mars/manalysis/MHillas.cc

Timestamp:

03/27/01 11:34:04 (24 years ago)

Author:

tbretz

Message:

*** empty log message ***

File:

• trunk/MagicSoft/Mars/manalysis/MHillas.cc (modified) (6 diffs)

trunk/MagicSoft/Mars/manalysis/MHillas.cc

@@ -20 +20 @@
 ClassImp(MHillas)
 
-    MHillas::MHillas(const char *name, const char *title) : fEllipse(NULL)
+MHillas::MHillas(const char *name, const char *title) :
+    fAlpha(0), fTheta(0), fWidth(0), fLength(0), fSize(0), fDist(0), fEllipse(NULL)
 {
     *fName  = name  ? name  : "MHillas";
@@ -36 +37 @@
 {
     *fLog << "Hillas Parameter:" << endl;
-    *fLog << " - Alpha  = " << fAlpha  << endl;
+    *fLog << " - Alpha  = " << fabs(fAlpha)  << endl;
     *fLog << " - Width  = " << fWidth  << endl;
     *fLog << " - Length = " << fLength << endl;
+    *fLog << " - Size   = " << fSize   << endl;
     *fLog << " - Dist   = " << fDist   << endl;
+}
+
+void MHillas::Paint(Option_t *)
+{
+    if (!fEllipse)
+        return;
+
+    fEllipse->Paint();
 }
 
@@ -60 +70 @@
     fEllipse->SetLineWidth(2);
     fEllipse->Draw();
+
+    /*
+     This is from TH1
+     TString opt = option;
+   opt.ToLower();
+   if (gPad && !opt.Contains("same")) {
+      //the following statement is necessary in case one attempts to draw
+      //a temporary histogram already in the current pad
+      if (TestBit(kCanDelete)) gPad->GetListOfPrimitives()->Remove(this);
+      gPad->Clear();
+   }
+   AppendPad(opt.Data());
+   */
 }
 
@@ -138 +161 @@
 
     //
-    // resolve four-fold ambiguity of solution
+    // calculate the length of the two axis
     //
     float axis1 =  2.0*c*s*sigmaxy + c*c*sigmaxx + s*s*sigmayy; // [mm^2]
@@ -148 +171 @@
     //
     // check for numerical negatives
+    // (very small number can get negative by chance)
     //
     if (axis1 < 0) axis1=0;
@@ -153 +177 @@
 
     //
-    // check the rotation of the axis (maybe turn by 90ø)
-    //
-    const int   rotation = axis1<axis2;
-
-    fLength = rotation ? sqrt(axis2) : sqrt(axis1);     // [mm]
-    fWidth  = rotation ? sqrt(axis1) : sqrt(axis2);     // [mm]
-
-    fAlpha  = rotation ?
-        fabs(atan((-xmean*c - ymean*s)/(c*ymean - s*xmean))) :
-        fabs(atan(( ymean*c - xmean*s)/(c*xmean + s*ymean))) ;
-
-        // [deg]
-
-    fAlpha *= kRad2Deg;                                  // [deg]
-
-    fDist   = sqrt(xmean*xmean + ymean*ymean);           // [mm]
-
-    fTheta  = atan(ymean/xmean);                         // [rad]
-    if (xmean<0) fTheta += kPI;                          // [rad]
-}
+    // calculate the main Hillas parameters
+    //
+    // fLength, fWidth describes the two axis of the ellipse
+    // fAlpha is the angle between the length-axis and the center
+    //    of the camera
+    // fDist is the distance between the center of the camera and the
+    //    denter of the ellipse
+    //
+    const int rotation = axis1<axis2;
+
+    fLength = rotation ? sqrt(axis2) : sqrt(axis1);  // [mm]
+    fWidth  = rotation ? sqrt(axis1) : sqrt(axis2);  // [mm]
+
+    const float a = c*xmean + s*ymean;
+    const float b = c*ymean - s*xmean;
+
+    fAlpha  = rotation ? atan(a/b) : atan(-b/a);     // [rad]
+    fAlpha *= kRad2Deg;                              // [deg]
+
+    fDist   = sqrt(xmean*xmean + ymean*ymean);       // [mm]
+
+    fTheta  = atan(ymean/xmean);                     // [rad]
+    if (xmean<0) fTheta += kPI;                      // [rad]
+}