Home My Page Projects Code Snippets Project Openings 3D graphics for Standard ML
Summary Activity SCM

SCM Repository

[sml3d] View of /src/common/matrix3f.sml
ViewVC logotype

View of /src/common/matrix3f.sml

Parent Directory Parent Directory | Revision Log Revision Log


Revision 210 - (download) (annotate)
Tue Jun 17 03:29:03 2008 UTC (11 years, 4 months ago) by jhr
File size: 6565 byte(s)
  Added rotations.
(* matrix3f.sml
 *
 * COPYRIGHT (c) 2008 John Reppy (http://cs.uchicago.edu/~jhr)
 * All rights reserved.
 *
 * Single-precision 3x3 matrices
 *)

structure Matrix3f : sig

    type mat3f

    type float = SML3dTypes.float
    type vec3f = SML3dTypes.vec3f

    val mat : {
	    m11 : float, m12 : float, m13 : float,
	    m21 : float, m22 : float, m23 : float,
	    m31 : float, m32 : float, m33 : float
	  } -> mat3f

  (* standard matrices *)
    val identity : mat3f			(* identity matrix *)
    val isoscale : float -> mat3f		(* isotropic scaling *)
    val scale    : vec3f -> mat3f		(* anisotropic scaling *)
    val rotateX  : float -> mat3f		(* rotation around X axis (in degrees) *)
    val rotateY  : float -> mat3f		(* rotation around Y axis (in degrees) *)
    val rotateZ  : float -> mat3f		(* rotation around Z axis (in degrees) *)
    val rotate   : float * vec3f -> mat3f	(* rotation around an axis (in degrees) *)

    val mxv : mat3f * vec3f -> vec3f	(* matrix times column vector *)
    val vxm : vec3f * mat3f -> vec3f	(* row vector times matrix *)
    val mxm : mat3f * mat3f -> mat3f	(* matrix times matrix *)

    val transpose : mat3f -> mat3f

    val det : mat3f -> float

    val inverse : mat3f -> mat3f option

  (* returns matrix as flat vector in column-major order *)
    val toVector : mat3f -> float Vector.vector

  end = struct

    type float = SML3dTypes.float
    type double = SML3dTypes.double
    type vec3f = SML3dTypes.vec3f

    val epsilonf : float = 1.0e~5
    val epsilond : double = 1.0e~7

  (* represented in column-major order, since that is what OpenGL LoadMatrix
   * expects.
   *)
    datatype mat3f = M of float Vector.vector

  (* index = 3*(c-1) + r - 1 = 3*c + r - 4 *)
    fun ! (a : float Vector.vector, (r, c)) = Unsafe.Vector.sub(a, 3*c + r - 4)
    infix !

    fun new {m11, m12, m13, m21, m22, m23, m31, m32, m33} =
	  M(Vector.fromList [m11, m21, m31, m12, m22, m32,  m13, m23, m33])

    val mat = new

    val identity = new {
	    m11 = 1.0, m12 = 0.0, m13 = 0.0,
	    m21 = 0.0, m22 = 1.0, m23 = 0.0,
	    m31 = 0.0, m32 = 0.0, m33 = 1.0
	  }

    fun isoscale s = new {
	    m11 = s,   m12 = 0.0, m13 = 0.0,
	    m21 = 0.0, m22 = s,   m23 = 0.0,
	    m31 = 0.0, m32 = 0.0, m33 = s
	  }

    fun scale {x, y, z} = new {
	    m11 = x,   m12 = 0.0, m13 = 0.0,
	    m21 = 0.0, m22 = y,   m23 = 0.0,
	    m31 = 0.0, m32 = 0.0, m33 = z
	  }

  (* rotation around X axis (in degrees) *)
    fun rotateX angle = let
	  val angle = Float.toRadians angle
	  val s = Float.sin angle
	  val c = Float.cos angle
	  in
	    new {
		m11 = 1.0, m12 = 0.0, m13 = 0.0,
		m21 = 0.0, m22 = c,   m23 = ~s,
		m31 = 0.0, m32 = s,   m33 = c
	      }
	  end

  (* rotation around Y axis (in degrees) *)
    fun rotateY angle = let
	  val angle = Float.toRadians angle
	  val s = Float.sin angle
	  val c = Float.cos angle
	  in
	    new {
		m11 = c,   m12 = 0.0, m13 = s,
		m21 = 0.0, m22 = 1.0, m23 = 0.0,
		m31 = ~s,  m32 = 0.0, m33 = c
	      }
	  end

  (* rotation around Z axis (in degrees) *)
    fun rotateZ angle = let
	  val angle = Float.toRadians angle
	  val s = Float.sin angle
	  val c = Float.cos angle
	  in
	    new {
		m11 = c,   m12 = ~s,  m13 = 0.0,
		m21 = s,   m22 = c,   m23 = 0.0,
		m31 = 0.0, m32 = 0.0, m33 = 1.0
	      }
	  end

  (* rotation around an axis (in degrees) *)
    fun rotate (angle, {x, y, z}) = let
	  val angle = Float.toRadians angle
	  val s = Float.sin angle
	  val c = Float.cos angle
	  val c' = 1.0 - c
	  in
	    new {
		m11 = x*x*c'+c,   m12 = x*y*c'-z*s, m13 = x*z*c'+y*s,
		m21 = y*x*c'+z*s, m22 = y*y*c'+c,   m23 = y*z*c'-x*s,
		m31 = x*z*c'-y*s, m32 = y*z*c'+x*s, m33 = z*z*c'+c
	      }
	  end

  (* multiply a matrix times a column vector *)
    fun mxv (M mat, {x, y, z}) = let
	  fun m arg = mat ! arg
	  in {
	    x = m(1,1) * x + m(1,2) * y + m(1,3) * z,
	    y = m(2,1) * x + m(2,2) * y + m(2,3) * z,
	    z = m(3,1) * x + m(3,2) * y + m(3,3) * z
	  } end
    
  (* multiply a row vector times a matrix *)
    fun vxm ({x, y, z}, M mat) = let
	  fun m arg = mat ! arg
	  in {
	    x = m(1,1) * x + m(2,1) * y + m(3,1) * z,
	    y = m(1,2) * x + m(2,2) * y + m(3,2) * z,
	    z = m(1,3) * x + m(2,3) * y + m(3,3) * z
	  } end

  (* matrix transpose *)
    fun transpose (M mat) = let
	  fun m arg = mat ! arg
	  in
	    new {
		m11 = m(1,1), m12 = m(2,1), m13 = m(3,1),
		m21 = m(1,2), m22 = m(2,2), m23 = m(3,2),
		m31 = m(1,3), m32 = m(2,3), m33 = m(3,3)
	      }
	  end

  (* matrix multiplication *)
    fun mxm (M aMat, M bMat) = let
	  fun a arg = aMat ! arg
	  fun b arg = bMat ! arg
	  in
	    new {
		m11 = a(1,1) * b(1,1) + a(1,2) * b(2,1) + a(1,3) * b(3,1),
		m12 = a(1,1) * b(1,2) + a(1,2) * b(2,2) + a(1,3) * b(3,2),
		m13 = a(1,1) * b(1,3) + a(1,2) * b(2,3) + a(1,3) * b(3,3),
		m21 = a(2,1) * b(1,1) + a(2,2) * b(2,1) + a(2,3) * b(3,1),
		m22 = a(2,1) * b(1,2) + a(2,2) * b(2,2) + a(2,3) * b(3,2),
		m23 = a(2,1) * b(1,3) + a(2,2) * b(2,3) + a(2,3) * b(3,3),
		m31 = a(3,1) * b(1,1) + a(3,2) * b(2,1) + a(3,3) * b(3,1),
		m32 = a(3,1) * b(1,2) + a(3,2) * b(2,2) + a(3,3) * b(3,2),
		m33 = a(3,1) * b(1,3) + a(3,2) * b(2,3) + a(3,3) * b(3,3)
	      }
	  end

    fun det (M mat) = let
	  fun m arg = mat ! arg
	  fun t (i,j,k) = m(1,i)*m(2,j)*m(3,k)
	  in
	    t(1,2,3) - t(1,3,2) - t(2,1,3) + t(2,3,1) + t(3,1,2) - t(3,2,1)
	  end

    fun inverse (M mat) = let
	  fun m arg = FP.ftod(mat ! arg)
	(* 2x2 determinant *)
	  fun det2 (a, b, c, d) = a*d - b*c
	(* 2x2 subdeterminants *)
	  val sd11 = det2 (m(2,2), m(2,3), m(3,2), m(3,3))
	  val sd12 = det2 (m(1,3), m(1,2), m(3,3), m(3,2))
	  val sd13 = det2 (m(1,2), m(1,3), m(2,2), m(2,3))
	  val sd21 = det2 (m(2,3), m(2,1), m(3,3), m(3,1))
	  val sd22 = det2 (m(1,1), m(1,3), m(3,1), m(3,3))
	  val sd23 = det2 (m(1,3), m(1,1), m(2,3), m(2,1))
	  val sd31 = det2 (m(2,1), m(2,2), m(3,1), m(3,2))
	  val sd32 = det2 (m(1,2), m(1,1), m(3,2), m(3,1))
	  val sd33 = det2 (m(1,1), m(1,2), m(2,1), m(2,2))
	(* determinant *)
	  val det = m(1,1)*sd11 - m(1,2)*sd12 + m(1,3)*sd13
	  in
	    if ((det > ~epsilond) andalso (det < epsilond))
	      then NONE
	      else let
		val detInv = 1.0 / det
		in
		  SOME(new {
		      m11 = FP.dtof(sd11*detInv),
		      m21 = FP.dtof(sd12*detInv),
		      m31 = FP.dtof(sd13*detInv),
		      m12 = FP.dtof(sd21*detInv),
		      m22 = FP.dtof(sd22*detInv),
		      m32 = FP.dtof(sd23*detInv),
		      m13 = FP.dtof(sd31*detInv),
		      m23 = FP.dtof(sd32*detInv),
		      m33 = FP.dtof(sd33*detInv)
		    })
		end
	  end

    fun toVector (M v) = v

  end

root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0