Merge pull request #74 from tlycken/instances

RFC: Use instances instead of types when specifying options

README.md

@@ -99,25 +99,25 @@ B-splines of quadratic or higher degree require solving an equation system to ob
 Some examples:
 ```jl
 # Nearest-neighbor interpolation
-itp = interpolate(a, BSpline{Constant}, OnCell)
+itp = interpolate(a, BSpline(Constant()), OnCell())
 v = itp[5.4]   # returns a[5]
 
 # (Multi)linear interpolation
-itp = interpolate(A, BSpline{Linear}, OnGrid)
+itp = interpolate(A, BSpline(Linear()), OnGrid())
 v = itp[3.2, 4.1]  # returns 0.9*(0.8*A[3,4]+0.2*A[4,4]) + 0.1*(0.8*A[3,5]+0.2*A[4,5])
 
 # Quadratic interpolation with reflecting boundary conditions
 # Quadratic is the lowest order that has continuous gradient
-itp = interpolate(A, BSpline{Quadratic{Reflect}}, OnCell)
+itp = interpolate(A, BSpline(Quadratic(Reflect())), OnCell())
 
 # Linear interpolation in the first dimension, and no interpolation (just lookup) in the second
-itp = interpolate(A, Tuple{BSpline{Linear}, BSpline{NoInterp}}, OnGrid)
+itp = interpolate(A, (BSpline(Linear()), NoInterp()), OnGrid())
 v = itp[3.65, 5]  # returns  0.35*A[3,5] + 0.65*A[4,5]
 ```
 There are more options available, for example:
 ```jl
 # In-place interpolation
-itp = interpolate!(A, BSpline{Quadratic{InPlace}}, OnCell)
+itp = interpolate!(A, BSpline(Quadratic(InPlace())), OnCell())
 ```
 which destroys the input `A` but also does not need to allocate as much memory.
 
@@ -133,7 +133,7 @@ In 1D
 A = rand(20)
 A_x = collect(1.:40.)
 knots = (a,)
-itp = interpolate(knots,A, Gridded{Linear})
+itp = interpolate(knots,A, Gridded(Linear()))
 itp[2.0] 
 ```
 
@@ -148,22 +148,22 @@ For example:
 ```jl
 A = rand(8,20)
 knots = ([x^2 for x = 1:8], [0.2y for y = 1:20])
-itp = interpolate(knots, A, Gridded{Linear})
+itp = interpolate(knots, A, Gridded(Linear()))
 itp[4,1.2]  # approximately A[2,6]
 ```
 One may also mix modes, by specifying a mode vector in the form of an explicit tuple:
 ```jl
-itp = interpolate(knots, A, Tuple{Gridded{Linear},Gridded{Constant}})
+itp = interpolate(knots, A, (Gridded(Linear()),Gridded(Constant())))
 ```
 
 
 Presently there are only three modes for gridded:
 ```jl
-Gridded{Linear}
+Gridded(Linear())
 ```
 whereby a linear interpolation is applied between knots,
 ```jl
-Gridded{Constant}
+Gridded(Constant())
 ```
 whereby nearest neighbor interpolation is used on the applied axis,
 ```jl
@@ -174,11 +174,7 @@ coordinates results in the throwing of an error.
 
 ## Extrapolation
 
-The call to `extrapolate` defines what happens if you try to index into the interpolation object with coordinates outside of `[1, size(data,d)]` in any dimension `d`. The implemented boundary conditions are `Throw` and `Flat`, with more options planned.
-
-## More examples
-
-There's an [IJulia notebook](http://nbviewer.ipython.org/github/tlycken/Interpolations.jl/blob/master/doc/Interpolations.jl.ipynb) that shows off some of the current functionality, and outlines where we're headed. I try to keep it up to date when I make any significant improvements and/or breaking changes, but if it's not, do file a PR.
+The call to `extrapolate` defines what happens if you try to index into the interpolation object with coordinates outside of `[1, size(data,d)]` in any dimension `d`. The implemented boundary conditions are `Throw`, `Flat`, `Linear`, `Periodic` and `Reflect`, with more options planned. `Periodic` and `Reflect` require that there is a method of `Base.mod` that can handle the indices used.
 
 ## Performance shootout
 

