packages feed

splines (empty) → 0.1

raw patch · 11 files changed

+656/−0 lines, 11 filesdep +basedep +containersdep +vector-spacesetup-changed

Dependencies added: base, containers, vector-space

Files

+ Setup.lhs view
@@ -0,0 +1,5 @@+#!/usr/bin/env runhaskell++> import Distribution.Simple+> main = defaultMain+
+ splines.cabal view
@@ -0,0 +1,36 @@+name:                   splines+version:                0.1+stability:              provisional++cabal-version:          >= 1.6+build-type:             Simple++author:                 James Cook <mokus@deepbondi.net>+maintainer:             James Cook <mokus@deepbondi.net>+license:                PublicDomain++category:               Graphics, Numerical, Math+synopsis:               B-Splines, other splines, and NURBS.+description:            This is a fairly simple implementation of a +                        general-purpose spline library, just to get the code+                        out there.  Its interface is still mildly unstable and+                        may change (hopefully not drastically) as new needs or+                        better style ideas come up.  Patches, suggestions+                        and/or feature requests are welcome.++source-repository head+    type: darcs+    location: http://code.haskell.org/~mokus/splines/++Library+  hs-source-dirs:       src+  exposed-modules:      Math.Spline+                        Math.Spline.BezierCurve+                        Math.Spline.BSpline+                        Math.Spline.Class+                        Math.Spline.ISpline+                        Math.Spline.Knots+                        Math.Spline.MSpline+                        Math.NURBS+  other-modules:        Math.Spline.BSpline.Internal+  build-depends:        base >= 3 && < 5, containers, vector-space
+ src/Math/NURBS.hs view
@@ -0,0 +1,79 @@+{-# LANGUAGE StandaloneDeriving, FlexibleContexts, UndecidableInstances, TypeFamilies #-}+module Math.NURBS+    ( NURBS+    , nurbs, toNURBS+    , evalNURBS, nurbsDomain+    , nurbsDegree, nurbsKnotVector, nurbsControlPoints+    , splitNURBS+    ) where++import Data.VectorSpace+import Math.Spline.Class (Spline, toBSpline)+import Math.Spline.BSpline.Internal+import Math.Spline.BSpline+import Math.Spline.Knots++newtype NURBS v = NURBS (BSpline (Scalar v, v))++deriving instance (Eq   v, Eq   (Scalar v), Eq   (Scalar (Scalar v))) => Eq   (NURBS v)+deriving instance (Ord  v, Ord  (Scalar v), Ord  (Scalar (Scalar v))) => Ord  (NURBS v)+instance (Show v, Show (Scalar v), Show (Scalar (Scalar v))) => Show (NURBS v) where+    showsPrec p (NURBS spline) = showParen (p>11)+        ( showString "nurbs "+        . showsPrec 11 spline+        )++toNURBS :: (Spline s v, Scalar v ~ Scalar (Scalar v)) => s v -> NURBS v+toNURBS = NURBS . mapControlPoints (\p -> (1,p)) . toBSpline++nurbs :: (VectorSpace v, Scalar v ~ w,+          VectorSpace w, Scalar w ~ w)+       => Knots (Scalar v) -> [(w, v)] -> NURBS v+nurbs kts cps = NURBS (bSpline kts cps)++-- |Constructs the homogeneous-coordinates B-spline that corresponds to this+-- NURBS curve+nurbsAsSpline (NURBS spline) = spline +    { controlPoints = map homogenize (controlPoints spline) }+    where+        homogenize (w,v) = (w, w *^ v)++-- |Constructs the NURBS curve corresponding to a homogeneous-coordinates B-spline+splineAsNURBS spline = NURBS spline +    { controlPoints = map unHomogenize (controlPoints spline) }+    where+        unHomogenize (w,v) = (w, recip w *^ v)+++evalNURBS+  :: (VectorSpace v, Scalar v ~ w,+      VectorSpace w, Scalar w ~ w,+      Fractional w, Ord w) =>+     NURBS v -> w -> v+evalNURBS nurbs = project . evalBSpline (nurbsAsSpline nurbs)+    where+        project (w,v) = recip w *^ v+++-- |Returns the domain of a NURBS - that is, the range of parameter values+-- over which a spline with this degree and knot vector has a full basis set.+nurbsDomain :: Scalar v ~ Scalar (Scalar v) => +    NURBS v -> Maybe (Scalar v, Scalar v)+nurbsDomain (NURBS spline) = knotDomain (knotVector spline) (degree spline)++nurbsDegree :: NURBS v -> Int+nurbsDegree (NURBS spline) = degree spline++nurbsKnotVector :: Scalar v ~ Scalar (Scalar v) => NURBS v -> Knots (Scalar v)+nurbsKnotVector (NURBS spline) = knotVector spline++nurbsControlPoints :: NURBS v -> [(Scalar v, v)]+nurbsControlPoints (NURBS spline) = controlPoints spline++splitNURBS :: (VectorSpace v, Scalar v ~ w,+               VectorSpace w, Scalar w ~ w,+               Ord w, Fractional w)+    => NURBS v -> Scalar v -> Maybe (NURBS v, NURBS v)+splitNURBS nurbs t = do+    (s0, s1) <- splitBSpline (nurbsAsSpline nurbs) t+    return (splineAsNURBS s0, splineAsNURBS s1)
+ src/Math/Spline.hs view
@@ -0,0 +1,17 @@+module Math.Spline+    ( Spline(..)+    +    , Knots, mkKnots, knots+    +    , BezierCurve, bezierCurve+    , BSpline, bSpline+    , MSpline, mSpline, toMSpline+    , ISpline, iSpline, toISpline+    ) where++import Math.Spline.Class+import Math.Spline.Knots+import Math.Spline.BezierCurve+import Math.Spline.BSpline+import Math.Spline.MSpline+import Math.Spline.ISpline
+ src/Math/Spline/BSpline.hs view
@@ -0,0 +1,91 @@+{-# LANGUAGE MultiParamTypeClasses, StandaloneDeriving, FlexibleContexts, UndecidableInstances, TypeFamilies, ParallelListComp #-}+module Math.Spline.BSpline+    ( BSpline+    , bSpline+    , evalBSpline+    , insertKnot+    , splitBSpline+    , differentiateBSpline, integrateBSpline+    ) where++import Math.Spline.Knots+import Math.Spline.BSpline.Internal++import Data.Maybe (fromMaybe)+import Data.VectorSpace++-- |@bSpline kts cps@ creates a B-spline with the given knot vector and control +-- points.  The degree is automatically inferred as the difference between the +-- number of spans in the knot vector (@numKnots kts - 1@) and the number of +-- control points (@length cps@).+bSpline :: Knots (Scalar a) -> [a] -> BSpline a+bSpline   _  [] = error "bSpline: no control points"+bSpline kts cps = fromMaybe (error "bSpline: too few knots") (maybeSpline kts cps)++maybeSpline :: Knots (Scalar a) -> [a] -> Maybe (BSpline a)+maybeSpline kts cps +    | n > m     = Nothing+    | otherwise = Just (Spline (m - n) kts cps)+    where+        n = length cps+        m = numKnots kts - 1++deriving instance (Eq   (Scalar v), Eq   v) => Eq   (BSpline v)+deriving instance (Ord  (Scalar v), Ord  v) => Ord  (BSpline v)+instance (Show (Scalar v), Show v) => Show (BSpline v) where+    showsPrec p (Spline _ kts cps) = showParen (p>10) +        ( showString "bSpline "+        . showsPrec 11 kts+        . showChar ' '+        . showsPrec 11 cps+        )++differentiateBSpline+  :: (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => BSpline v -> BSpline v+differentiateBSpline spline+    | numKnots ks  < 2  = spline+    | numKnots ks == 2  = bSpline ks [zeroV]+    | otherwise         = bSpline ks' ds'+    where+        ks' = mkKnots . init . tail $ ts+        ds' = zipWith (*^) (tail cs) (zipWith (^-^) (tail ds) ds)+        +        ks = knotVector spline; ts = knots ks+        ds = controlPoints spline+        +        p  = degree spline+        cs = [fromIntegral p / (t1 - t0) | (t0,t1) <- spans p ts]++integrateBSpline+  :: (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => BSpline v -> BSpline v+integrateBSpline spline = bSpline (mkKnots ts') (scanl (^+^) zeroV ds')+    where+        ds' = zipWith (*^) cs (controlPoints spline)+        ts = knots (knotVector spline)+        ts' = head ts : ts ++ [last ts]+        p = degree spline + 1+        cs = [(t1 - t0) / fromIntegral p | (t0,t1) <- spans p ts]++spans n xs = zip xs (drop n xs)++splitBSpline+  :: (VectorSpace v, Ord (Scalar v), Fractional (Scalar v)) =>+     BSpline v -> Scalar v -> Maybe (BSpline v, BSpline v)+splitBSpline spline@(Spline p kv ds) t +    | inDomain  = Just (Spline p (mkKnots us0) ds0, Spline p (mkKnots us1) ds1)+    | otherwise = Nothing+    where+        inDomain = case knotDomain kv p of+            Nothing         -> False+            Just (t0, t1)   -> t >= t0 || t <= t1+        +        us = knots kv+        dss = deBoor spline t+        +        us0 = takeWhile (<t) us ++ replicate (p+1) t+        ds0 = trimTo (drop (p+1) us0) (map head dss)+        +        us1 = replicate (p+1) t ++ dropWhile (<=t) us+        ds1 = reverse (trimTo (drop (p+1) us1) (map last dss))++        trimTo list  xs = zipWith const xs list
+ src/Math/Spline/BSpline/Internal.hs view
@@ -0,0 +1,70 @@+{-# LANGUAGE FlexibleContexts #-}+module Math.Spline.BSpline.Internal+    (BSpline(..), mapControlPoints, evalBSpline, insertKnot, deBoor) where++import Math.Spline.Knots++import Data.List (zipWith4)+import Data.Monoid+import Data.VectorSpace++data BSpline v = Spline+    { degree        :: !Int+    , knotVector    :: Knots (Scalar v)+    , controlPoints :: [v]+    }++mapControlPoints f spline = spline+    { controlPoints = map f (controlPoints spline)+    , knotVector = knotVector spline+    }++evalBSpline spline = head . last . deBoor spline++-- |Insert one knot into a 'BSpline'+insertKnot+  :: (VectorSpace a, Ord (Scalar a), Fractional (Scalar a)) =>+     BSpline a -> Scalar a -> BSpline a+insertKnot spline x = spline+    { knotVector    = knotVector spline `mappend` knot x+    , controlPoints = zipWith4 (interp x) us (drop p us) ds (tail ds)+    }+    where+        us = knots (knotVector spline)+        p  = degree spline+        ds = extend (controlPoints spline)+++-- duplicate the endpoints of a list; for example,+-- extend [1..5] -> [1,1,2,3,4,5,5]+extend []       = []+extend (x:xs)   = x : extend' x xs+    where   extend' x []      = [x,x]+            extend' x (x':xs) = x:   extend' x' xs++deBoor spline x = go us (controlPoints spline)+    where+        us = knots (knotVector spline)+        +        -- Upper endpoints of the intervals are the same for+        -- each row in the table (they just line up differently+        -- with the lower endpoints):+        uHi = drop (degree spline + 1) us+        +        -- On each pass, the lower endpoints of the +        -- interpolation intervals advance and the new +        -- coefficients are given by linear interpolation+        -- on the current intervals:+        go       _ [] = []+        go (_:uLo) ds = ds : go uLo ds'+            where+                ds' = zipWith4 (interp x) uLo uHi+                                          ds (tail ds)++interp x x0 x1 y0 y1+    |  x <  x0  = y0+    |  x >= x1  = y1+    | otherwise = lerp y0 y1 a+    where+        a = (x - x0) / (x1 - x0)+
+ src/Math/Spline/BezierCurve.hs view
@@ -0,0 +1,57 @@+{-# LANGUAGE MultiParamTypeClasses, FlexibleInstances, UndecidableInstances #-}+module Math.Spline.BezierCurve+    ( BezierCurve, bezierCurve, splitBezierCurve+    , evalSpline+    ) where++import Math.Spline.BSpline+import Math.Spline.Class+import Math.Spline.Knots++import Control.Applicative+import Data.VectorSpace++-- |A BezierCurve curve on @0 <= x <= 1@.+data BezierCurve v = BezierCurve !Int [v] deriving (Eq, Ord)++-- |Construct a Bezier curve from a list of control points.  The degree+-- of the curve is one less than the number of control points.+bezierCurve :: [v] -> BezierCurve v+bezierCurve cs+    | null cs   = error "bezierCurve: no control points given"+    | otherwise = BezierCurve (length cs - 1) cs++instance Show v => Show (BezierCurve v) where+    showsPrec p (BezierCurve _ cs) = showParen (p>10)+        ( showString "bezierCurve "+        . showsPrec 11 cs+        )++instance (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => Spline BezierCurve v where+    splineDomain (BezierCurve _  _) = Just (0,1)+    evalSpline   (BezierCurve _ cs) = head . last . deCasteljau cs+    splineDegree (BezierCurve p  _) = p+    knotVector   (BezierCurve p  _) = fromList [(0, p+1), (1, p+1)]+    toBSpline = bSpline <$> knotVector <*> controlPoints++instance Spline BezierCurve v => ControlPoints BezierCurve v where+    controlPoints (BezierCurve _ cs) = cs++deCasteljau :: VectorSpace v => [v] -> Scalar v -> [[v]]+deCasteljau [] t = []+deCasteljau cs t = cs : deCasteljau (zipWith interp cs (tail cs)) t+    where+        interp x0 x1 = lerp x0 x1 t++-- |Split and rescale a Bezier curve.  Given a 'BezierCurve' @b@ and a point +-- @t@, @splitBezierCurve b t@ creates 2 curves @(b1, b2)@ such that (up to +-- reasonable numerical accuracy expectations):+-- +-- > evalSpline b1  x    == evalSpline b (x * t)+-- > evalSpline b2 (x-t) == evalSpline b (x * (1-t))+-- +splitBezierCurve :: VectorSpace v => BezierCurve v -> Scalar v -> (BezierCurve v, BezierCurve v)+splitBezierCurve (BezierCurve n cs) t = +    ( BezierCurve n (map head css)+    , BezierCurve n (reverse (map last css))+    ) where css = deCasteljau cs t
+ src/Math/Spline/Class.hs view
@@ -0,0 +1,44 @@+{-# LANGUAGE MultiParamTypeClasses, FlexibleInstances, FlexibleContexts, UndecidableInstances #-}+module Math.Spline.Class where++import Control.Applicative+import Math.Spline.Knots+import qualified Math.Spline.BSpline.Internal as BSpline++import Data.VectorSpace++-- |A spline is a piecewise polynomial vector-valued function.  The necessary+-- and sufficient instance definition is 'toBSpline'.+class (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => Spline s v where+    -- |Returns the domain of a spline.  In the case of B-splines, this is+    -- the domain on which a spline with this degree and knot vector has a +    -- full basis set.  In other cases, it should be no larger than +    -- @splineDomain . toBSpline@, but may be smaller.  Within this domain,+    -- 'evalSpline' should agree with @'evalSpline' . 'toBSpline'@ (not +    -- necessarily exactly, but up to reasonable expectations of numerical +    -- accuracy).+    splineDomain :: s v -> Maybe (Scalar v, Scalar v)+    splineDomain = knotDomain <$> knotVector <*> splineDegree+    +    evalSpline :: s v -> Scalar v -> v+    evalSpline = evalSpline . toBSpline+    +    splineDegree :: s v -> Int+    splineDegree = splineDegree . toBSpline+    +    knotVector :: s v -> Knots (Scalar v)+    knotVector = knotVector . toBSpline+    +    toBSpline :: s v -> BSpline.BSpline v++class Spline s v => ControlPoints s v where+    controlPoints :: s v -> [v]++instance (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => Spline BSpline.BSpline v where+    evalSpline spline = head . last . BSpline.deBoor spline+    splineDegree = BSpline.degree+    knotVector = BSpline.knotVector+    toBSpline = id++instance Spline BSpline.BSpline v => ControlPoints BSpline.BSpline v where+    controlPoints = BSpline.controlPoints
+ src/Math/Spline/ISpline.hs view
@@ -0,0 +1,74 @@+{-# LANGUAGE+        MultiParamTypeClasses,+        FlexibleInstances, FlexibleContexts, UndecidableInstances,+        ParallelListComp,+        StandaloneDeriving+  #-}+module Math.Spline.ISpline+    ( ISpline, iSpline, toISpline+    , evalSpline+    ) where++import Math.Spline.BSpline+import Math.Spline.Class+import Math.Spline.Knots++import Data.VectorSpace++-- |The I-Spline basis functions are the integrals of the M-splines, or+-- alternatively the integrals of the B-splines normalized to the range+-- [0,1].  Every I-spline basis function increases monotonically from 0 to 1,+-- thus it is useful as a basis for monotone functions.  An I-Spline curve+-- is monotone if and only if every non-zero control point has the same sign.+data ISpline v = ISpline+    { iSplineDegree        :: !Int+    , iSplineKnotVector    :: Knots (Scalar v)+    , iSplineControlPoints :: [v]+    }++deriving instance (Eq   (Scalar v), Eq   v) => Eq   (ISpline v)+deriving instance (Ord  (Scalar v), Ord  v) => Ord  (ISpline v)+instance (Show (Scalar v), Show v) => Show (ISpline v) where+    showsPrec p (ISpline _ kts cps) = showParen (p>10) +        ( showString "iSpline "+        . showsPrec 11 kts+        . showChar ' '+        . showsPrec 11 cps+        )+++-- |@iSpline kts cps@ creates an I-spline with the given knot vector and control +-- points.  The degree is automatically inferred as the difference between the +-- number of spans in the knot vector (@numKnots kts - 1@) and the number of +-- control points (@length cps@).+iSpline :: Knots (Scalar a) -> [a] -> ISpline a+iSpline kts cps +    | n > m     = error "iSpline: too few knots"+    | otherwise = ISpline (m - n) kts cps+    where+        n = length cps+        m = numKnots kts - 1++instance (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => Spline ISpline v where+    splineDegree = (1 +) . iSplineDegree+    knotVector spline = mkKnots (head ts : ts ++ [last ts])+        where ts = knots (iSplineKnotVector spline)+    toBSpline spline = bSpline (knotVector spline) (scanl (^+^) zeroV cs)+        where cs = iSplineControlPoints spline++instance Spline ISpline v => ControlPoints ISpline v where+    controlPoints      (ISpline _  _ cs) = cs++toISpline :: (Spline s v, Eq v) => s v -> ISpline v+toISpline = fromBSpline . toBSpline++fromBSpline spline+    | head ds == zeroV +    && numKnots ks >= 2 = iSpline (mkKnots (init (tail ts))) (tail ds')+    | otherwise         = iSpline (mkKnots (init       ts )) ds'+    where+        ks = knotVector spline+        ts = knots ks+        ds = controlPoints spline+        +        ds' = zipWith (^-^) ds (zeroV:ds)
+ src/Math/Spline/Knots.hs view
@@ -0,0 +1,108 @@+module Math.Spline.Knots+    ( Knots+    , knot, multipleKnot+    , mkKnots, fromList+    +    , knots, numKnots+    , toList, distinctKnots, numDistinctKnots+    +    , knotMultiplicity, setKnotMultiplicity+    +    , knotDomain+    ) where++import Prelude hiding (sum)+import Data.Foldable (Foldable(foldMap), sum)+import qualified Data.Map as M+import Data.Monoid (Monoid(..))+import Data.Maybe (fromMaybe)++-- |Knot vectors - multisets of points in a 1-dimensional space.+data Knots a = Knots !Int (M.Map a Int) deriving (Eq, Ord)++instance Show a => Show (Knots a) where+    showsPrec p ks@(Knots 1 _) = showParen (p > 10)+        ( showString "knot "+        . showsPrec 11 (head $ knots ks)+        )+    showsPrec p ks = showParen (p > 10)+        ( showString "mkKnots "+        . showsPrec 11 (knots ks)+        )++instance (Ord a) => Monoid (Knots a) where+    mempty = Knots 0 M.empty+    mappend (Knots n1 v1) (Knots n2 v2) =+        Knots (n1 + n2) (M.filter (/=0) (M.unionWith (+) v1 v2))++instance Foldable Knots where+    foldMap f = foldMap f . knots+++-- |Create a knot vector consisting of one knot.+knot :: Ord a => a -> Knots a+knot x = multipleKnot x 1++-- |Create a knot vector consisting of one knot with the specified multiplicity.+multipleKnot :: Ord a => a -> Int -> Knots a+multipleKnot k n +    | n <= 0    = Knots 0 (M.empty)+    | otherwise = Knots n (M.singleton k n)++-- |Create a knot vector consisting of all the knots in a list.+mkKnots :: (Ord a) => [a] -> Knots a+mkKnots ks = fromList (map (\k -> (k,1)) ks)++-- |Create a knot vector consisting of all the knots and corresponding +-- multiplicities in a list.+fromList :: (Ord k) => [(k, Int)] -> Knots k+fromList ks = Knots (sum kMap) kMap+    where kMap = M.fromListWith (+) (filter ((>0).snd) ks)++-- |Returns a list of all distinct knots in ascending order along with+-- their multiplicities.+toList :: Knots k -> [(k, Int)]+toList (Knots _ ks) = M.toList ks++-- |Returns the number of knots (not necessarily distinct) in a knot vector.+numKnots :: Knots t -> Int+numKnots (Knots n _) = n++-- |Returns the number of distinct knots in a knot vector.+numDistinctKnots :: Knots t -> Int+numDistinctKnots (Knots _ ks) = M.size ks++-- |Returns a list of all knots (not necessarily distinct) of a knot vector in ascending order+knots :: Knots t -> [t]+knots (Knots _ ks) = concat [replicate n k | (k,n) <- M.toAscList ks]++-- |Returns a list of all distinct knots of a knot vector in ascending order+distinctKnots :: Knots t -> [t]+distinctKnots (Knots _ ks) = M.keys ks++-- |Looks up the multiplicity of a knot (which is 0 if the point is not a knot)+knotMultiplicity :: (Ord k) => k -> Knots k -> Int+knotMultiplicity k (Knots _ ks) = fromMaybe 0 (M.lookup k ks)++-- |Returns a new knot vector with the given knot set to the specified +-- multiplicity and all other knots unchanged.+setKnotMultiplicity :: Ord k => k -> Int -> Knots k -> Knots k+setKnotMultiplicity k n (Knots m ks)+    | n <= 0    = Knots (m     - n') (M.delete k ks)+    | otherwise = Knots (m + n - n') (M.insert k n ks)+    where+        n' = knotMultiplicity k (Knots m ks)++-- |@knotDomain kts p@ return the domain of a B-spline or NURBS with knot+-- vector @kts@ and degree @p@.  This is the subrange spanned by all+-- except the first and last @p@ knots.  Outside this domain, the spline+-- does not have a complete basis set.  De Boor's algorithm assumes that+-- the basis functions sum to 1, which is only true on this range, and so+-- this is also precisely the domain on which de Boor's algorithm is valid.+knotDomain :: Knots a -> Int -> Maybe (a,a)+knotDomain ks@(Knots n _) p +    | n > 2*p   = Just (head (drop p kts), head (drop p (reverse kts)))+    | otherwise = Nothing+    where+        kts = knots ks+
+ src/Math/Spline/MSpline.hs view
@@ -0,0 +1,75 @@+{-# LANGUAGE+        MultiParamTypeClasses,+        FlexibleInstances, FlexibleContexts, UndecidableInstances,+        ParallelListComp,+        StandaloneDeriving+  #-}+module Math.Spline.MSpline+    ( MSpline, mSpline, toMSpline+    , evalSpline+    ) where++import Math.Spline.BSpline+import Math.Spline.Class+import Math.Spline.Knots++import Data.VectorSpace++-- |M-Splines are B-splines normalized so that the integral of each basis +-- function over the spline domain is 1.+data MSpline v = MSpline+    { mSplineDegree        :: !Int+    , mSplineKnotVector    :: Knots (Scalar v)+    , mSplineControlPoints :: [v]+    }++deriving instance (Eq   (Scalar v), Eq   v) => Eq   (MSpline v)+deriving instance (Ord  (Scalar v), Ord  v) => Ord  (MSpline v)+instance (Show (Scalar v), Show v) => Show (MSpline v) where+    showsPrec p (MSpline _ kts cps) = showParen (p>10) +        ( showString "mSpline "+        . showsPrec 11 kts+        . showChar ' '+        . showsPrec 11 cps+        )+++-- |@mSpline kts cps@ creates a M-spline with the given knot vector and control +-- points.  The degree is automatically inferred as the difference between the +-- number of spans in the knot vector (@numKnots kts - 1@) and the number of +-- control points (@length cps@).+mSpline :: Knots (Scalar a) -> [a] -> MSpline a+mSpline kts cps+    | n > m     = error "mSpline: too few knots"+    | otherwise = MSpline (m - n) kts cps+    where+        n = length cps+        m = numKnots kts - 1++spans n xs = zip xs (drop n xs)++instance (VectorSpace v, Fractional (Scalar v), Ord (Scalar v)) => Spline MSpline v where+    splineDegree = mSplineDegree+    knotVector   = mSplineKnotVector+    toBSpline (MSpline p ks cs) = bSpline ks cs'+        where+            n = p + 1; n' = fromIntegral n+            cs' = [ (n' / (t1 - t0)) *^ c +                  | c <- cs+                  | (t0, t1) <- spans n (knots ks)+                  ]++instance Spline MSpline v => ControlPoints MSpline v where+    controlPoints = mSplineControlPoints++toMSpline :: Spline s v => s v -> MSpline v+toMSpline = fromBSpline . toBSpline++fromBSpline spline = mSpline ks cs+    where+        n = splineDegree spline + 1; n' = fromIntegral n+        ks = knotVector spline+        cs =  [ ((t1 - t0) / n') *^ c+              | c <- controlPoints spline+              | (t0, t1) <- spans n (knots ks)+              ]