loadMatrix, saveMatrix

4 files changed, 117 insertions, 7 deletions

diff --git a/packages/base/src/C/vector-aux.c b/packages/base/src/C/vector-aux.c
index 7cdc750..f1bb371 100644
--- a/packages/base/src/C/vector-aux.c
+++ b/packages/base/src/C/vector-aux.c
@@ -12,6 +12,7 @@ typedef float  complex TCF;
 #include <string.h>
 #include <math.h>
 #include <stdio.h>
+#include <stdlib.h>
 
 #define MACRO(B) do {B} while (0)
 #define ERROR(CODE) MACRO(return CODE;)
@@ -648,4 +649,49 @@ int zipQ(int code, KQVEC(a), KQVEC(b), QVEC(r)) {
     }
 }
 
+////////////////////////////////////////////////////////////////////////////////
+
+int vectorScan(char * file, int* n, double**pp){
+    FILE * fp;
+    fp = fopen (file, "r");
+    int nbuf = 100*100;
+    double * p = (double*)malloc(nbuf*sizeof(double));
+    int k=0;
+    double d;
+    int ok;
+    for (;;) {
+        ok = fscanf(fp,"%lf",&d);
+        if (ok<1) {
+            break;
+        }
+        if (k==nbuf) {
+            nbuf = nbuf * 2;
+            p = (double*)realloc(p,nbuf*sizeof(double));
+            //printf("R\n");
+        }
+        p[k++] = d;
+    }
+    *n = k;
+    *pp = p;
+    fclose(fp);
+    OK
+}    
+
+int saveMatrix(char * file, char * format, KDMAT(a)){
+    FILE * fp;
+    fp = fopen (file, "w");
+    int r, c;
+    for (r=0;r<ar; r++) {
+        for (c=0; c<ac; c++) {
+            fprintf(fp,format,ap[r*ac+c]);
+            if (c<ac-1) {
+                fprintf(fp," ");
+            } else {
+                fprintf(fp,"\n");
+            }
+        }
+    }
+    fclose(fp);
+    OK
+}
 
diff --git a/packages/base/src/Data/Packed/IO.hs b/packages/base/src/Data/Packed/IO.hs
index dbb2943..db03d5f 100644
--- a/packages/base/src/Data/Packed/IO.hs
+++ b/packages/base/src/Data/Packed/IO.hs
@@ -14,7 +14,7 @@
 
 module Data.Packed.IO (
     dispf, disps, dispcf, vecdisp, latexFormat, format,
-    readMatrix, fromArray2D
+    readMatrix, fromArray2D, loadMatrix, saveMatrix
 ) where
 
 import Data.Packed
@@ -22,6 +22,8 @@ import Data.Packed.Development
 import Text.Printf(printf)
 import Data.List(intersperse)
 import Data.Complex
+import Numeric.Vectorized(vectorScan,saveMatrix)
+import Control.Applicative((<$>))
 
 {- | Creates a string from a matrix given a separator and a function to show each entry. Using
 this function the user can easily define any desired display function:
@@ -139,3 +141,17 @@ dispcf d m = sdims m ++ "\n" ++ format "  " (showComplex d) m
 readMatrix :: String -> Matrix Double
 readMatrix = fromLists . map (map read). map words . filter (not.null) . lines
 
+--------------------------------------------------------------------------------
+
+apparentCols :: FilePath -> IO Int
+apparentCols s = f . dropWhile null . map words . lines <$> readFile s
+  where
+    f [] = 0
+    f (x:_) = length x
+
+loadMatrix :: FilePath -> IO (Matrix Double)
+loadMatrix f = do
+    v <- vectorScan f
+    c <- apparentCols f
+    return (reshape c v)
+
diff --git a/packages/base/src/Numeric/LinearAlgebra/Data.hs b/packages/base/src/Numeric/LinearAlgebra/Data.hs
index 2754576..45fc00c 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Data.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Data.hs
@@ -43,9 +43,10 @@ module Numeric.LinearAlgebra.Data(
     find, maxIndex, minIndex, maxElement, minElement, atIndex,
 
     -- * IO
-    disp, dispf, disps, dispcf, latexFormat, format, readMatrix,
-
-    -- | loadMatrix, saveMatrix, fromFile, fileDimensions, fscanfVector, fprintfVector, freadVector, fwriteVector
+    disp,
+    loadMatrix, saveMatrix,
+    latexFormat,
+    dispf, disps, dispcf, format,
 
 -- * Conversion
     Convert(..),
diff --git a/packages/base/src/Numeric/Vectorized.hs b/packages/base/src/Numeric/Vectorized.hs
index 3814579..a2d7f70 100644
--- a/packages/base/src/Numeric/Vectorized.hs
+++ b/packages/base/src/Numeric/Vectorized.hs
@@ -16,20 +16,29 @@ module Numeric.Vectorized (
     FunCodeS(..), toScalarR, toScalarF, toScalarC, toScalarQ,
     FunCodeV(..), vectorMapR, vectorMapC, vectorMapF, vectorMapQ,
     FunCodeSV(..), vectorMapValR, vectorMapValC, vectorMapValF, vectorMapValQ,
-    FunCodeVV(..), vectorZipR, vectorZipC, vectorZipF, vectorZipQ
+    FunCodeVV(..), vectorZipR, vectorZipC, vectorZipF, vectorZipQ,
+    vectorScan, saveMatrix
 ) where
 
 import Data.Packed.Internal.Common
 import Data.Packed.Internal.Signatures
 import Data.Packed.Internal.Vector
+import Data.Packed.Internal.Matrix
 
 import Data.Complex
-import Foreign.Marshal.Alloc(free)
-import Foreign.Marshal.Array(newArray)
+import Foreign.Marshal.Alloc(free,malloc)
+import Foreign.Marshal.Array(newArray,copyArray)
 import Foreign.Ptr(Ptr)
+import Foreign.Storable(peek)
 import Foreign.C.Types
+import Foreign.C.String
 import System.IO.Unsafe(unsafePerformIO)
 
+import Control.Monad(when)
+import Control.Applicative((<$>))
+
+
+
 fromei x = fromIntegral (fromEnum x) :: CInt
 
 data FunCodeV = Sin
@@ -271,3 +280,41 @@ vectorZipQ = vectorZipAux c_vectorZipQ
 
 foreign import ccall unsafe "zipQ" c_vectorZipQ :: CInt -> TQVQVQV
 
+--------------------------------------------------------------------------------
+
+foreign import ccall unsafe "vectorScan" c_vectorScan
+    :: CString -> Ptr CInt -> Ptr (Ptr Double) -> IO CInt
+
+vectorScan :: FilePath -> IO (Vector Double)
+vectorScan s = do
+    pp <- malloc
+    pn <- malloc
+    cs <- newCString s
+    ok <- c_vectorScan cs pn pp
+    when (not (ok == 0)) $
+        error ("vectorScan \"" ++ s ++"\"")
+    n <- fromIntegral <$> peek pn
+    p <- peek pp
+    v <- createVector n
+    free pn
+    free cs
+    unsafeWith v $ \pv -> copyArray pv p n
+    free p
+    free pp
+    return v
+
+--------------------------------------------------------------------------------
+
+foreign import ccall unsafe "saveMatrix" c_saveMatrix
+    :: CString -> CString -> TM
+
+saveMatrix :: FilePath -> String -> Matrix Double -> IO ()
+saveMatrix name format m = do
+    cname   <- newCString name
+    cformat <- newCString format
+    app1 (c_saveMatrix cname cformat) mat m "saveMatrix"
+    free cname
+    free cformat
+    return ()
+
+