doc/devdocs.md

@@ -18,7 +18,7 @@ First let's create an interpolation object:
     0.134746
     0.430876
 
-    julia> yitp = interpolate(A, BSpline{Linear}, OnGrid)
+    julia> yitp = interpolate(A, BSpline(Linear()), OnGrid())
     5-element Interpolations.BSplineInterpolation{Float64,1,Float64,Interpolations.BSpline{Interpolations.Linear},Interpolations.OnGrid}:
     0.74838
     0.995383

src/Interpolations.jl

@@ -32,9 +32,10 @@ using WoodburyMatrices, Ratios, AxisAlgorithms
 
 import Base: convert, size, getindex, gradient, scale, promote_rule
 
-abstract InterpolationType
+abstract Flag
+abstract InterpolationType <: Flag
 immutable NoInterp <: InterpolationType end
-abstract GridType
+abstract GridType <: Flag
 immutable OnGrid <: GridType end
 immutable OnCell <: GridType end
 
@@ -44,16 +45,16 @@ abstract AbstractInterpolation{T,N,IT<:DimSpec{InterpolationType},GT<:DimSpec{Gr
 abstract AbstractInterpolationWrapper{T,N,ITPT,IT,GT} <: AbstractInterpolation{T,N,IT,GT}
 abstract AbstractExtrapolation{T,N,ITPT,IT,GT} <: AbstractInterpolationWrapper{T,N,ITPT,IT,GT}
 
-abstract BoundaryCondition
-immutable Flat <: BoundaryCondition end
-immutable Line <: BoundaryCondition end
-immutable Free <: BoundaryCondition end
-immutable Periodic <: BoundaryCondition end
-immutable Reflect <: BoundaryCondition end
-immutable InPlace <: BoundaryCondition end
+immutable Throw <: Flag end
+immutable Flat <: Flag end
+immutable Line <: Flag end
+immutable Free <: Flag end
+immutable Periodic <: Flag end
+immutable Reflect <: Flag end
+immutable InPlace <: Flag end
 # InPlaceQ is exact for an underlying quadratic. This is nice for ground-truth testing
 # of in-place (unpadded) interpolation.
-immutable InPlaceQ <: BoundaryCondition end
+immutable InPlaceQ <: Flag end
 typealias Natural Line
 
 # TODO: size might have to be faster?

src/b-splines/b-splines.jl

@@ -6,17 +6,17 @@ export
     Linear,
     Quadratic
 
-abstract Degree{N}
+abstract Degree{N} <: Flag
 
 immutable BSpline{D<:Degree} <: InterpolationType end
-BSpline{D<:Degree}(::Type{D}) = BSpline{D}
+BSpline{D<:Degree}(::D) = BSpline{D}()
 
 bsplinetype{D<:Degree}(::Type{BSpline{D}}) = D
 
 immutable BSplineInterpolation{T,N,TCoefs,IT<:DimSpec{BSpline},GT<:DimSpec{GridType},pad} <: AbstractInterpolation{T,N,IT,GT}
     coefs::Array{TCoefs,N}
 end
-function BSplineInterpolation{N,TCoefs,TWeights<:Real,IT<:DimSpec{BSpline},GT<:DimSpec{GridType},pad}(::Type{TWeights}, A::AbstractArray{TCoefs,N}, ::Type{IT}, ::Type{GT}, ::Val{pad})
+function BSplineInterpolation{N,TCoefs,TWeights<:Real,IT<:DimSpec{BSpline},GT<:DimSpec{GridType},pad}(::Type{TWeights}, A::AbstractArray{TCoefs,N}, ::IT, ::GT, ::Val{pad})
     isleaftype(IT) || error("The b-spline type must be a leaf type (was $IT)")
     isleaftype(TCoefs) || warn("For performance reasons, consider using an array of a concrete type (eltype(A) == $(eltype(A)))")
 
@@ -47,12 +47,12 @@ function size{T,N,TCoefs,IT,GT,pad}(itp::BSplineInterpolation{T,N,TCoefs,IT,GT,p
     d <= N ? size(itp.coefs, d) - 2*padextract(pad, d) : 1
 end
 
-function interpolate{TWeights,TCoefs,IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(::Type{TWeights}, ::Type{TCoefs}, A, ::Type{IT}, ::Type{GT})
+function interpolate{TWeights,TCoefs,IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(::Type{TWeights}, ::Type{TCoefs}, A, it::IT, gt::GT)
     Apad, Pad = prefilter(TWeights, TCoefs, A, IT, GT)
-    BSplineInterpolation(TWeights, Apad, IT, GT, Pad)
+    BSplineInterpolation(TWeights, Apad, it, gt, Pad)
 end
-function interpolate{IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(A::AbstractArray, ::Type{IT}, ::Type{GT})
-    interpolate(tweight(A), tcoef(A), A, IT, GT)
+function interpolate{IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(A::AbstractArray, it::IT, gt::GT)
+    interpolate(tweight(A), tcoef(A), A, it, gt)
 end
 
 # We can't just return a tuple-of-types due to julia #12500
@@ -66,9 +66,9 @@ tcoef(A::AbstractArray{Float32}) = Float32
 tcoef(A::AbstractArray{Rational{Int}}) = Rational{Int}
 tcoef{T<:Integer}(A::AbstractArray{T}) = typeof(float(zero(T)))
 
-interpolate!{TWeights,IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(::Type{TWeights}, A, ::Type{IT}, ::Type{GT}) = BSplineInterpolation(TWeights, prefilter!(TWeights, A, IT, GT), IT, GT, Val{0}())
-function interpolate!{IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(A::AbstractArray, ::Type{IT}, ::Type{GT})
-    interpolate!(tweight(A), A, IT, GT)
+interpolate!{TWeights,IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(::Type{TWeights}, A, it::IT, gt::GT) = BSplineInterpolation(TWeights, prefilter!(TWeights, A, IT, GT), it, gt, Val{0}())
+function interpolate!{IT<:DimSpec{BSpline},GT<:DimSpec{GridType}}(A::AbstractArray, it::IT, gt::GT)
+    interpolate!(tweight(A), A, it, gt)
 end
 
 include("constant.jl")

src/b-splines/quadratic.jl

@@ -1,5 +1,5 @@
-immutable Quadratic{BC<:BoundaryCondition} <: Degree{2} end
-Quadratic{BC<:BoundaryCondition}(::Type{BC}) = Quadratic{BC}
+immutable Quadratic{BC<:Flag} <: Degree{2} end
+Quadratic{BC<:Flag}(::BC) = Quadratic{BC}()
 
 function define_indices_d{BC}(::Type{BSpline{Quadratic{BC}}}, d, pad)
     symix, symixm, symixp = symbol("ix_",d), symbol("ixm_",d), symbol("ixp_",d)
@@ -73,7 +73,7 @@ function index_gen{Q<:Quadratic,IT<:DimSpec{BSpline}}(::Type{BSpline{Q}}, ::Type
     end
 end
 
-padding{BC<:BoundaryCondition}(::Type{BSpline{Quadratic{BC}}}) = Val{1}()
+padding{BC<:Flag}(::Type{BSpline{Quadratic{BC}}}) = Val{1}()
 padding(::Type{BSpline{Quadratic{Periodic}}}) = Val{0}()
 
 function inner_system_diags{T,Q<:Quadratic}(::Type{T}, n::Int, ::Type{Q})

src/extrapolation/extrapolation.jl

@@ -1,11 +1,15 @@
-export Throw
-
 type Extrapolation{T,N,ITPT,IT,GT,ET} <: AbstractInterpolationWrapper{T,N,ITPT,IT,GT}
     itp::ITPT
 end
-Extrapolation{T,ITPT,IT,GT,ET}(::Type{T}, N, itp::ITPT, ::Type{IT}, ::Type{GT}, ::Type{ET}) =
+Extrapolation{T,ITPT,IT,GT,ET}(::Type{T}, N, itp::ITPT, ::Type{IT}, ::Type{GT}, et::ET) =
     Extrapolation{T,N,ITPT,IT,GT,ET}(itp)
 
+# DimSpec{Flag} is not enough for extrapolation dispatch, since we allow nested tuples
+# However, no tuples should be nested deeper than this; the first level is for different
+# schemes in different dimensions, and the second level is for different schemes in
+# different directions.
+typealias ExtrapDimSpec Union{Flag,Tuple{Vararg{Union{Flag,NTuple{2,Flag}}}}}
+
 """
 `extrapolate(itp, scheme)` adds extrapolation behavior to an interpolation object, according to the provided scheme.
 
@@ -14,15 +18,15 @@ The scheme can take any of these values:
 * `Throw` - throws a BoundsError for out-of-bounds indices
 * `Flat` - for constant extrapolation, taking the closest in-bounds value
 * `Linear` - linear extrapolation (the wrapped interpolation object must support gradient)
-* `Reflect` - reflecting extrapolation
-* `Periodic` - periodic extrapolation
+* `Reflect` - reflecting extrapolation (indices must support `mod`)
+* `Periodic` - periodic extrapolation (indices must support `mod`)
 
-You can also combine schemes in tuples. For example, the scheme `Tuple{Linear, Flat}` will use linear extrapolation in the first dimension, and constant in the second.
+You can also combine schemes in tuples. For example, the scheme `(Linear(), Flat())` will use linear extrapolation in the first dimension, and constant in the second.
 
-Finally, you can specify different extrapolation behavior in different direction. `Tuple{Tuple{Linear,Flat}, Flat}` will extrapolate linearly in the first dimension if the index is too small, but use constant etrapolation if it is too large, and always use constant extrapolation in the second dimension.
+Finally, you can specify different extrapolation behavior in different direction. `((Linear(),Flat()), Flat())` will extrapolate linearly in the first dimension if the index is too small, but use constant etrapolation if it is too large, and always use constant extrapolation in the second dimension.
 """
-extrapolate{T,N,IT,GT,ET}(itp::AbstractInterpolation{T,N,IT,GT}, ::Type{ET}) =
-    Extrapolation(T,N,itp,IT,GT,ET)
+extrapolate{T,N,IT,GT,ET<:ExtrapDimSpec}(itp::AbstractInterpolation{T,N,IT,GT}, et::ET) =
+    Extrapolation(T,N,itp,IT,GT,et)
 
 include("throw.jl")
 include("flat.jl")

src/extrapolation/filled.jl

@@ -1,6 +1,5 @@
 nindexes(N::Int) = N == 1 ? "1 index" : "$N indexes"
 
-
 type FilledExtrapolation{T,N,ITP<:AbstractInterpolation,IT,GT,FT} <: AbstractExtrapolation{T,N,ITP,IT,GT}
     itp::ITP
     fillvalue::FT

src/extrapolation/throw.jl

@@ -1,5 +1,3 @@
-immutable Throw end
-
 function extrap_prep{N,d}(::Type{Throw}, ::Val{N}, ::Val{d})
     xsym = symbol("xs_", d)
     :(lbound(etp, $d) <= $xsym <= ubound(etp, $d) || throw(BoundsError()))

src/gridded/gridded.jl

@@ -1,7 +1,7 @@
 export Gridded
 
 immutable Gridded{D<:Degree} <: InterpolationType end
-Gridded{D<:Degree}(::Type{D}) = Gridded{D}
+Gridded{D<:Degree}(::D) = Gridded{D}()
 
 griddedtype{D<:Degree}(::Type{Gridded{D}}) = D
 
@@ -13,7 +13,7 @@ immutable GriddedInterpolation{T,N,TCoefs,IT<:DimSpec{Gridded},K<:Tuple{Vararg{V
     knots::K
     coefs::Array{TCoefs,N}
 end
-function GriddedInterpolation{N,TCoefs,TWeights<:Real,IT<:DimSpec{Gridded},pad}(::Type{TWeights}, knots::NTuple{N,GridIndex}, A::AbstractArray{TCoefs,N}, ::Type{IT}, ::Val{pad})
+function GriddedInterpolation{N,TCoefs,TWeights<:Real,IT<:DimSpec{Gridded},pad}(::Type{TWeights}, knots::NTuple{N,GridIndex}, A::AbstractArray{TCoefs,N}, ::IT, ::Val{pad})
     isleaftype(IT) || error("The b-spline type must be a leaf type (was $IT)")
     isleaftype(TCoefs) || warn("For performance reasons, consider using an array of a concrete type (eltype(A) == $(eltype(A)))")
 
@@ -48,16 +48,16 @@ iextract{T<:GridType}(::Type{T}, d) = T
     end
 end
 
-function interpolate{TWeights,TCoefs,Tel,N,IT<:DimSpec{Gridded}}(::Type{TWeights}, ::Type{TCoefs}, knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, ::Type{IT})
-    GriddedInterpolation(TWeights, knots, A, IT, Val{0}())
+function interpolate{TWeights,TCoefs,Tel,N,IT<:DimSpec{Gridded}}(::Type{TWeights}, ::Type{TCoefs}, knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, it::IT)
+    GriddedInterpolation(TWeights, knots, A, it, Val{0}())
 end
-function interpolate{Tel,N,IT<:DimSpec{Gridded}}(knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, ::Type{IT})
-    interpolate(tweight(A), tcoef(A), knots, A, IT)
+function interpolate{Tel,N,IT<:DimSpec{Gridded}}(knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, it::IT)
+    interpolate(tweight(A), tcoef(A), knots, A, it)
 end
 
-interpolate!{TWeights,Tel,N,IT<:DimSpec{Gridded}}(::Type{TWeights}, knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, ::Type{IT}) = GriddedInterpolation(TWeights, knots, A, IT, Val{0}())
-function interpolate!{Tel,N,IT<:DimSpec{Gridded}}(knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, ::Type{IT})
-    interpolate!(tweight(A), tcoef(A), knots, A, IT)
+interpolate!{TWeights,Tel,N,IT<:DimSpec{Gridded}}(::Type{TWeights}, knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, it::IT) = GriddedInterpolation(TWeights, knots, A, it, Val{0}())
+function interpolate!{Tel,N,IT<:DimSpec{Gridded}}(knots::NTuple{N,GridIndex}, A::AbstractArray{Tel,N}, it::IT)
+    interpolate!(tweight(A), tcoef(A), knots, A, it)
 end
 
 include("constant.jl")

test/b-splines/constant.jl

@@ -10,12 +10,12 @@ A2 = rand(Float64, N1, N1) * 100
 A3 = rand(Float64, N1, N1, N1) * 100
 
 for (constructor, copier) in ((interpolate, x->x), (interpolate!, copy))
-    itp1c = @inferred(constructor(copier(A1), BSpline(Constant), OnCell))
-    itp1g = @inferred(constructor(copier(A1), BSpline(Constant), OnGrid))
-    itp2c = @inferred(constructor(copier(A2), BSpline(Constant), OnCell))
-    itp2g = @inferred(constructor(copier(A2), BSpline(Constant), OnGrid))
-    itp3c = @inferred(constructor(copier(A3), BSpline(Constant), OnCell))
-    itp3g = @inferred(constructor(copier(A3), BSpline(Constant), OnGrid))
+    itp1c = @inferred(constructor(copier(A1), BSpline(Constant()), OnCell()))
+    itp1g = @inferred(constructor(copier(A1), BSpline(Constant()), OnGrid()))
+    itp2c = @inferred(constructor(copier(A2), BSpline(Constant()), OnCell()))
+    itp2g = @inferred(constructor(copier(A2), BSpline(Constant()), OnGrid()))
+    itp3c = @inferred(constructor(copier(A3), BSpline(Constant()), OnCell()))
+    itp3g = @inferred(constructor(copier(A3), BSpline(Constant()), OnGrid()))
 
     # Evaluation on provided data points
     # 1D

test/b-splines/linear.jl

@@ -10,7 +10,7 @@ f(x) = g1(x)
 A1 = Float64[f(x) for x in 1:xmax]
 
 for (constructor, copier) in ((interpolate, x->x), (interpolate!, copy))
-    itp1c = @inferred(constructor(copier(A1), BSpline(Linear), OnCell))
+    itp1c = @inferred(constructor(copier(A1), BSpline(Linear()), OnCell()))
 
     # Just interpolation
     for x in 1:.2:xmax
@@ -20,7 +20,7 @@ for (constructor, copier) in ((interpolate, x->x), (interpolate!, copy))
     # Rational element types
     fr(x) = (x^2) // 40 + 2
     A1R = Rational{Int}[fr(x) for x in 1:10]
-    itp1r = @inferred(constructor(copier(A1R), BSpline(Linear), OnGrid))
+    itp1r = @inferred(constructor(copier(A1R), BSpline(Linear()), OnGrid()))
     @test @inferred(size(itp1r)) == size(A1R)
     @test_approx_eq_eps itp1r[23//10] fr(23//10) abs(.1*fr(23//10))
     @test typeof(itp1r[23//10]) == Rational{Int}
@@ -31,7 +31,7 @@ for (constructor, copier) in ((interpolate, x->x), (interpolate!, copy))
     f(x,y) = g1(x)*g2(y)
     ymax = 10
     A2 = Float64[f(x,y) for x in 1:xmax, y in 1:ymax]
-    itp2 = @inferred(constructor(copier(A2), BSpline(Linear), OnGrid))
+    itp2 = @inferred(constructor(copier(A2), BSpline(Linear()), OnGrid()))
     @test @inferred(size(itp2)) == size(A2)
 
     for x in 2.1:.2:xmax-1, y in 1.9:.2:ymax-.9

test/b-splines/mixed.jl

@@ -7,14 +7,14 @@ N = 10
 for (constructor, copier) in ((interpolate, x->x), (interpolate!, copy))
     A2 = rand(Float64, N, N) * 100
     for BC in (Flat,Line,Free,Periodic,Reflect,Natural), GT in (OnGrid, OnCell)
-        itp_a = @inferred(constructor(copier(A2), Tuple{BSpline(Linear), BSpline(Quadratic(BC))}, GT))
-        itp_b = @inferred(constructor(copier(A2), Tuple{BSpline(Quadratic(BC)), BSpline(Linear)}, GT))
+        itp_a = @inferred(constructor(copier(A2), (BSpline(Linear()), BSpline(Quadratic(BC()))), GT()))
+        itp_b = @inferred(constructor(copier(A2), (BSpline(Quadratic(BC())), BSpline(Linear())), GT()))
         @test @inferred(size(itp_a)) == size(A2)
         @test @inferred(size(itp_b)) == size(A2)
 
         for j = 2:N-1, i = 2:N-1
-            @test_approx_eq itp_a[i,j] A2[i,j]
-            @test_approx_eq itp_b[i,j] A2[i,j]
+            @test_approx_eq_eps itp_a[i,j] A2[i,j] sqrt(eps(A2[i,j]))
+            @test_approx_eq_eps itp_b[i,j] A2[i,j] sqrt(eps(A2[i,j]))
         end
 
         for i = 1:10

test/b-splines/multivalued.jl

@@ -23,20 +23,20 @@ Base.promote_rule{T1,T2<:Number}(::Type{MyPair{T1}}, ::Type{T2}) = MyPair{promot
 
 # 1d
 A = reinterpret(MyPair{Float64}, rand(2, 10), (10,))
-itp = interpolate(A, BSpline(Constant), OnGrid)
+itp = interpolate(A, BSpline(Constant()), OnGrid())
 itp[3.2]
-itp = interpolate(A, BSpline(Linear), OnGrid)
+itp = interpolate(A, BSpline(Linear()), OnGrid())
 itp[3.2]
-itp = interpolate(A, BSpline(Quadratic(Flat)), OnGrid)
+itp = interpolate(A, BSpline(Quadratic(Flat())), OnGrid())
 itp[3.2]
 
 # 2d
 A = reinterpret(MyPair{Float64}, rand(2, 10, 5), (10,5))
-itp = interpolate(A, BSpline(Constant), OnGrid)
+itp = interpolate(A, BSpline(Constant()), OnGrid())
 itp[3.2,1.8]
-itp = interpolate(A, BSpline(Linear), OnGrid)
+itp = interpolate(A, BSpline(Linear()), OnGrid())
 itp[3.2,1.8]
-itp = interpolate(A, BSpline(Quadratic(Flat)), OnGrid)
+itp = interpolate(A, BSpline(Quadratic(Flat())), OnGrid())
 itp[3.2,1.8]
 
 end