Skip to content

Instantly share code, notes, and snippets.

@leftaroundabout
Last active September 27, 2020 19:42
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save leftaroundabout/4fd6ef8642029607e1b222783b9d1c1e to your computer and use it in GitHub Desktop.
Save leftaroundabout/4fd6ef8642029607e1b222783b9d1c1e to your computer and use it in GitHub Desktop.
{-# LANGUAGE OverloadedLists #-}
import qualified Data.Vector.Unboxed as VU
import Data.Vector.Unboxed (Vector, (!))
import Data.MemoTrie (memo2)
import Graphics.Dynamic.Plot.R2
import Text.Printf
newtype Knots = Knots {getIncreasingKnotsSeq :: Vector Double}
baseFuncs :: Knots -- ^ \(\{u_i\}_i\)
-> Int -- ^ Index \(i\) at which to evaluate
-> Int -- ^ Spline degree \(p\)
-> Double -- ^ Position \(u\) at which to evaluate
-> Double
baseFuncs (Knots us) i₀ p₀ u = VU.unsafeHead $ gol i₀ p₀
where gol i 0 = VU.generate (p₀+1) $ \j ->
if u >= us!(i+j)
&& (i+j>=VU.length us || u < us!(i+j+1))
then 1 else 0
gol i p = case gol i (p-1) of
res' -> VU.izipWith
(\j l r -> let i' = i+j
in (u - us!i')/(us!(i'+p) - us!i') * l
+ (us!(i'+p+1) - u)/(us!(i'+p+1) - us!(i'+1)) * r)
res' (VU.unsafeTail res')
main :: IO ()
main = do
plotWindow
[ plotLatest
$ plotDelay 4 mempty :
[ legendName (printf "𝑖=%i, 𝑝=%i" i p)
. plotDelay 1
. continFnPlot $ baseFuncs (Knots us) i p
| p <- [0..4]
, i <- [0..n-p-2]
]
, dynamicAxes
]
return ()
where us = [0, 0.2, 0.3, 0.6, 0.9, 1]
n = VU.length us
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment