Implement a few more BC:s + whitespace and test fixes

src/Interpolations.jl

@@ -87,6 +87,7 @@ for IT in (
         Linear{OnCell},
         Quadratic{Flat,OnCell},
         Quadratic{Flat,OnGrid},
+        Quadratic{LinearBC,OnGrid},
         Quadratic{LinearBC,OnCell}
     )
     for EB in (
@@ -137,7 +138,8 @@ end
 for IT in (
         Quadratic{Flat,OnCell},
         Quadratic{Flat,OnGrid},
-        Quadratic{LinearBC,OnCell}
+        Quadratic{LinearBC,OnGrid},
+        Quadratic{LinearBC,OnCell},
     )
     @ngenerate N promote_type_grid(T, x...) function prefilter{T,N}(A::Array{T,N},it::IT)
         ret, pad = similar_with_padding(A,it)

src/extrapolation.jl

@@ -1,27 +1,26 @@
 abstract ExtrapolationBehavior
-type ExtrapError <: ExtrapolationBehavior end
 
+type ExtrapError <: ExtrapolationBehavior end
 function extrap_gen(::OnGrid, ::ExtrapError, N)
     quote
         @nexprs $N d->(1 <= x_d <= size(itp,d) || throw(BoundsError()))
     end
 end
 function extrap_gen(::OnCell, ::ExtrapError, N)
     quote
-        @nexprs $N d->(.5 <= x_d <= size(itp,d)+.5 || throw(BoundsError()))
+        @nexprs $N d->(.5 <= x_d <= size(itp,d) + .5 || throw(BoundsError()))
     end
 end
 
 type ExtrapNaN <: ExtrapolationBehavior end
-
 function extrap_gen(::OnGrid, ::ExtrapNaN, N)
     quote
         @nexprs $N d->(1 <= x_d <= size(itp,d) || return convert(T, NaN))
     end
 end
 function extrap_gen(::OnCell, ::ExtrapNaN, N)
     quote
-        @nexprs $N d->(.5 <= x_d <= size(itp,d)+.5 || return convert(T, NaN))
+        @nexprs $N d->(.5 <= x_d <= size(itp,d) + .5 || return convert(T, NaN))
     end
 end
 
@@ -89,7 +88,7 @@ function extrap_gen(::OnGrid, ::ExtrapLinear, N)
             if x_d < 1
                 fx_d = x_d - convert(typeof(x_d), 1)
 
-                k = itp[1] - itp[2]
+                k = -4*itp.coefs[1]
                 return itp[1] - k * fx_d
             end
             if x_d > size(itp, d)
@@ -102,4 +101,4 @@ function extrap_gen(::OnGrid, ::ExtrapLinear, N)
         end
     end
 end
-#extrap_gen(::OnCell, e::ExtrapLinear, N) = extrap_gen(OnGrid(), e, N)
+extrap_gen(::OnCell, e::ExtrapLinear, N) = extrap_gen(OnGrid(), e, N)

src/linear.jl

@@ -10,10 +10,7 @@ end
 function bc_gen(::Linear{OnCell}, N)
     quote
         @nexprs $N d->begin
-            ix_d = ifloor(x_d)
-            if .5 <= x_d < 1
-                ix_d = 1
-            end
+            ix_d = clamp(ifloor(x_d), 1, size(itp, d))
         end
     end
 end

src/quadratic.jl

@@ -25,9 +25,9 @@ end
 function coefficients(::Quadratic, N)
     quote
         @nexprs $N d->begin
-            cm_d = (fx_d-.5)^2 / 2
+            cm_d = .5 * (fx_d-.5)^2
             c_d = .75 - fx_d^2
-            cp_d = (fx_d+.5)^2 / 2
+            cp_d = .5 * (fx_d+.5)^2
         end
     end
 end

test/linear.jl

@@ -52,7 +52,7 @@ f(x,y) = sin(x/10)*cos(y/6)
 
 xmax, ymax = 30,10
 
-A = [f(x,y) for x in 1:xmax, y in 1:ymax]
+A = Float64[f(x,y) for x in 1:xmax, y in 1:ymax]
 
 itp1 = Interpolation(A, Linear(OnGrid()), ExtrapError())
 

test/quadratic.jl

@@ -10,7 +10,7 @@ A = Float64[f(x) for x in 1:xmax]
 
 ## ExtrapError
 
-itp1 = Interpolation(A, Quadratic(ExtendInner(),OnCell()), ExtrapError())
+itp1 = Interpolation(A, Quadratic(Flat(),OnCell()), ExtrapError())
 
 for x in [3.1:.2:4.3]  
     @test_approx_eq_eps f(x) itp1[x] abs(.1*f(x))
@@ -21,19 +21,22 @@ end
 
 ## ExtrapNaN
 
-itp2 = Interpolation(A, Quadratic(ExtendInner(),OnCell()), ExtrapNaN())
+itp2 = Interpolation(A, Quadratic(Flat(),OnCell()), ExtrapNaN())
 
-for x in [3.1:.2:4.3]  
-    @test_approx_eq_eps f(x) itp1[x] abs(.1*f(x))
+for x in [3.1:.2:4.3]
+    @test_approx_eq_eps f(x) itp2[x] abs(.1*f(x))
 end
 
-xlo, xhi = itp2[.9], itp2[xmax+.2]
-@test isnan(xlo)
-@test isnan(xhi)
+@test isnan(itp2[.4])
+@test !isnan(itp2[.5])
+@test !isnan(itp2[.6])
+@test !isnan(itp2[xmax+.4])
+@test !isnan(itp2[xmax+.5])
+@test isnan(itp2[xmax+.6])
 
 # Flat/ExtrapConstant
 
-itp3 = Interpolation(A, Quadratic(Flat(),OnCell()), ExtrapConstant())
+itp3 = Interpolation(A, Quadratic(Flat(),OnGrid()), ExtrapConstant())
 
 # Check inbounds and extrap values
 
@@ -42,8 +45,8 @@ for x in [3.1:.2:4.3]
 end
 
 xlo, xhi = itp3[.9], itp3[xmax+.2]
-@test xlo == A[1]
-@test xhi == A[end]
+@test_approx_eq xlo A[1]
+@test_approx_eq xhi A[end]
 
 # Check continuity
 xs = [0:.1:length(A)+1]

test/visual.jl

@@ -1,69 +1,35 @@
 module VisualTests
 
+print("loading packages...")
 using Gadfly, Interpolations
-
+println("done!")
 const nx = 10
 xg = 1:nx
-xf = -2nx:.1:3nx
+xf = -.2nx:.1:4.3nx
 
 f1(x) = sin(2pi/nx * (x-1))
 
 y1 = map(f1,xg)
 
-for EB in (
-        ExtrapNaN,
-        ExtrapConstant,
-        ExtrapReflect
+for (IT, EB) in (
+        (Constant{OnCell}, ExtrapConstant),
+        (Constant{OnGrid}, ExtrapConstant),
+        (Linear{OnCell}, ExtrapNaN),
+        (Linear{OnGrid}, ExtrapNaN),
+        (Quadratic{Flat,OnCell}, ExtrapConstant),
+        (Quadratic{Flat,OnGrid}, ExtrapConstant),
+        (Quadratic{LinearBC,OnCell}, ExtrapLinear),
+        (Quadratic{LinearBC,OnGrid}, ExtrapLinear),
+        (Quadratic{LinearBC,OnGrid}, ExtrapReflect),
+        (Quadratic{Flat,OnCell}, ExtrapReflect)
     )
+    itp = Interpolation(y1, IT(), EB())
 
-    for IT in (
-            Constant{OnCell},
-            Constant{OnGrid},
-            Linear{OnCell},
-            Linear{OnGrid},
-            Quadratic{ExtendInner,OnCell},
-            Quadratic{ExtendInner,OnGrid},
-            Quadratic{Flat,OnCell},
-            Quadratic{Flat,OnGrid},
-        )
-
-        itp = Interpolation(y1, IT(), EB())
-
-        display(plot(
-            layer(x=xg,y=y1,Geom.point),
-            layer(x=xf,y=[itp[x] for x in xf],Geom.path),
-            Guide.title("$(typeof(Interpolations.degree(IT()))), $(typeof(Interpolations.gridrepresentation(IT()))), $EB")
-        ))
-    end
+    display(plot(
+        layer(x=xg,y=y1,Geom.point),
+        layer(x=xf,y=[itp[x] for x in xf],Geom.path),
+        Guide.title("$IT, $EB")
+    ))
 end
 
-for IT in (
-        Constant{OnCell},
-        Constant{OnGrid},
-        Linear{OnCell},
-        Linear{OnGrid},
-        Quadratic{ExtendInner,OnCell},
-        Quadratic{ExtendInner,OnGrid},
-        Quadratic{Flat,OnCell},
-        Quadratic{Flat,OnGrid},
-    )
-
-    # Treat ExtrapError specially, since it will throw bounds
-    # errors for any x outside 1:nx
-    itperr = Interpolation(y1, IT(), ExtrapError())
-    if Interpolations.gridrepresentation(IT()) == OnCell()
-        display(plot(
-            layer(x=xg,y=y1,Geom.point),
-            layer(x=.5:.1:nx+.5,y=[itperr[x] for x in .5:.1:nx+.5],Geom.path),
-            Guide.title("$(Interpolations.degree(IT())), $(Interpolations.gridrepresentation(IT())), ExtrapError")
-        ))
-    else
-        display(plot(
-            layer(x=xg,y=y1,Geom.point),
-            layer(x=1:.1:nx,y=[itperr[x] for x in 1:.1:nx],Geom.path),
-            Guide.title("$(Interpolations.degree(IT())), $(Interpolations.gridrepresentation(IT())), ExtrapError")
-        ))
-    end
 end
-
-end