Changeset 4160 for trunk/psLib/src/math/psMatrix.c

Timestamp:

Jun 8, 2005, 12:28:08 PM (21 years ago)

Author:

desonia

Message:

Bug #403 -- implemented and tested non-square matrix multiply.

File:

• trunk/psLib/src/math/psMatrix.c (modified) (6 diffs)

trunk/psLib/src/math/psMatrix.c

@@ -21 +21 @@
  *  @author Robert DeSonia, MHPCC
  *
- *  @version $Revision: 1.30 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-06-08 02:38:02 $
+ *  @version $Revision: 1.31 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-06-08 22:28:07 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -33 +33 @@
 #include <gsl/gsl_matrix.h>
 #include <gsl/gsl_linalg.h>
+#include <gsl/gsl_blas.h>
 #include <gsl/gsl_permutation.h>
 #include <gsl/gsl_eigen.h>
@@ -138 +139 @@
     numRows = inImage->numRows;
     numCols = inImage->numCols;
-    for(i=0; i<numRows; i++) {
-        for(j=0; j<numCols; j++) {
-            if(inImage->type.type == PS_TYPE_F32) {
-                outGslMatrix->data[i*numCols+j] = (psF64)inImage->data.F32[i][j];
-            } else {
+    if(inImage->type.type == PS_TYPE_F32) {
+        for(i=0; i<numRows; i++) {
+            for(j=0; j<numCols; j++) {
+                outGslMatrix->data[i*numCols+j] = inImage->data.F32[i][j];
+            }
+        }
+    } else {
+        for(i=0; i<numRows; i++) {
+            for(j=0; j<numCols; j++) {
                 outGslMatrix->data[i*numCols+j] = inImage->data.F64[i][j];
             }
@@ -160 +165 @@
     numRows = outImage->numRows;
     numCols = outImage->numCols;
-    for(i=0; i<numRows; i++) {
-        for(j=0; j<numCols; j++) {
-            if(outImage->type.type == PS_TYPE_F32) {
-                outImage->data.F32[i][j] = (psF32)inGslMatrix->data[i*numCols+j];
-            } else {
+    if(outImage->type.type == PS_TYPE_F32) {
+        for(i=0; i<numRows; i++) {
+            for(j=0; j<numCols; j++) {
+                outImage->data.F32[i][j] = inGslMatrix->data[i*numCols+j];
+            }
+        }
+    } else {
+        for(i=0; i<numRows; i++) {
+            for(j=0; j<numCols; j++) {
                 outImage->data.F64[i][j] = inGslMatrix->data[i*numCols+j];
             }
@@ -401 +410 @@
     }
 
-    outImage = psImageRecycle(outImage, inImage1->numCols, inImage2->numRows, inImage2->type.type);
+    outImage = psImageRecycle(outImage, inImage2->numCols, inImage1->numRows, inImage2->type.type);
 
     // Initialize GSL data
@@ -411 +420 @@
 
     // Perform multiplication
+    // gsl_blas_dgemm(CblasNoTrans, CblasNoTrans, 1.0, m1, m2, 0.0, m3);
     gsl_linalg_matmult(m1, m2, m3);
-    // N.B., gsl_blas_sgemm/dgemm could have been used instead, e.g.,
-    // gsl_blas_dgemm(CblasNoTrans, CblasNoTrans, 1.0, m1, m2, m3)
 
     // Copy GSL matrix data to psImage data