interpolation 0.0 → 0.1

raw patch · 14 files changed

+430/−54 lines, 14 filesdep +arraydep +containersdep +hmatrix-bandedPVP ok

version bump matches the API change (PVP)

Dependencies added: array, containers, hmatrix-banded

API changes (from Hackage documentation)

Files

• ChangeLog +9/−0
• example/Fit.hs +104/−10
• example/Plot.hs +3/−3
• interpolation.cabal +15/−4
• private/Numeric/Interpolation/Private/Basis.hs +9/−1
• private/Numeric/Interpolation/Private/List.hs +27/−0
• src/Numeric/Interpolation/Basis.hs +1/−8
• src/Numeric/Interpolation/Basis/Compact.hs +8/−20
• src/Numeric/Interpolation/NodeList.hs +26/−0
• src/Numeric/Interpolation/Sample.hs +113/−0
• src/Numeric/Interpolation/Type.hs +23/−7
• test/Test.hs +5/−1
• test/Test/Overlap.hs +25/−0
• test/Test/Sample.hs +62/−0

@@ -0,0 +1,9 @@+0.1:++* Hermite1 interpolation: different order of coefficients++  Interleave node values and derivatives+  in order to get a narrow banded Gramian matrix+  from a sampled interpolation basis.++  Rename from 'cubic' to 'hermite1'.

