Changeset 3714

Timestamp:

Apr 18, 2005, 6:16:02 PM (21 years ago)

Author:

gusciora

Message:

Created tests for inverting linear psPlaneTransforms. Fixed a bug with the
linear transform inversion code.

Location:

trunk/psLib

Files:

• src/astro/psCoord.c (modified) (4 diffs)
• src/astronomy/psCoord.c (modified) (4 diffs)
• src/dataManip/psFunctions.c (modified) (2 diffs)
• src/image/psImageStats.c (modified) (2 diffs)
• src/imageops/psImageStats.c (modified) (2 diffs)
• src/math/psPolynomial.c (modified) (2 diffs)
• src/math/psSpline.c (modified) (2 diffs)
• test/astronomy/tst_psAstrometry01.c (modified) (4 diffs)
• test/image/Makefile.am (modified) (2 diffs)

trunk/psLib/src/astro/psCoord.c

@@ -10 +10 @@
 *  @author GLG, MHPCC
 *
-*  @version $Revision: 1.63 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2005-04-19 01:22:59 $
+*  @version $Revision: 1.64 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2005-04-19 04:16:02 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -84 +84 @@
  X1 = (-D + ((F*A)/C)) / (E - ((F*B)/C)) +
       (X2 * -((F/C) / (E - ((F*B)/C)))) +
-             (Y2 * (1.0 / (E - ((F*B)/C))));
-        Y1 = (-D + ((E*A)/B))/(F - ((C*E)/B)) +
-             (X2 * -((E/B) / (F - ((C*E)/B)))) +
-             (Y2 * (1.0 / (F - ((C*E)/B))));
+      (Y2 * (1.0 / (E - ((F*B)/C))));
+ Y1 = (-D + ((E*A)/B))/(F - ((C*E)/B)) +
+      (X2 * -((E/B) / (F - ((C*E)/B)))) +
+      (Y2 * (1.0 / (F - ((C*E)/B))));
  
 XXX: Since thre is now a general psPlaneTransformInvert() function, we
 should rename this.
