Changeset 13063 for trunk/psLib/src/math/psEllipse.c

Timestamp:

Apr 27, 2007, 12:13:33 PM (19 years ago)

Author:

magnier

Message:

adding polarization terms

File:

• trunk/psLib/src/math/psEllipse.c (modified) (1 diff)

trunk/psLib/src/math/psEllipse.c

@@ -7 +7 @@
 #include "psConstants.h"
 #include "psEllipse.h"
+
+// ellipse rotation (major, minor, theta) -> (x2, y2, xy)
+psEllipsePol psEllipseAxesToPol(psEllipseAxes axes)
+{
+    psEllipsePol pol;
+
+    pol.e0 = PS_SQR(axes.major) + PS_SQR(axes.minor);
+    double ds = PS_SQR(axes.major) - PS_SQR(axes.minor);
+
+    pol.e1 = ds*cos(2*axes.theta);
+    pol.e2 = ds*sin(2*axes.theta);
+
+    assert (isfinite(pol.e0));
+    assert (isfinite(pol.e1));
+    assert (isfinite(pol.e2));
+
+    return pol;
+}
+
+// ellipse rotation (major, minor, theta) -> (x2, y2, xy)
+psEllipsePol psEllipseShapeToPol(psEllipseShape shape)
+{
+    psEllipsePol pol;
+
+    double r = 1.0 / (1.0 - PS_SQR(shape.sx)*PS_SQR(shape.sy)*PS_SQR(shape.sxy));
+
+    pol.e0 = r*(PS_SQR(shape.sx) + PS_SQR(shape.sy));
+    pol.e1 = r*(PS_SQR(shape.sx) - PS_SQR(shape.sy));
+    pol.e2 = r*(2.0*PS_SQR(shape.sx)*PS_SQR(shape.sy)*shape.sxy);
+
+    assert (isfinite(pol.e0));
+    assert (isfinite(pol.e1));
+    assert (isfinite(pol.e2));
+
+    return pol;
+}
+
+// ellipse rotation (major, minor, theta) -> (x2, y2, xy)
+psEllipseAxes psEllipsePolToAxes(psEllipsePol pol, float minMinorAxis)
+{
+    psEllipseAxes axes;
+
+    axes.theta = +0.5 * atan2 (pol.e2, pol.e1); // theta in radians
+
+    double cs = cos(2*axes.theta);
+    double sn = sin(2*axes.theta);
+    double ds = 0;
+
+    if ((cs > 0.707) || (cs < -0.707)) {
+        ds = pol.e1 / cs;
+    } else {
+        ds = pol.e2 / sn;
+    }
+
+    float LIMIT = PS_SQR(minMinorAxis);
+    if (2.0*(pol.e0 - ds) < LIMIT) {
+        axes.major = sqrt(pol.e0 - LIMIT);
+        axes.minor = sqrt(LIMIT);
+    } else {
+        axes.major = sqrt(0.5*(pol.e0 + ds));
+        axes.minor = sqrt(0.5*(pol.e0 - ds));
+    }
+
+    assert (isfinite(axes.major));
+    assert (isfinite(axes.minor));
+    assert (isfinite(axes.theta));
+
+    return axes;
+}
 
 // ellipse rotation (major, minor, theta) -> (x2, y2, xy)