@@ -4,8 +4,13 @@ import qualified Numeric.Interpolation.Piecewise as Piecewise import qualified Numeric.Interpolation.Type as Type +import qualified Data.Packed.ST as PackST import qualified Data.Packed.Matrix as Matrix import qualified Data.Packed.Vector as Vector+import Data.Packed.Matrix (Matrix)+import Data.Packed.Vector (Vector)++import qualified Numeric.LinearAlgebra.Banded as Banded import qualified Numeric.Container as Container import Numeric.Container ((<\>)) @@ -15,6 +20,9 @@  import System.Random (randomRs, mkStdGen) import Control.Monad.HT (void)+import Control.Monad (when, zipWithM_, forM_)++import qualified Data.Foldable as Fold import Data.Monoid ((<>))  @@ -26,19 +34,93 @@       (randomRs (0,7) (mkStdGen 23))       (randomRs (-0.2,0.2) (mkStdGen 42)) +basisMatrixFull ::+   Type.T Double Double ny -> [Double] -> [Double] -> Matrix Double+basisMatrixFull typ xs txs0 =+   let txs = Vector.fromList txs0+   in  Matrix.fromColumns $+       map (flip Container.cmap txs . Piecewise.interpolateConstantExt typ) $+       Type.basisFunctions typ xs++basisMatrixSparse ::+   Type.T Double Double ny -> [Double] -> [Double] -> Matrix Double+basisMatrixSparse typ xs txs = PackST.runSTMatrix $ do+   mat <- PackST.newMatrix 0 (length txs) (length $ Type.basisFunctions typ xs)+   zipWithM_+      (\k -> mapM_ (uncurry (PackST.writeMatrix mat k))) [0..] $+      map (Type.sampleBasisFunctions typ xs) txs+   return mat+ fit ::    Type.T Double Double ny ->    [Double] -> [(Double, Double)] -> Nodes.T Double ny fit typ xs target =-   let txs = Vector.fromList $ map fst target-       tys = Vector.fromList $ map snd target-       matrix =-          Matrix.fromColumns $-          map (flip Container.cmap txs . Piecewise.interpolateConstantExt typ) $-          Type.basisFunctions typ xs+   let (txs, tys) = unzip target+       matrix = basisMatrixSparse typ xs txs    in  Type.coefficientsToInterpolator typ xs $-       Vector.toList $ matrix <\> tys+       Vector.toList $ matrix <\> Vector.fromList tys +matrixDiff ::+   Type.T Double Double ny ->+   [Double] -> [(Double, Double)] -> Double+matrixDiff typ xs target =+   let txs = map fst target+   in  Container.maxElement $ Container.cmap abs $+       Container.sub+          (basisMatrixFull typ xs txs)+          (basisMatrixSparse typ xs txs)+++mulSparseMatrixVector ::+   Int -> [[(Int, Double)]] -> [Double] -> Vector Double+mulSparseMatrixVector size samples tys = PackST.runSTVector $ do+   vec <- PackST.newVector 0 size+   forM_ (zip samples tys) $ \(row,ty) ->+      forM_ row $ \(k,y) ->+         PackST.modifyVector vec k (+y*ty)+   return vec++bandedGramian ::+   Int -> Int -> [[(Int, Double)]] -> Banded.SymmetricMatrix Double+bandedGramian size width samples =+      Banded.SymmetricMatrix $ PackST.runSTMatrix $ do+   mat <- PackST.newMatrix 0 size width+   forM_ samples $ \row ->+      forM_ row $ \(k,yk) ->+      forM_ row $ \(j,yj) ->+         when (k<=j) $ PackST.modifyMatrix mat k (j-k) (+yk*yj)+   return mat++fitBanded ::+   Type.T Double Double ny ->+   [Double] -> [(Double, Double)] -> Nodes.T Double ny+fitBanded typ xs target =+   let size = length $ Type.basisFunctions typ xs+       (txs, tys) = unzip target+       samples = map (Type.sampleBasisFunctions typ xs) txs+       matrix =+          Banded.choleskyDecompose $+          bandedGramian size (Type.basisOverlap typ) samples+   in  Type.coefficientsToInterpolator typ xs $ Vector.toList $+       Banded.choleskySolve matrix $ mulSparseMatrixVector size samples tys++bandedDiff ::+   (ny -> ny -> Double) ->+   Type.T Double Double ny ->+   [Double] -> [(Double, Double)] -> Double+bandedDiff absDiff typ xs target =+   maximum $+   zipWith absDiff+      (Fold.toList $ fit typ xs target)+      (Fold.toList $ fitBanded typ xs target)++absDiffSingle :: Double -> Double -> Double+absDiffSingle x y = abs (x-y)++absDiffPair :: (Double,Double) -> (Double,Double) -> Double+absDiffPair (x,dx) (y,dy) = max (abs (x-y)) (abs (dx-dy))++ plotBasisFunctions ::    Type.T Double Double ny -> [Double] -> Plot2D.T Double Double plotBasisFunctions nodeType xs =@@ -53,11 +135,11 @@    let xs = [0, 1, 3, 4, 6, 7]        exs = (-1) : xs ++ [8]    void $ GP.plotDefault $ plotBasisFunctions Type.linear xs-   void $ GP.plotDefault $ plotBasisFunctions Type.cubic xs+   void $ GP.plotDefault $ plotBasisFunctions Type.hermite1 xs    void $ GP.plotDefault $ plotBasisFunctions Type.cubicLinear exs    void $ GP.plotDefault $ plotBasisFunctions Type.cubicParabola exs    let linearNodes = fit Type.linear xs noisy-       hermite1Nodes = fit Type.cubic xs noisy+       hermite1Nodes = fit Type.hermite1 xs noisy        cubicLinearNodes = fit Type.cubicLinear exs noisy        cubicParabolaNodes = fit Type.cubicParabola exs noisy    void $ GP.plotDefault $@@ -65,7 +147,19 @@       <>       (Plot2D.functions Graph2D.lines (Plot2D.linearScale 1000 (-2,10)) $        Piecewise.interpolateConstantExt Type.linear linearNodes :-       Piecewise.interpolateConstantExt Type.cubic hermite1Nodes :+       Piecewise.interpolateConstantExt Type.hermite1 hermite1Nodes :        Piecewise.interpolateConstantExt Type.cubicLinear cubicLinearNodes :        Piecewise.interpolateConstantExt Type.cubicParabola cubicParabolaNodes :        [])++   putStrLn "differences between matrices should be almost zero:"+   putStrLn $ "linear: " ++ show (matrixDiff Type.linear xs noisy)+   putStrLn $ "hermite1: " ++ show (matrixDiff Type.hermite1 xs noisy)+   putStrLn $ "cubicLinear: " ++ show (matrixDiff Type.cubicLinear exs noisy)+   putStrLn $ "cubicParabola: " ++ show (matrixDiff Type.cubicParabola exs noisy)++   putStrLn "differences between samples should be almost zero:"+   putStrLn $ "linear: " ++ show (bandedDiff absDiffSingle Type.linear xs noisy)+   putStrLn $ "hermite1: " ++ show (bandedDiff absDiffPair Type.hermite1 xs noisy)+   putStrLn $ "cubicLinear: " ++ show (bandedDiff absDiffPair Type.cubicLinear exs noisy)+   putStrLn $ "cubicParabola: " ++ show (bandedDiff absDiffPair Type.cubicParabola exs noisy)

