\section{Points and Vectors: RSAGL.Vector} \begin{code}
{-# LANGUAGE FlexibleInstances, MultiParamTypeClasses, FunctionalDependencies, PatternGuards #-}
module RSAGL.Math.Vector
(Point3D(..),
origin_point_3d,
Vector3D(..),
SurfaceVertex3D(..),
zero_vector,
point3d,
points3d,
point2d,
points2d,
vector3d,
dotProduct,
angleBetween,
crossProduct,
distanceBetween,
distanceBetweenSquared,
aNonZeroVector,
aLargeVector,
vectorSum,
vectorScale,
vectorScaleTo,
vectorToFrom,
vectorNormalize,
vectorAverage,
vectorLength,
vectorLengthSquared,
newell,
Xyz(..),
XYZ,
vectorString,
randomXYZ,
fixOrtho,
fixOrtho2,
fixOrtho2Left,
orthos)
where

import Control.DeepSeq
import Control.Parallel.Strategies
import RSAGL.Math.Angle
import System.Random
import RSAGL.Math.AbstractVector
import RSAGL.Math.Types
import RSAGL.Math.ListUtils

\end{code} \subsection{Generic 3-dimensional types and operations} \begin{code}
type XYZ = (RSdouble,RSdouble,RSdouble)

class Xyz a where
toXYZ :: a -> XYZ
fromXYZ :: XYZ -> a

instance Xyz (RSdouble,RSdouble,RSdouble) where
toXYZ = id
fromXYZ = id

vectorString :: Xyz a => a -> String
vectorString xyz = let (x,y,z) = toXYZ xyz
in (show x) ++ "," ++ (show y) ++ "," ++ (show z)

uncurry3d :: (RSdouble -> RSdouble -> RSdouble -> a) -> XYZ -> a
uncurry3d fn (x,y,z) = fn x y z

\end{code} \subsection{Points in 3-space} \begin{code}
data Point3D = Point3D {-# UNPACK #-} !RSdouble {-# UNPACK #-} !RSdouble {-# UNPACK #-} !RSdouble

origin_point_3d :: Point3D
origin_point_3d = Point3D 0 0 0

point3d :: (RSdouble,RSdouble,RSdouble) -> Point3D
point3d = uncurry3d Point3D

point2d :: (RSdouble,RSdouble) -> Point3D
point2d (x,y) = point3d (x,y,0)

points3d :: [(RSdouble,RSdouble,RSdouble)] -> [Point3D]
points3d = map point3d

points2d :: [(RSdouble,RSdouble)] -> [Point3D]
points2d = map point2d

instance Xyz Point3D where
toXYZ (Point3D x y z) = (x,y,z)
fromXYZ (x,y,z) = Point3D x y z

instance AbstractZero Point3D where
zero = origin_point_3d

instance AbstractSubtract Point3D Vector3D where
sub = vectorToFrom

instance NFData Point3D

\end{code} \subsection{Vectors in 3-space} \begin{code}
data Vector3D = Vector3D {-# UNPACK #-} !RSdouble {-# UNPACK #-} !RSdouble {-# UNPACK #-} !RSdouble

zero_vector :: Vector3D
zero_vector = Vector3D 0 0 0

vector3d :: (RSdouble,RSdouble,RSdouble) -> Vector3D
vector3d = uncurry3d Vector3D

instance Xyz Vector3D where
toXYZ (Vector3D x y z) = (x,y,z)
fromXYZ (x,y,z) = Vector3D x y z

instance NFData Vector3D

instance AbstractZero Vector3D where
zero = zero_vector

instance AbstractSubtract Vector3D Vector3D where
sub = vectorToFrom

instance AbstractScale Vector3D where
scalarMultiply = vectorScale

instance AbstractMagnitude Vector3D where
magnitude = vectorLength

instance AbstractVector Vector3D

\end{code} A \texttt{SurfaceVertex3D} is a a point on an orientable surface, having a position and a normal vector. \subsection{Surface Vertices in 3-space} \begin{code}
data SurfaceVertex3D = SurfaceVertex3D { sv3d_position :: Point3D,
sv3d_normal :: Vector3D }

instance NFData SurfaceVertex3D where
rnf (SurfaceVertex3D p v) = rnf (p,v)

\end{code} \subsection{Vector Arithmetic} \begin{code}
aNonZeroVector :: Vector3D -> Maybe Vector3D
aNonZeroVector v = case vectorLength v of
x | x <= 0 -> Nothing
x | isDenormalized x -> Nothing
x | isNaN x -> Nothing
x | isInfinite x -> Nothing
_ | otherwise -> Just v

aLargeVector :: RSdouble -> Vector3D -> Maybe Vector3D
aLargeVector x v_ =
case aNonZeroVector v_ of
Just v | vectorLength v > x -> Just v
_ | otherwise -> Nothing

dotProduct :: Vector3D -> Vector3D -> RSdouble
dotProduct (Vector3D ax ay az) (Vector3D bx by bz) =
(ax*bx) + (ay*by) + (az*bz)

angleBetween :: Vector3D -> Vector3D -> Angle
angleBetween a b = arcCosine $dotProduct (vectorNormalize a) (vectorNormalize b) crossProduct :: Vector3D -> Vector3D -> Vector3D crossProduct (Vector3D ax ay az) (Vector3D bx by bz) = Vector3D (ay*bz - az*by) (az*bx - ax*bz) (ax*by - ay*bx) distanceBetween :: (Xyz xyz) => xyz -> xyz -> RSdouble distanceBetween a b = vectorLength$ vectorToFrom a b

distanceBetweenSquared :: (Xyz xyz) => xyz -> xyz -> RSdouble
distanceBetweenSquared a b = vectorLengthSquared $vectorToFrom a b displace :: (Xyz xyz) => xyz -> Vector3D -> xyz displace xyz (Vector3D x2 y2 z2) = let (x1,y1,z1) = toXYZ xyz in fromXYZ (x1+x2,y1+y2,z1+z2) vectorAdd :: Vector3D -> Vector3D -> Vector3D vectorAdd = displace vectorSum :: [Vector3D] -> Vector3D vectorSum vectors = foldr vectorAdd zero_vector vectors vectorToFrom :: (Xyz xyz) => xyz -> xyz -> Vector3D vectorToFrom a b = let (ax,ay,az) = toXYZ a (bx,by,bz) = toXYZ b in Vector3D (ax - bx) (ay - by) (az - bz) vectorLength :: Vector3D -> RSdouble vectorLength = sqrt . vectorLengthSquared vectorLengthSquared :: Vector3D -> RSdouble vectorLengthSquared (Vector3D x y z) = (x*x + y*y + z*z) vectorScale :: RSdouble -> Vector3D -> Vector3D vectorScale s (Vector3D x y z) = Vector3D (x*s) (y*s) (z*s)  \end{code} vectorScaleTo forces the length of a vector to a certain value, without changing the vector's direction. \begin{code} vectorScaleTo :: RSdouble -> Vector3D -> Vector3D vectorScaleTo new_length vector = vectorScale new_length$ vectorNormalize vector

\end{code} vectorNormalize forces the length of a vector to 1, without changing the vector's direction. \begin{code}
vectorNormalize :: Vector3D -> Vector3D
vectorNormalize v =
let l = vectorLength v
in maybe (Vector3D 0 0 0) (vectorScale (1/l)) $aNonZeroVector v  \end{code} vectorAverage takes the average of a list of vectors. The result is a normalized vector where the length of the element vectors is not reflected in the result. \begin{code} vectorAverage :: [Vector3D] -> Vector3D vectorAverage vects = vectorNormalize$ vectorSum $map vectorNormalize vects  \end{code} \subsection{Generating normal vectors} \texttt{newell} calculates the normal vector of an arbitrary polygon. If the points specified are non-coplanar, \texttt{newell} often calculates a reasonable result; if they are colinear or singular \texttt{newell} will fail. The result is a normalized vector. \begin{code} newell :: [Point3D] -> Maybe Vector3D newell points = fmap vectorNormalize$ aNonZeroVector $vectorSum$ map newell_ $loopedDoubles points where newell_ (Point3D x0 y0 z0,Point3D x1 y1 z1) = (Vector3D ((y0 - y1)*(z0 + z1)) ((z0 - z1)*(x0 + x1)) ((x0 - x1)*(y0 + y1)))  \end{code} \subsection{Randomly Generated Coordinates} \texttt{randomXYZ} can generate random coordinates within the cube where x, y, and z are each in the range (lo,hi). \begin{code} randomXYZ :: (RandomGen g,Xyz p) => (RSdouble,RSdouble) -> g -> (p,g) randomXYZ (lo,hi) g = (fromXYZ (f2f x,f2f y,f2f z),g') where (g_,g') = split g (x:y:z:_) = randomRs (f2f lo,f2f hi) g_ :: [Double]  \end{code} \subsection{Orthagonal Vectors} \texttt{fixOrtho a v} finds the vector, orthogonal to a, that has the least angle to v. \texttt{fixOrtho2 right up} yields \texttt{(up,forward)}. \texttt{fixOrtho2Left right up} yields \texttt{(up,backward)}. \texttt{orthos} finds two arbitrary vectors orthogonal to the parameter. \begin{code} fixOrtho :: Vector3D -> Vector3D -> Vector3D fixOrtho a = fst . fixOrtho2 a fixOrtho2 :: Vector3D -> Vector3D -> (Vector3D,Vector3D) fixOrtho2 a v = (vectorNormalize$ crossProduct a $vectorScale (-1) b,vectorNormalize b) where b = crossProduct a v fixOrtho2Left :: Vector3D -> Vector3D -> (Vector3D,Vector3D) fixOrtho2Left a v = (vectorNormalize$ crossProduct a b,vectorNormalize b)
where b = vectorScale (-1) $crossProduct a v orthos :: Vector3D -> (Vector3D,Vector3D) orthos v@(Vector3D x y z) | abs y >= abs x && abs z >= abs x = fixOrtho2 v (Vector3D (abs x + abs y + abs z) y z) orthos v@(Vector3D x y z) | abs x >= abs y && abs z >= abs y = fixOrtho2 v (Vector3D x (abs x + abs y + abs z) z) orthos v@(Vector3D x y z) | abs x >= abs z && abs y >= abs z = fixOrtho2 v (Vector3D x y (abs x + abs y + abs z)) orthos v = error$ "orthos: (" ++ show v ++ ")"

\end{code}