; "matfile.scm", Read MAT-File Format version 4 (MATLAB)
; Copyright (C) 2001, 2002, 2003 Aubrey Jaffer
;
;Permission to copy this software, to modify it, to redistribute it,
;to distribute modified versions, and to use it for any purpose is
;granted, subject to the following restrictions and understandings.
;
;1. Any copy made of this software must include this copyright notice
;in full.
;
;2. I have made no warranty or representation that the operation of
;this software will be error-free, and I am under no obligation to
;provide any services, by way of maintenance, update, or otherwise.
;
;3. In conjunction with products arising from the use of this
;material, there shall be no use of my name in any advertising,
;promotional, or sales literature without prior written consent in
;each case.
(require 'array)
(require 'byte)
(require 'byte-number)
(require-if 'compiling 'string-case) ; string-ci->symbol used by matfile:load
(define (unwritten-stubber name)
(lambda (arg) (slib:error 'name 'not 'written "matfile.scm")))
(define bytes->vax-d-double (unwritten-stubber 'bytes->vax-d-double))
(define bytes->vax-g-double (unwritten-stubber 'bytes->vax-g-double))
(define bytes->cray-double (unwritten-stubber 'bytes->cray-double))
(define bytes->vax-d-float (unwritten-stubber 'bytes->vax-d-float))
(define bytes->vax-g-float (unwritten-stubber 'bytes->vax-g-float))
(define bytes->cray-float (unwritten-stubber 'bytes->cray-float))
;;@code{(require 'matfile)}
;;@ftindex matfile
;;@ftindex matlab
;;
;;@uref{http://www.mathworks.com/access/helpdesk/help/pdf_doc/matlab/matfile_format.pdf}
;;
;;@noindent
;;This package reads MAT-File Format version 4 (MATLAB) binary data
;;files. MAT-files written from big-endian or little-endian computers
;;having IEEE format numbers are currently supported. Support for files
;;written from VAX or Cray machines could also be added.
;;
;;@noindent
;;The numeric and text matrix types handled; support for @dfn{sparse}
;;matrices awaits a sample file.
(define (bytes->long lst)
(bytes->integer lst -4))
(define (bytes->short lst)
(bytes->integer lst -2))
(define (bytes->ushort lst)
(bytes->integer lst 2))
;;Version 4 MAT-file endianness cannot be detected solely from the
;;first word; it is ambiguous when 0.
(define (matfile:read-matrix port)
(define null (integer->char 0))
(define (read1 endian type mrows ncols imagf namlen)
(set! type (bytes->long type))
(set! mrows (bytes->long mrows))
(set! ncols (bytes->long ncols))
(set! imagf (bytes->long imagf))
(set! namlen (+ -1 (bytes->long namlen)))
(let ((d-prot (modulo (quotient type 10) 10))
(d-endn (case (quotient type 1000)
((0 ;ieee-little-endian
2 ;vax-d-float
3) endian) ;vag-g-float
((1 ;ieee-big-endian
4) (- endian)) ;cray
(else #f)))
(m-type (case (modulo type 10)
((0) 'numeric)
((1) 'text)
((2) 'sparse)
(else #f))))
(define d-leng (case d-prot
((0) 8)
((1 2) 4)
((3 4) 2)
((5) 1)
(else #f)))
(define d-conv (case d-prot
((0) (case (quotient type 1000)
((0 1) bytes->ieee-double)
((2) bytes->vax-d-double)
((3) bytes->vax-g-double)
((4) bytes->cray-double)))
((1) (case (quotient type 1000)
((0 1) bytes->ieee-float)
((2) bytes->vax-d-float)
((3) bytes->vax-g-float)
((4) bytes->cray-float)))
((2) bytes->long)
((3) bytes->short)
((4) bytes->ushort)
((5) (if (eqv? 'text m-type)
(lambda (lst) (integer->char (byte-ref lst 0)))
(lambda (lst) (byte-ref lst 0))))
(else #f)))
;;(@print d-leng d-endn m-type type mrows ncols imagf namlen d-conv)
(cond ((and (= 0 (modulo (quotient type 100) 10) (quotient type 65536))
d-leng d-endn m-type
(<= 0 imagf 1)
(< 0 mrows #xFFFFFF)
(< 0 ncols #xFFFFFF)
(< 0 namlen #xFFFF))
(set! imagf (case imagf ((0) #f) ((1) #t)))
(let ((namstr (make-string namlen))
(mat (case m-type
((numeric) (make-array
(case d-prot
((0) ((if imagf A:floC64b A:floR64b)))
((1) ((if imagf A:floC32b A:floR32b)))
((2) (A:fixZ32b))
((3) (A:fixZ16b))
((4) (A:fixN16b))
((5) (A:fixN8b))
(else (slib:error 'p 'type d-prot)))
mrows ncols))
((text) (make-array "." mrows ncols))
((sparse) (slib:error 'sparse '?))))
(d-endn*leng (* -1 d-endn d-leng)))
(do ((idx 0 (+ 1 idx)))
((>= idx namlen))
(string-set! namstr idx (read-char port)))
;;(@print namstr)
(if (not (eqv? null (read-char port)))
(slib:error 'matfile 'string 'missing null))
(do ((jdx 0 (+ 1 jdx)))
((>= jdx ncols))
(do ((idx 0 (+ 1 idx)))
((>= idx mrows))
(array-set! mat (d-conv (read-bytes d-endn*leng port))
idx jdx)))
(if imagf
(do ((jdx 0 (+ 1 jdx)))
((>= jdx ncols))
(do ((idx 0 (+ 1 idx)))
((>= idx mrows))
(array-set! mat
(+ (* (d-conv (read-bytes d-endn*leng port))
+i)
(array-ref mat idx jdx))
idx jdx))))
(list namstr mat)))
(else #f))))
;;(trace read1)
(let* ((type (read-bytes 4 port))
(mrows (read-bytes 4 port))
(ncols (read-bytes 4 port))
(imagf (read-bytes 4 port))
(namlen (read-bytes 4 port)))
;;Try it with either endianness:
(or (read1 1 type mrows ncols imagf namlen)
(read1 -1
(bytes-reverse type)
(bytes-reverse mrows)
(bytes-reverse ncols)
(bytes-reverse imagf)
(bytes-reverse namlen)))))
;;@body @1 should be a string naming an existing file containing a
;;MATLAB Version 4 MAT-File. The @0 procedure reads matrices from the
;;file and returns a list of the results; a list of the name string and
;;array for each matrix.
(define (matfile:read filename)
(call-with-open-ports
(open-file filename 'rb)
(lambda (port)
(do ((mat (matfile:read-matrix port) (matfile:read-matrix port))
(mats '() (cons mat mats)))
((or (not mat) (eof-object? (peek-char port)))
(if (and (null? mats) (not mat))
'()
(reverse (cons mat mats))))))))
;;@body @1 should be a string naming an existing file containing a
;;MATLAB Version 4 MAT-File. The @0 procedure reads matrices from the
;;file and defines the @code{string-ci->symbol} for each matrix to its
;;corresponding array. @0 returns a list of the symbols defined.
(define (matfile:load filename)
(require 'string-case)
(let ((mats (matfile:read filename)))
(for-each (lambda (nam-mat)
(and nam-mat
(slib:eval
(list 'define
(string-ci->symbol (car nam-mat))
(list 'quote (cadr nam-mat))))))
mats)
(map string-ci->symbol (map car mats))))
;;(trace-all "/home/jaffer/slib/matfile.scm")