Changeset 3105 for trunk/Ohana/src/opihi/mana/fitcontour.c

Timestamp:

Jan 28, 2005, 5:17:26 AM (21 years ago)

Author:

eugene

Message:

minor changes

File:

• trunk/Ohana/src/opihi/mana/fitcontour.c (modified) (3 diffs)

trunk/Ohana/src/opihi/mana/fitcontour.c

@@ -36 +36 @@
     r = hypot (x, y);
 
-    /* 
-    t = 2*atan2(y,x);
-    sn = sin(t);
-    cs = cos(t);
-    */
+    /* calculate sin(2t), cos(2t) using 1/2 angle tri relationship above */
     sn1 = y / r;
     cs1 = x / r;
@@ -63 +59 @@
   gaussj (C, NTERM, B, 1);
   
+  /** this is somewhat weak: if the object is too elongated, Rmin can be < 0 **/
   dR = hypot (B[1][0], B[2][0]);
   Rmaj = B[0][0] + dR;
@@ -101 +98 @@
 }
 
+
+/* this routine fits a single sine and cosine to the value of dr as a function of the angle around the central
+ * point, theta.  this is NOT an exact fit, an is increasingly in error for more eccentric ellipses.  however,
+ * it does a reasonable job of approximating the value of the major and minor axes, and it gets the angle of
+ * orientation right at well.  We then assume the values of major and minor axis and angle are correct to get
+ * the parameters for the elliptical gaussian fit:
+ 
+ z = (x^2) / (2 Rx^2) + (y^2) / (2 Ry^2) + Rxy x y
+
+ the functions above give the correct relationship between these:
+
+   R1 = SQ(Rmaj) + SQ(Rmin);
+   R2 = SQ(Rmaj) - SQ(Rmin);
+   R3 = Rmaj*Rmin;
+
+   Rx = R3 / sqrt (R1 - R2*cs);
+   Ry = R3 / sqrt (R1 + R2*cs);
+   Rxy = -sn*R2 / SQ(R3);
+
+ to derive these relationships, write the equation for an ellipse with major axis in the x-dir and minor in
+ the y-dir:
+
+ z = (x^2) / (2 Rmaj^2) + (y^2) / (2 Rmin^2)
+
+ apply the rotation matrix:
+
+ (x)' = (cos, -sin)(x,y)
+ (y)' = (sin, +cos)(x,y)  (modulo the sign on the sine)
+
+ then group the x^2 terms, the y^2 terms, and the xy terms to find the three coeffs.
+
+ NOTE: I spent a while having trouble deriving this from the polar form of an ellipse:
+
+ x = Rmaj*cos(theta)
+ y = Rmin*sin(theta)
+
+ it turns out that 'theta' above is NOT the angle (0,1)-(0,0)-(x,y). rather, it should be viewed as a
+ paraterization of the ellipse.  
+
+*/