@@ -27,7 +27,7 @@       Plot2D.functions Graph2D.lines          (Plot2D.linearScale 1000 (-2,15))          [Piecewise.interpolateConstantExt Type.linear linearNodes,-          Piecewise.interpolateConstantExt Type.cubic hermite1Nodes,-          Piecewise.interpolateConstantExt Type.cubic cubicLinearNodes,-          Piecewise.interpolateConstantExt Type.cubic cubicParabolaNodes,+          Piecewise.interpolateConstantExt Type.hermite1 hermite1Nodes,+          Piecewise.interpolateConstantExt Type.hermite1 cubicLinearNodes,+          Piecewise.interpolateConstantExt Type.hermite1 cubicParabolaNodes,           sin]

@@ -1,5 +1,5 @@ Name:             interpolation-Version:          0.0+Version:          0.1 License:          BSD3 License-File:     LICENSE Author:           Henning Thielemann@@ -8,7 +8,7 @@ Category:         Math Synopsis:         piecewise linear and cubic Hermite interpolation Description:-  Represent real functions by linear or cubic segments.+  Represent real functions by linear or cubic polynomial segments.   The package provides both data structures   for efficient lookup of interpolation intervals,   and computation of basis functions.@@ -30,9 +30,13 @@     We use a distorted sinus as target.   .   The package needs only Haskell 98.-Tested-With:      GHC==7.4.2, GHC==7.6.3, GHC==7.8.2+  Most of the package dependencies are only needed for the examples+  and are only installed if you enable to build them.+Tested-With:      GHC==7.4.2, GHC==7.6.3, GHC==7.8.4, GHC==7.10.1 Cabal-Version:    >=1.8 Build-Type:       Simple+Extra-Source-Files:+  ChangeLog  Flag buildExamples   description: Build example executables@@ -43,7 +47,7 @@   default:     True  Source-Repository this-  Tag:         0.0+  Tag:         0.1   Type:        darcs   Location:    http://code.haskell.org/~thielema/interpolation/ @@ -66,9 +70,11 @@     Numeric.Interpolation.Basis     Numeric.Interpolation.Basis.Compact     Numeric.Interpolation.Basis.Full+    Numeric.Interpolation.Sample   Other-Modules:     Numeric.Interpolation.Private.Piece     Numeric.Interpolation.Private.Basis+    Numeric.Interpolation.Private.List  Executable interpolation-plot   Main-Is:             Plot.hs@@ -90,6 +96,7 @@   If flag(buildExamples) && flag(hmatrix)     Build-Depends:       interpolation,+      hmatrix-banded >=0.0 && <0.1,       hmatrix >=0.15 && <0.16,       random >=1.0 && <1.1,       gnuplot >=0.5.2 && <0.6,@@ -103,10 +110,14 @@   Main-Is:             Test.hs   Other-Modules:     Test.Piece+    Test.Sample+    Test.Overlap   Hs-Source-Dirs:      test, private   GHC-Options:         -Wall   Build-Depends:     interpolation,     QuickCheck >=2.4 && <2.8,     utility-ht >=0.0.9 && <0.1,+    array >=0.4 && <0.6,+    containers >=0.4 && <0.6,     base >=4.5 && <4.8

@@ -5,8 +5,16 @@ import qualified Data.List.Match as Match  +_hermite1Split :: [a] -> [b] -> [(b, b)]+_hermite1Split xs = uncurry zip . Match.splitAt xs+ hermite1Split :: [a] -> [b] -> [(b, b)]-hermite1Split xs = uncurry zip . Match.splitAt xs+hermite1Split _ = pairs++pairs :: [a] -> [(a,a)]+pairs (x0:x1:xs) = (x0,x1) : pairs xs+pairs [] = []+pairs _ = error "pairs: odd number of elements"   parabolaDerivative ::

