(* kernel.sml
*
* COPYRIGHT (c) 2010 The Diderot Project (http://diderot.cs.uchicago.edu)
* All rights reserved.
*
* QUESTION: should we
*)
structure Kernel : sig
type coefficient = Rational.rat
(* polynomial represented as list of coefficients, where ith element is
* coefficient for x^i.
*)
type polynomial = coefficient list
type kernel = {
support : int, (* number of samples to left/right *)
segs : polynomial list, (* piece-wise polynomial that defines the kernel. Since *)
(* a kernel is symmetric, we only store the positive *)
(* segments, with segs[i] covering the domain [i,i+1) *)
}
val differentiate : polynomial -> polynomial
val evaluate : polynomial * int -> Rational.rat
end = struct
structure R = Rational
type coefficient = R.rat
val zero = R.fromInt 0
val one = R.fromInt 1
(* polynomial represented as list of coefficients, where ith element is
* coefficient for x^i.
*)
type polynomial = coefficient list
type kernel = {
support : int, (* number of samples to left/right *)
segs : polynomial list, (* piece-wise polynomial that defines the *)
(* kernel; there are 2*support segments *)
}
fun differentiate [] = raise Fail "invalid polynomial"
| differentiate [_] = [zero]
| differentiate (_::coeffs) = let
fun lp (_, []) = []
| lp (i, c::r) = R.mul(F.fromInt i, c) :: lp(i+1, r)
in
lp (1, coeffs)
end
(* evaluate a polynomial at an integer coordinate (used to test continuity) *)
fun evaluate (poly, x) = let
val x = F.fromInt x
fun eval (sum, [], xn) = sum
| eval (sum, c::r, xn) = eval(F.add(sum, R.mul(c, xn)), r, R.mul(x, xn))
in
eval (zero, poly, one)
end
(* some standard kernels *)
local
val op / = R./
fun r i = R.fromInt i
in
val tent : kernel = { (* linear interpolation *)
support = 1,
segs = [[r 1, r ~1]]
}
val ctmr : kernel = { (* Catmull-Rom interpolation *)
support = 2,
segs = [
[r 1, r 0, ~5/2, 3/2],
[r 2, r ~4, 5/2, ~1/2]
]
}
val bspl3 : kernel = { (* cubic bspline reconstruction, doesn't interpolate *)
support = 2,
segs = [
[ 2/3, r 0, r ~1, 1/2 ],
[ 4/3, r ~2, r 1, ~1/6 ]'
]
}
val bspl5 : kernel = { (* quintic bspline reconstruction, doesn't interpolate *)
support = 3,
segs = [
[ 11/20, r 0, ~1/2, r 0, 1/4, ~1/12 ],
[ 17/40, 5/8, ~7/4, 5/4, ~3/8, 1/24 ],
[ 81/40, ~27/8, 9/4, ~3/4, 1/8, ~1/120 ]
]
}
end
end