+ 
+X1, Y1 = (2, 3)
+ 
+    X2 = 1 + 4 + 9   = 14
+    Y2 = 4 + 10 + 18 = 32
+ 
+in inverse
+ 
+Y1 = (32 - ((5/2) * 14) - 4 + ((5)/2)) / (6 - ((15)/2));
+Y1 = (32 - 35 - 4 + (5/2)) / (6 - 15/2);
+Y1 = (-7 + (5/2)) / (12/2 - 15/2);
+Y1 = (-7 + (5/2)) / (-3/2);
+Y1 = (-14/2 + 5/2) / (-3/2);
+Y1 = (-9/2) / (-3/2);
+Y1 = (9/2) / (3/2);
+Y1 = 3
+ 
+X1 = (Y2 - ((F/C) * X2) - D + ((F*A)/C)) / (E - ((F*B)/C));
+X1 = (32 - ((6/3) * 14) - 4 + ((6)/3)) / (5 - ((12)/3));
+X1 = (32 - (2 * 14) - 4 + 2 / (5 - 4);
+X1 = (32 - 28 - 4 + 2) / 1;
+X1 = 2;
  
  *****************************************************************************/
@@ -125 +147 @@
     psPlaneTransform *out = psPlaneTransformAlloc(2, 2);
 
-    out->x->coeff[0][0] = -D + ((F*A)/C) / (E - ((F*B)/C));
+    out->x->coeff[0][0] = (-D + ((F*A)/C)) / (E - ((F*B)/C));
     out->x->coeff[1][0] = -(F/C) / (E - ((F*B)/C));
     out->x->coeff[0][1] =  1.0 / (E - ((F*B)/C));
-    out->y->coeff[0][0] = -D + ((E*A)/B) / (F - ((C*E)/B));
+    out->y->coeff[0][0] = (-D + ((E*A)/B)) / (F - ((C*E)/B));
     out->y->coeff[1][0] = -(E/B) / (F - ((C*E)/B));
     out->y->coeff[0][1] =  1.0 / (F - ((C*E)/B));
@@ -1140 +1162 @@
     //
     if (p_psIsProjectionLinear((psPlaneTransform *) in)) {
-        printf("COOL: is linear\n");
         return(p_psPlaneTransformLinearInvert((psPlaneTransform *) in));
     }

trunk/psLib/src/astronomy/psCoord.c

@@ -10 +10 @@
 *  @author GLG, MHPCC
 *
-*  @version $Revision: 1.63 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2005-04-19 01:22:59 $
+*  @version $Revision: 1.64 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2005-04-19 04:16:02 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -84 +84 @@
  X1 = (-D + ((F*A)/C)) / (E - ((F*B)/C)) +
       (X2 * -((F/C) / (E - ((F*B)/C)))) +
-             (Y2 * (1.0 / (E - ((F*B)/C))));
-        Y1 = (-D + ((E*A)/B))/(F - ((C*E)/B)) +
-             (X2 * -((E/B) / (F - ((C*E)/B)))) +
-             (Y2 * (1.0 / (F - ((C*E)/B))));
+      (Y2 * (1.0 / (E - ((F*B)/C))));
+ Y1 = (-D + ((E*A)/B))/(F - ((C*E)/B)) +
+      (X2 * -((E/B) / (F - ((C*E)/B)))) +
+      (Y2 * (1.0 / (F - ((C*E)/B))));
  
 XXX: Since thre is now a general psPlaneTransformInvert() function, we
 should rename this.
+ 
+X1, Y1 = (2, 3)
+ 
+    X2 = 1 + 4 + 9   = 14
+    Y2 = 4 + 10 + 18 = 32
+ 
+in inverse
+ 
+Y1 = (32 - ((5/2) * 14) - 4 + ((5)/2)) / (6 - ((15)/2));
+Y1 = (32 - 35 - 4 + (5/2)) / (6 - 15/2);
+Y1 = (-7 + (5/2)) / (12/2 - 15/2);
+Y1 = (-7 + (5/2)) / (-3/2);
+Y1 = (-14/2 + 5/2) / (-3/2);
+Y1 = (-9/2) / (-3/2);
+Y1 = (9/2) / (3/2);
+Y1 = 3
+ 
+X1 = (Y2 - ((F/C) * X2) - D + ((F*A)/C)) / (E - ((F*B)/C));
+X1 = (32 - ((6/3) * 14) - 4 + ((6)/3)) / (5 - ((12)/3));
+X1 = (32 - (2 * 14) - 4 + 2 / (5 - 4);
+X1 = (32 - 28 - 4 + 2) / 1;
+X1 = 2;
  
  *****************************************************************************/
@@ -125 +147 @@
     psPlaneTransform *out = psPlaneTransformAlloc(2, 2);
 
-    out->x->coeff[0][0] = -D + ((F*A)/C) / (E - ((F*B)/C));
+    out->x->coeff[0][0] = (-D + ((F*A)/C)) / (E - ((F*B)/C));
     out->x->coeff[1][0] = -(F/C) / (E - ((F*B)/C));
     out->x->coeff[0][1] =  1.0 / (E - ((F*B)/C));
-    out->y->coeff[0][0] = -D + ((E*A)/B) / (F - ((C*E)/B));
+    out->y->coeff[0][0] = (-D + ((E*A)/B)) / (F - ((C*E)/B));
     out->y->coeff[1][0] = -(E/B) / (F - ((C*E)/B));
     out->y->coeff[0][1] =  1.0 / (F - ((C*E)/B));
@@ -1140 +1162 @@
     //
     if (p_psIsProjectionLinear((psPlaneTransform *) in)) {
-        printf("COOL: is linear\n");
         return(p_psPlaneTransformLinearInvert((psPlaneTransform *) in));
     }

trunk/psLib/src/dataManip/psFunctions.c

@@ -7 +7 @@
  *  polynomials.  It also contains a Gaussian functions.
  *
- *  @version $Revision: 1.97 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-04-19 00:34:04 $
+ *  @version $Revision: 1.98 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-04-19 04:16:02 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
  *
- *  XXX: What happens if the polyEVal functions are called with data of the wrong
+ *  XXX: What happens if the polyEval functions are called with data of the wrong
  *       type?
  *  XXX: Should the "coeffErr[]" be used as well?  Bug ???.  Ignore coeffErr
  *
+ *  XXX: In the various polyAlloc(n) functions, n is really the order of the
+ *  polynomial plus 1.  To create a 2nd-order polynomial, n == 3.
  */
 /*****************************************************************************/
@@ -175 +177 @@
 static void dPolynomial2DFree(psDPolynomial2D* myPoly)
 {
-    psS32 x = 0;
-
-    for (x = 0; x < myPoly->nX; x++) {
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nX is %d\n", myPoly->nX);
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nY is %d\n", myPoly->nY);
+    for (psS32 x = 0; x < myPoly->nX; x++) {
         psFree(myPoly->coeff[x]);
         psFree(myPoly->coeffErr[x]);

trunk/psLib/src/image/psImageStats.c

@@ -9 +9 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.71 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-04-08 17:58:57 $
+ *  @version $Revision: 1.72 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-04-19 04:16:02 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -40 +40 @@
  
 XXX: use static psVectors
+ 
+XXX: optimize this.  2k vs 4k, sample mean, takes8 seconds on Gene's machine.
+Should take .2.
  *****************************************************************************/
 psStats* psImageStats(psStats* stats,

trunk/psLib/src/imageops/psImageStats.c

@@ -9 +9 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.71 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-04-08 17:58:57 $
+ *  @version $Revision: 1.72 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-04-19 04:16:02 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -40 +40 @@
  
 XXX: use static psVectors
+ 
+XXX: optimize this.  2k vs 4k, sample mean, takes8 seconds on Gene's machine.
+Should take .2.
  *****************************************************************************/
 psStats* psImageStats(psStats* stats,

trunk/psLib/src/math/psPolynomial.c

@@ -7 +7 @@
  *  polynomials.  It also contains a Gaussian functions.
  *
- *  @version $Revision: 1.97 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-04-19 00:34:04 $
+ *  @version $Revision: 1.98 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-04-19 04:16:02 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
  *
- *  XXX: What happens if the polyEVal functions are called with data of the wrong
+ *  XXX: What happens if the polyEval functions are called with data of the wrong
  *       type?
  *  XXX: Should the "coeffErr[]" be used as well?  Bug ???.  Ignore coeffErr
  *
+ *  XXX: In the various polyAlloc(n) functions, n is really the order of the
+ *  polynomial plus 1.  To create a 2nd-order polynomial, n == 3.
  */
 /*****************************************************************************/
@@ -175 +177 @@
 static void dPolynomial2DFree(psDPolynomial2D* myPoly)
 {
-    psS32 x = 0;
-
-    for (x = 0; x < myPoly->nX; x++) {
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nX is %d\n", myPoly->nX);
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nY is %d\n", myPoly->nY);
+    for (psS32 x = 0; x < myPoly->nX; x++) {
         psFree(myPoly->coeff[x]);
         psFree(myPoly->coeffErr[x]);

trunk/psLib/src/math/psSpline.c

@@ -7 +7 @@
  *  polynomials.  It also contains a Gaussian functions.
  *
- *  @version $Revision: 1.97 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-04-19 00:34:04 $
+ *  @version $Revision: 1.98 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-04-19 04:16:02 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
  *
- *  XXX: What happens if the polyEVal functions are called with data of the wrong
+ *  XXX: What happens if the polyEval functions are called with data of the wrong
  *       type?
  *  XXX: Should the "coeffErr[]" be used as well?  Bug ???.  Ignore coeffErr
  *
+ *  XXX: In the various polyAlloc(n) functions, n is really the order of the
+ *  polynomial plus 1.  To create a 2nd-order polynomial, n == 3.
  */
 /*****************************************************************************/
@@ -175 +177 @@
 static void dPolynomial2DFree(psDPolynomial2D* myPoly)
 {
-    psS32 x = 0;
-
-    for (x = 0; x < myPoly->nX; x++) {
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nX is %d\n", myPoly->nX);
+    //printf("dPolynomial2DFree(): HMMM: myPoly->nY is %d\n", myPoly->nY);
+    for (psS32 x = 0; x < myPoly->nX; x++) {
         psFree(myPoly->coeff[x]);
         psFree(myPoly->coeffErr[x]);

trunk/psLib/test/astronomy/tst_psAstrometry01.c

@@ -5 +5 @@
 *  @author GLG, MHPCC
 *
-*  @version $Revision: 1.21 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2005-04-19 01:22:23 $
+*  @version $Revision: 1.22 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2005-04-19 04:16:02 $
 *
 * XXX: Must test
@@ -828 +828 @@
 
 /******************************************************************************
-XXX: This is only a rudimentary test of the psPlaneTransform()
+XXX: This is only a rudimentary test of the psPlaneTransformInvert()
 function.  It tests a few NULL input parameter conditions.
  
@@ -844 +844 @@
     //
     trans->x->coeff[0][0] = 1.0;
-    trans->x->coeff[0][1] = 1.0;
-    trans->x->coeff[1][0] = 1.0;
+    trans->x->coeff[0][1] = 3.0;
+    trans->x->coeff[1][0] = 2.0;
     trans->x->coeff[1][1] = 1.0;
-    trans->y->coeff[0][0] = 1.0;
-    trans->y->coeff[0][1] = 1.0;
-    trans->y->coeff[1][0] = 1.0;
+    trans->y->coeff[0][0] = 4.0;
+    trans->y->coeff[0][1] = 6.0;
+    trans->y->coeff[1][0] = 5.0;
     trans->y->coeff[1][1] = 1.0;
 
@@ -878 +878 @@
     }
 
-    printf("----------------------------------------------------------------------------------\n");
-    printf("Calling psPlaneTransformInvert with acceptable data.\n");
-    transOut = psPlaneTransformInvert(NULL, trans, myRegion, 10);
-    if (transOut == NULL) {
-        printf("TEST ERROR: psPlaneTransformInvert() returned a NULL psPlaneTransform.\n");
-        testStatus = false;
+    printf("Calling psPlaneTransformInvert with an acceptable linear transformation.\n");
+    for (psS32 n = 0 ; n < 100 ; n++) {
+        //
+        // We create a random linear transformation and hope it is invertible.
+        //
+        trans->x->coeff[0][0] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->x->coeff[0][1] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->x->coeff[1][0] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->x->coeff[1][1] = 0.1;
+        trans->y->coeff[0][0] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->y->coeff[0][1] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->y->coeff[1][0] = ((psF32) (random() % 1000000)) / 100000.0;
+        trans->y->coeff[1][1] = 0.1;
+
+        transOut = psPlaneTransformInvert(NULL, trans, myRegion, 10);
+        if (transOut == NULL) {
+            printf("TEST ERROR: psPlaneTransformInvert() returned a NULL psPlaneTransform.\n");
+            testStatus = false;
+        }
+        psPlane *in = psPlaneAlloc();
+        in->x = 2.0;
+        in->y = 3.0;
+
+        psPlane *mid = psPlaneTransformApply(NULL, trans, in);
+        psPlane *out = psPlaneTransformApply(NULL, transOut, mid);
+
+        if (((fabs(in->x - out->x) > FLT_EPSILON)) ||
+                ((fabs(in->y - out->y) > FLT_EPSILON))) {
+            printf("TEST ERROR: in coords were (%f, %f), out coords were (%f, %f).\n",
+                   in->x, in->y, out->x, out->y);
+            testStatus = false;
+        }
+
+        psFree(in);
+        psFree(mid);
+        psFree(out);
+        psFree(transOut);
     }
 
     psFree(trans);
-    psFree(transOut);
     psFree(myRegion);
 
     return(!testStatus);
 }
+//This code will

trunk/psLib/test/image/Makefile.am

@@ -20 +20 @@
         tst_psImageManip \
         tst_psImageStats \
+        tst_psImageStatsSpeed \
         tst_psImageExtraction \
         tst_psImageConvolve \
@@ -43 +44 @@
 tst_psImageManip_SOURCES =  tst_psImageManip.c
 tst_psImageStats_SOURCES =  tst_psImageStats.c
+tst_psImageStatsSpeed_SOURCES =  tst_psImageStatsSpeed.c
 tst_psImageExtraction_SOURCES =  tst_psImageExtraction.c
 tst_psImageConvolve_SOURCES =  tst_psImageConvolve.c