@@ -0,0 +1,27 @@+module Numeric.Interpolation.Private.List where++import qualified Data.List as List+++mapAdjacent3 :: (a -> a -> a -> b) -> [a] -> [b]+mapAdjacent3 f xs0 =+   let xs1 = drop 1 xs0+       xs2 = drop 1 xs1+   in  List.zipWith3 f xs0 xs1 xs2++mapAdjacentMaybe3 :: (Maybe a -> a -> Maybe a -> b) -> [a] -> [b]+mapAdjacentMaybe3 f xs =+   let jxs = map Just xs+   in  zipWith3 f (Nothing : jxs) xs (drop 1 jxs ++ [Nothing])++mapAdjacentMaybe5 ::+   (Maybe a -> Maybe a -> a -> Maybe a -> Maybe a -> b) ->+   [a] -> [b]+mapAdjacentMaybe5 f xs =+   let jxs = map Just xs+       lxs1 = Nothing : jxs+       lxs2 = Nothing : lxs1+       rxs1 = drop 1 $ jxs ++ repeat Nothing+       rxs2 = drop 1 $ rxs1+   in  List.zipWith5 f lxs2 lxs1 xs rxs1 rxs2+

@@ -13,8 +13,7 @@ import qualified Numeric.Interpolation.NodeList as Nodes import Numeric.Interpolation.Private.Basis           (parabolaDerivativeCenterNode, hermite1Split)--import qualified Data.List as List+import Numeric.Interpolation.Private.List (mapAdjacent3, )   {- |@@ -34,12 +33,6 @@    Nodes.fromList . zip xs . hermite1Split xs  --mapAdjacent3 :: (a -> a -> a -> b) -> [a] -> [b]-mapAdjacent3 f xs0 =-   let xs1 = drop 1 xs0-       xs2 = drop 1 xs1-   in  List.zipWith3 f xs0 xs1 xs2  {- | Cf. 'coefficientsToLinear'

@@ -8,18 +8,16 @@    parabolaBasisDerivativeCenter,    parabolaBasisDerivativeLeft,    )+import Numeric.Interpolation.Private.List (+   mapAdjacentMaybe3,+   mapAdjacentMaybe5,+   )  import Control.Monad (liftM, liftM2) -import qualified Data.List as List import Data.Maybe (catMaybes)  -mapAdjacentMaybe3 :: (Maybe a -> a -> Maybe a -> b) -> [a] -> [b]-mapAdjacentMaybe3 f xs =-   let jxs = map Just xs-   in  zipWith3 f (Nothing : jxs) xs (drop 1 jxs ++ [Nothing])- generic :: ny -> ny -> [x] -> [Nodes.T x ny] generic nz ny =    mapAdjacentMaybe3@@ -34,23 +32,13 @@  hermite1 :: (Num b) => [a] -> [Nodes.T a (b, b)] hermite1 xs =-   generic (0,0) (1,0) xs-   ++-   generic (0,0) (0,1) xs-+   concat $+   zipWith (\f df -> [f,df])+      (generic (0,0) (1,0) xs)+      (generic (0,0) (0,1) xs)   -mapAdjacentMaybe5 ::-   (Maybe a -> Maybe a -> a -> Maybe a -> Maybe a -> b) ->-   [a] -> [b]-mapAdjacentMaybe5 f xs =-   let jxs = map Just xs-       lxs1 = Nothing : jxs-       lxs2 = Nothing : lxs1-       rxs1 = drop 1 $ jxs ++ repeat Nothing-       rxs2 = drop 1 $ rxs1-   in  List.zipWith5 f lxs2 lxs1 xs rxs1 rxs2  cubicAutoGeneric ::    (Num b) =>

@@ -8,11 +8,37 @@  import Data.Tuple.HT (mapFst) +import Data.Traversable (Traversable, traverse)+import Data.Foldable (Foldable, foldMap)+import Data.Monoid (mempty, (<>))++import Control.Applicative (liftA3, pure)+ import Prelude hiding (lookup)   data T x y = Interval | Node (x, y) (T x y) (T x y)    deriving (Eq, Ord, Show)++instance Functor (T x) where+   fmap f =+      let go Interval = Interval+          go (Node (x,y) l r) = Node (x, f y) (go l) (go r)+      in  go++instance Foldable (T x) where+   foldMap f =+      let go Interval = mempty+          go (Node (_x,y) l r) = go l <> f y <> go r+      in  go++instance Traversable (T x) where+   traverse f =+      let go Interval = pure Interval+          go (Node (x,y) l0 r0) =+             liftA3 (\l m r -> Node (x,m) l r) (go l0) (f y) (go r0)+      in  go+  {- | list must be sorted with respect to first element

@@ -0,0 +1,113 @@+module Numeric.Interpolation.Sample (+   T,+   linear,+   hermite1,+   cubicLinear,+   cubicParabola,+   ) where++import qualified Numeric.Interpolation.NodeList as Nodes+import qualified Numeric.Interpolation.Piece as Piece+import Numeric.Interpolation.Private.List (mapAdjacentMaybe3, )+import Numeric.Interpolation.Private.Basis (+   parabolaBasisDerivativeRight,+   parabolaBasisDerivativeCenter,+   parabolaBasisDerivativeLeft,+   )+++type T x y = [x] -> x -> [(Int, y)]++linear :: (Fractional a, Ord a) => T a a+linear nodeXs =+   let nodes = Nodes.fromList $ zip nodeXs [0..]+   in  \x ->+          case Nodes.lookup nodes x of+             (Just (l,nl), Just (r,nr)) ->+                [(nl, Piece.linear (l,1) (r,0) x),+                 (nr, Piece.linear (l,0) (r,1) x)]+             (Just (_l,nl), Nothing) -> [(nl, 1)]+             (Nothing, Just (_r,nr)) -> [(nr, 1)]+             (Nothing, Nothing) -> []++hermite1 :: (Fractional a, Ord a) => T a a+hermite1 nodeXs =+   let nodes = Nodes.fromList $ zip nodeXs [0..]+   in  \x ->+          case Nodes.lookup nodes x of+             (Just (l,nl), Just (r,nr)) ->+                [(2*nl+0, Piece.hermite1 (l,(1,0)) (r,(0,0)) x),+                 (2*nl+1, Piece.hermite1 (l,(0,1)) (r,(0,0)) x),+                 (2*nr+0, Piece.hermite1 (l,(0,0)) (r,(1,0)) x),+                 (2*nr+1, Piece.hermite1 (l,(0,0)) (r,(0,1)) x)]+             (Just (_l,nl), Nothing) -> [(2*nl, 1)]+             (Nothing, Just (_r,nr)) -> [(2*nr, 1)]+             (Nothing, Nothing) -> []++cubicLinear :: (Fractional a, Ord a) => T a a+cubicLinear nodeXs =+   let nodes =+          Nodes.fromList $ zip nodeXs $ zip [0..] $+          mapAdjacentMaybe3 (\l _ r -> (l,r)) nodeXs+   in  \x ->+          case Nodes.lookup nodes x of+             (Nothing, Nothing) -> []+             (Just (_l,(nl,_)), Nothing) -> [(nl-1, 1)]+             (Nothing, Just (_r,(nr,_))) -> [(nr+1, 1)]+             (Just (l,(nl,(mll,_))), Just (r,(nr,(_,mrr)))) ->+                let interL ll =+                       (nl-1, Piece.hermite1 (l,(0,recip(ll-r))) (r,(0,0)) x)+                    interR rr =+                       (nr+1, Piece.hermite1 (l,(0,0)) (r,(0,recip(rr-l))) x)+                in  case (mll,mrr) of+                       (Just ll, Just rr) ->+                          interL ll :+                          (nl, Piece.hermite1 (l,(1,0)) (r,(0,recip(l-rr))) x) :+                          (nr, Piece.hermite1 (l,(0,recip(r-ll))) (r,(1,0)) x) :+                          interR rr :+                          []+                       (Just ll, Nothing) -> interL ll : [(nl, 1)]+                       (Nothing, Just rr) -> interR rr : [(nr, 1)]+                       (Nothing, Nothing) -> []++cubicParabola :: (Fractional a, Ord a) => T a a+cubicParabola nodeXs =+   let nodes =+          Nodes.fromList $ zip nodeXs $ zip [0..] $+          mapAdjacentMaybe3 (\l _ r -> (l,r)) nodeXs+   in  \x ->+          case Nodes.lookup nodes x of+             (Nothing, Nothing) -> []+             (Just (_l,(nl,_)), Nothing) -> [(nl-1, 1)]+             (Nothing, Just (_r,(nr,_))) -> [(nr+1, 1)]+             (Just (l,(nl,(mll,_))), Just (r,(nr,(_,mrr)))) ->+                let interL ll =+                       (nl-1,+                        Piece.hermite1+                           (l,(0, parabolaBasisDerivativeLeft ll l r))+                           (r,(0, 0))+                           x)+                    interR rr =+                       (nr+1,+                        Piece.hermite1+                           (l,(0, 0))+                           (r,(0, parabolaBasisDerivativeRight l r rr))+                           x)+                in  case (mll,mrr) of+                       (Just ll, Just rr) ->+                          interL ll :+                          (nl,+                           Piece.hermite1+                              (l, (1, parabolaBasisDerivativeCenter ll l r))+                              (r, (0, parabolaBasisDerivativeLeft l r rr))+                              x) :+                          (nr,+                           Piece.hermite1+                              (l, (0, parabolaBasisDerivativeRight ll l r))+                              (r, (1, parabolaBasisDerivativeCenter l r rr))+                              x) :+                          interR rr :+                          []+                       (Just ll, Nothing) -> interL ll : [(nl, 1)]+                       (Nothing, Just rr) -> interR rr : [(nr, 1)]+                       (Nothing, Nothing) -> []

@@ -1,7 +1,7 @@ module Numeric.Interpolation.Type (    T(..),    linear,-   cubic,+   hermite1,    cubicLinear,    cubicParabola,    ) where@@ -9,6 +9,7 @@ import qualified Numeric.Interpolation.NodeList as Nodes import qualified Numeric.Interpolation.Piece as Piece import qualified Numeric.Interpolation.Basis as Basis+import qualified Numeric.Interpolation.Sample as Sample import Numeric.Interpolation.Private.Basis (hermite1Split)  @@ -16,24 +17,31 @@    Cons {       ssvFromNodes :: [x] -> [y] -> String,       interpolatePiece :: Piece.T x y ny,+      basisOverlap :: Int+         {- ^+         maximum difference of indices of basis functions that overlap plus one+         -},       basisFunctions :: [x] -> [Nodes.T x ny],+      sampleBasisFunctions :: [x] -> x -> [(Int, y)],       coefficientsToInterpolator :: [x] -> [y] -> Nodes.T x ny,       valueFromNode :: ny -> y    } -linear :: T Double Double Double+linear :: (Fractional a, Ord a, Show a) => T a a a linear =    Cons {       ssvFromNodes =          \xs ys -> unlines $ zipWith (\x y -> show x ++ " " ++ show y) xs ys,       interpolatePiece = Piece.linear,+      basisOverlap = 2,       basisFunctions = Basis.linear,+      sampleBasisFunctions = Sample.linear,       coefficientsToInterpolator = Basis.coefficientsToLinear,       valueFromNode = id    } -cubic :: T Double Double (Double, Double)-cubic =+hermite1 :: (Fractional a, Ord a, Show a) => T a a (a, a)+hermite1 =    Cons {       ssvFromNodes =          \xs ys ->@@ -41,41 +49,49 @@             zipWith (\x (y,dy) -> show x ++ " " ++ show y ++ " " ++ show dy) xs $             hermite1Split xs ys,       interpolatePiece = Piece.hermite1,+      basisOverlap = 4,       basisFunctions = Basis.hermite1,+      sampleBasisFunctions = Sample.hermite1,       coefficientsToInterpolator = Basis.coefficientsToHermite1,       valueFromNode = fst    } -cubicLinear :: T Double Double (Double, Double)+cubicLinear :: (Fractional a, Ord a, Show a) => T a a (a, a) cubicLinear =    Cons {       ssvFromNodes =          \xs ys -> unlines $ zipWith (\x y -> show x ++ " " ++ show y) xs ys,       interpolatePiece = Piece.hermite1,+      basisOverlap = 4,       basisFunctions = Basis.cubicLinear,+      sampleBasisFunctions = Sample.cubicLinear,       coefficientsToInterpolator = Basis.coefficientsToCubicLinear,       valueFromNode = fst    } -cubicParabola :: T Double Double (Double, Double)+cubicParabola :: (Fractional a, Ord a, Show a) => T a a (a, a) cubicParabola =    Cons {       ssvFromNodes =          \xs ys -> unlines $ zipWith (\x y -> show x ++ " " ++ show y) xs ys,       interpolatePiece = Piece.hermite1,+      basisOverlap = 4,       basisFunctions = Basis.cubicParabola,+      sampleBasisFunctions = Sample.cubicParabola,       coefficientsToInterpolator = Basis.coefficientsToCubicParabola,       valueFromNode = fst    }  -_cubicMean :: T Double Double (Double, Double)+_cubicMean :: (Fractional a, Ord a, Show a) => T a a (a, a) _cubicMean =    Cons {       ssvFromNodes =          \xs ys -> unlines $ zipWith (\x y -> show x ++ " " ++ show y) xs ys,       interpolatePiece = Piece.hermite1,+      basisOverlap = 4,       basisFunctions = Basis.cubicParabola, -- Basis.cubicMean,+      sampleBasisFunctions = Sample.cubicParabola, -- Sample.cubicMean,       coefficientsToInterpolator = Basis.coefficientsToCubicParabola, -- not correct       valueFromNode = fst    }

@@ -1,6 +1,8 @@ module Main where  import qualified Test.Piece as Piece+import qualified Test.Sample as Sample+import qualified Test.Overlap as Overlap   run :: String -> [(String, IO ())] -> IO ()@@ -8,5 +10,7 @@    mapM_ (\(msg,act) -> putStr (prefix ++ '.' : msg ++ ": ") >> act)  main :: IO ()-main =+main = do    run "Piece" Piece.tests+   run "Sample" Sample.tests+   run "Overlap" Overlap.tests

@@ -0,0 +1,25 @@+module Test.Overlap where++import qualified Numeric.Interpolation.Type as Type++import Test.QuickCheck (quickCheck, )++++test :: Type.T Double y ny -> IO ()+test typ =+   quickCheck $ \xs xi ->+      let samples = map fst $ Type.sampleBasisFunctions typ xs xi+          {- not total:+          maximum samples - minimum samples < Type.basisOverlap typ+          -}+      in  all (< minimum samples + Type.basisOverlap typ) samples+++tests :: [(String, IO ())]+tests =+   ("linear", test Type.linear) :+   ("hermite1", test Type.hermite1) :+   ("cubicLinear", test Type.cubicLinear) :+   ("cubicParabola", test Type.cubicParabola) :+   []

@@ -0,0 +1,62 @@+module Test.Sample where++import qualified Numeric.Interpolation.Type as Type+import qualified Numeric.Interpolation.Piecewise as Piecewise++import qualified Data.Set as Set+import Data.Array (accumArray, listArray, )+import Data.List.HT (lengthAtLeast, )++import Test.QuickCheck (Property, quickCheck, (==>), )++++withSortedRatios ::+   ([Rational] -> Rational -> a) ->+   ([Integer] -> Integer -> a)+withSortedRatios f nodeXs x =+   f (map fromInteger $ Set.toAscList $ Set.fromList nodeXs) (fromInteger x)++checkEq ::+   (Ord x, Eq y, Num y) =>+   Type.T x y ny -> [x] -> x -> Bool+checkEq typ nodeXs x =+   let ys =+          map+             (flip (Piecewise.interpolateConstantExt typ) x)+             (Type.basisFunctions typ nodeXs)+       bounds = (0, length ys - 1)+   in  listArray bounds ys+       ==+       accumArray (flip const) 0 bounds+          (Type.sampleBasisFunctions typ nodeXs x)+++linear :: [Integer] -> Integer -> Bool+linear = withSortedRatios $ checkEq Type.linear++hermite1 :: [Integer] -> Integer -> Bool+hermite1 = withSortedRatios $ checkEq Type.hermite1+++derivativeFree ::+   Type.T Rational Rational ny ->+   [Integer] -> Integer -> Property+derivativeFree typ =+   withSortedRatios $ \nodeXs x ->+      lengthAtLeast 4 nodeXs ==> checkEq typ nodeXs x++cubicLinear :: [Integer] -> Integer -> Property+cubicLinear = derivativeFree Type.cubicLinear++cubicParabola :: [Integer] -> Integer -> Property+cubicParabola = derivativeFree Type.cubicParabola+++tests :: [(String, IO ())]+tests =+   ("linear", quickCheck linear) :+   ("hermite1", quickCheck hermite1) :+   ("cubicLinear", quickCheck cubicLinear) :+   ("cubicParabola", quickCheck cubicParabola) :+   []