Changeset 26418

Timestamp:

Dec 15, 2009, 6:15:25 PM (16 years ago)

Author:

eugene

Message:

add some matrix math operations and fix svd

Location:

branches/eam_branches/20091201/Ohana/src/opihi

Files:

• cmd.data/Makefile (modified) (1 diff)
• cmd.data/init.c (modified) (2 diffs)
• cmd.data/matrix.c (added)
• cmd.data/svd.c (modified) (1 diff)
• cmd.data/test/matrix.sh (added)
• cmd.data/test/svd.sh (added)
• doc/svdcmp_eispack.f (added)
• include/data.h (modified) (1 diff)

branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/Makefile

@@ -69 +69 @@
 $(SRC)/load.$(ARCH).o           \
 $(SRC)/lookup.$(ARCH).o \
+$(SRC)/matrix.$(ARCH).o \
 $(SRC)/mkrgb.$(ARCH).o  \
 $(SRC)/mcreate.$(ARCH).o        \

branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/init.c

@@ -58 +58 @@
 int load             PROTO((int, char **));
 int lookup           PROTO((int, char **));
+int matrix           PROTO((int, char **));
 int mkrgb            PROTO((int, char **));
 int mcreate          PROTO((int, char **));
@@ -195 +196 @@
   {1, "minterp",      minterp,          "interpolate image pixels"},
   {1, "iminterp",     minterp,          "interpolate image pixels"},
+  {1, "matrix",       matrix,           "matrix math operations"},
   {1, "mkrgb",        mkrgb,            "convert 3 images to rgb jpeg (use Kapa for better control)"},
   {1, "mset",         mset,             "insert a vector in an image"},

branches/eam_branches/20091201/Ohana/src/opihi/cmd.data/svd.c

@@ -57 +57 @@
   /* use a bcopy instead? */
 
+  // try C.R. Bond's version -- requires matrices in the form A[row][col] not A[row*Ncol + col]
+  if (1) { 
+      int j;
+      double **a, **u, **v, *q;
+      ALLOCATE (a, double *, Ny);
+      ALLOCATE (u, double *, Ny);
+      ALLOCATE (v, double *, Ny);
+      for (i = 0; i < Ny; i++) {
+          ALLOCATE (a[i], double, Nx);
+          ALLOCATE (u[i], double, Nx);
+          ALLOCATE (v[i], double, Nx);
+      }   
+      ALLOCATE (q, double, Nx);
+
+      for (j = 0; j < Ny; j++) {
+          for (i = 0; i < Nx; i++) {
+              a[j][i] = A[j*Nx + i];
+              u[j][i] = 0;
+              v[j][i] = 0;
+          }
+      }
+      for (i = 0; i < Nx; i++) {
+          q[i] = 0;
+      }
+
+      status = svdcmp_bond_raw (Ny, Nx, 1, 1, FLT_EPSILON, 1e-6, a, q, u, v);
+      fprintf (stderr, "status: %d\n", status);
+
+      // copy u q v back to U W V:
+      for (j = 0; j < Ny; j++) {
+          for (i = 0; i < Nx; i++) {
+              U[j*Nx + i] = u[j][i];
+              V[j*Nx + i] = v[j][i];
+          }
+      }
+      for (i = 0; i < Nx; i++) {
+          W[i] = q[i];
+      }
+
+      for (j = 0; j < Ny; j++) {
+          free(a[j]);
+          free(u[j]);
+          free(v[j]);
+      }
+      free (a);
+      free (u);
+      free (v);
+      free (q);
+
+      return TRUE;
+  }
+
   status = svdcmp (U, W, V, Nx, Ny);
   if (!status) {

branches/eam_branches/20091201/Ohana/src/opihi/include/data.h

@@ -89 +89 @@
 int svdcmp (float *a, opihi_flt *w, float *v, int Nx, int Ny);
 
+/* in svdcmp_bond_raw.c */
+int svdcmp_bond_raw(int m, int n, int withu, int withv, double eps, double tol, double **a, double *q, double **u, double **v);
+
 /* mrqmin.c */
 opihi_flt mrqcof (opihi_flt *x, opihi_flt *y, opihi_flt *dy, int Npts,