Improve scaled testing

src/scaling/scaling.jl

@@ -5,26 +5,38 @@ end
 ScaledInterpolation{T,ITPT,IT,GT,RT}(::Type{T}, N, itp::ITPT, ::Type{IT}, ::Type{GT}, ranges::RT) =
     ScaledInterpolation{T,N,ITPT,IT,GT,RT}(itp, ranges)
 function scale{T,N,IT,GT}(itp::AbstractInterpolation{T,N,IT,GT}, ranges::Range...)
-    length(ranges) == N || error(string("Must scale $N-dimensional interpolation object with exactly $N ranges (you used $(length(ranges)))"))
+    length(ranges) == N || error("Must scale $N-dimensional interpolation object with exactly $N ranges (you used $(length(ranges)))")
+    for d in 1:N
+        iextract(IT,d) != NoInterp || ranges[d] == 1:size(itp,d) || error("NoInterp dimension $d must be scaled with unit range 1:$(size(itp,d))")
+    end
+
     ScaledInterpolation(T,N,itp,IT,GT,ranges)
 end
 
-@generated function getindex{T,N}(sitp::ScaledInterpolation{T,N}, xs...)
-    length(xs) == N || error(string("Must index into ScaledInterpolation{", N, "} with exactly ", N, " indices (you used ", length(xs), ")"))
-    :(getindex(sitp.itp, coordlookup(sitp.ranges, tuple(xs...))...))
+@generated function getindex{T,N,ITPT,IT<:DimSpec}(sitp::ScaledInterpolation{T,N,ITPT,IT}, xs...)
+    length(xs) == N || error("Must index into $N-dimensional scaled interpolation object with exactly $N indices (you used $(length(xs))")
+    if IT <: Tuple
+        # handle dimspecs
+        interp_indices = map(it -> IT.parameters[it] != NoInterp, 1:length(IT.parameters))
+        :(getindex(sitp.itp, map(coordlookup, $interp_indices, sitp.ranges, xs)...))
+    else
+        :(getindex(sitp.itp, coordlookup(sitp.ranges, tuple(xs...))...))
+    end
 end
 
 size(sitp::ScaledInterpolation, d) = size(sitp.itp, d)
 
 coordlookup(r::FloatRange, x) = (r.divisor * x - r.start) ./ r.step + one(eltype(r))
 coordlookup(r::StepRange, x) = (x - r.start) ./ r.step + one(eltype(r))
 coordlookup(r::UnitRange, x) = x - r.start + one(eltype(r))
+coordlookup(i::Bool, r::Range, x) = i ? coordlookup(r, x) : x
 coordlookup{N}(r::NTuple{N}, x::NTuple{N}) = map(coordlookup,r,x)
+coordlookup{N}(i::NTuple{N}, r::NTuple{N}, x::NTuple{N}) = map(coordlookup, i, r, x)
 
 gradient{T,N}(sitp::ScaledInterpolation{T,N}, xs...) = gradient!(Array(T,N), sitp, xs...)
-@generated function gradient!{T,N}(g, sitp::ScaledInterpolation{T,N}, xs...)
+@generated function gradient!{T,N,ITPT,IT}(g, sitp::ScaledInterpolation{T,N,ITPT,IT}, xs...)
     ndims(g) == 1 || error(string("g must be a vector (but had ", ndims(g), " dimensions)"))
-    length(xs) == N || error(string("Must index into ScaledInterpolation{", N, "} with exactly ", N, " indices (you used ", length(xs), ")"))
+    length(xs) == count_interp_dims(IT, N) || error("Must index into ScaledInterpolation{$N} with exactly $N indices (you used $(length(xs))")
 
     quote
         length(g) == N || error(string("g must be a vector of length ", N, " (was ", length(g), ")"))

test/scaling/dimspecs.jl

@@ -0,0 +1,26 @@
+module ScalingDimspecTests
+
+using Interpolations, DualNumbers, Base.Test
+
+xs = -pi:(2pi/10):pi-2pi/10
+ys = -2:.1:2
+f(x,y) = sin(x) * y^2
+
+itp = interpolate(Float64[f(x,y) for x in xs, y in ys], Tuple{BSpline(Quadratic(Periodic)), BSpline(Linear)}, OnGrid)
+sitp = scale(itp, xs, ys)
+
+# Don't test too near the edges until issue #64 is resolved
+for (ix,x0) in enumerate(xs[5:end-5]), (iy,y0) in enumerate(ys[2:end-1])
+    x, y = x0 + 2pi/20, y0 + .05
+    @test_approx_eq sitp[x0, y0] f(x0,y0)
+    @test_approx_eq_eps sitp[x0, y0] f(x0,y0) 0.05
+
+    g = gradient(sitp, x, y)
+    fx = epsilon(f(dual(x,1), y))
+    fy = (f(x, ys[iy+2]) - f(x, ys[iy+1])) / (ys[iy+2] - ys[iy+1])
+
+    @test_approx_eq_eps g[1] fx 0.15
+    @test_approx_eq_eps g[2] fy 0.05 # gradients for linear interpolation is "easy"
+end
+
+end

test/scaling/nointerp.jl

@@ -0,0 +1,34 @@
+module ScalingNoInterpTests
+
+using Interpolations, Base.Test
+
+xs = -pi:2pi/10:pi
+f1(x) = sin(x)
+f2(x) = cos(x)
+f3(x) = sin(x) .* cos(x)
+f(x,y) = y == 1 ? f1(x) : (y == 2 ? f2(x) : (y == 3 ? f3(x) : error("invalid value for y (must be 1, 2 or 3, you used $y)")))
+ys = 1:3
+
+A = hcat(f1(xs), f2(xs), f3(xs))
+
+itp = interpolate(A, Tuple{BSpline(Quadratic(Periodic)), NoInterp}, OnGrid)
+sitp = scale(itp, xs, ys)
+
+for (ix,x0) in enumerate(xs[1:end-1]), y0 in ys
+    x,y = x0, y0
+    @test_approx_eq_eps sitp[x,y] f(x,y) .05
+end
+
+# Test error messages for incorrect initialization
+function message_is(message)
+	r -> r.err.msg == message || error("Incorrect error message: expected '$message' but was '$(r.msg)'")
+end
+Test.with_handler(message_is("Must scale 2-dimensional interpolation object with exactly 2 ranges (you used 1)")) do
+	@test scale(itp, xs)
+end
+Test.with_handler(message_is("NoInterp dimension 2 must be scaled with unit range 1:3", handler)) do
+	@test scale(itp, xs, -1:1)
+end
+
+end
+

test/scaling/runtests.jl

@@ -1,71 +1,4 @@
-module ScalingTests
-
-using Interpolations
-using Base.Test
-
-# Model linear interpolation of y = -3 + .5x by interpolating y=x
-# and then scaling to the new x range
-
-itp = interpolate(1:1.0:10, BSpline(Linear), OnGrid)
-
-sitp = scale(itp, -3:.5:1.5)
-
-for (x,y) in zip(-3:.05:1.5, 1:.1:10)
-    @test_approx_eq sitp[x] y
-end
-
-# Verify that it works in >1D, with different types of ranges
-
-gauss(phi, mu, sigma) = exp(-(phi-mu)^2 / (2sigma)^2)
-f(x,y) = gauss(x, 0.5, 4) * gauss(y, -.5, 2)
-
-xs = -5:.5:5
-ys = -4:.2:4
-zs = Float64[f(x,y) for x in xs, y in ys]
-
-itp2 = interpolate(zs, BSpline(Quadratic(Flat)), OnGrid)
-sitp2 = scale(itp2, xs, ys)
-
-for x in xs, y in ys
-    @test_approx_eq f(x,y) sitp2[x,y]
-end
-
-# Test gradients of scaled grids
-xs = -pi:.1:pi
-ys = sin(xs)
-itp = interpolate(ys, BSpline(Linear), OnGrid)
-sitp = scale(itp, xs)
-
-for x in -pi:.1:pi
-	g = gradient(sitp, x)[1]
-	@test_approx_eq_eps cos(x) g .05
-end
-
-end
-
-module ScalingDimspecTests
-
-using Interpolations, DualNumbers, Base.Test
-
-xs = -pi:(2pi/10):pi
-ys = -2:.1:2
-f(x,y) = sin(x) * y^2
-
-itp = interpolate(Float64[f(x,y) for x in xs, y in ys], Tuple{BSpline(Quadratic(Periodic)), BSpline(Linear)}, OnGrid)
-sitp = scale(itp, xs, ys)
-
-# Don't test too near the edges until issue #64 is resolved
-for (ix,x0) in enumerate(xs[5:end-5]), (iy,y0) in enumerate(ys[2:end-1])
-    x, y = x0 + 2pi/20, y0 + .05
-    @test_approx_eq sitp[x0, y0] f(x0,y0)
-    @test_approx_eq_eps sitp[x0, y0] f(x0,y0) 0.05
-
-    g = gradient(sitp, x, y)
-    fx = epsilon(f(dual(x,1), y))
-    fy = (f(x, ys[iy+2]) - f(x, ys[iy+1])) / (ys[iy+2] - ys[iy+1])
-
-    @test_approx_eq_eps g[1] fx 0.15
-    @test_approx_eq_eps g[2] fy 0.05 # gradients for linear interpolation is "easy"
-end
-
-end
+include("scaling.jl")
+include("dimspecs.jl")
+include("nointerp.jl")
+include("withextrap.jl")

test/scaling/scaling.jl

@@ -0,0 +1,44 @@
+module ScalingTests
+
+using Interpolations
+using Base.Test
+
+# Model linear interpolation of y = -3 + .5x by interpolating y=x
+# and then scaling to the new x range
+
+itp = interpolate(1:1.0:10, BSpline(Linear), OnGrid)
+
+sitp = scale(itp, -3:.5:1.5)
+
+for (x,y) in zip(-3:.05:1.5, 1:.1:10)
+    @test_approx_eq sitp[x] y
+end
+
+# Verify that it works in >1D, with different types of ranges
+
+gauss(phi, mu, sigma) = exp(-(phi-mu)^2 / (2sigma)^2)
+f(x,y) = gauss(x, 0.5, 4) * gauss(y, -.5, 2)
+
+xs = -5:.5:5
+ys = -4:.2:4
+zs = Float64[f(x,y) for x in xs, y in ys]
+
+itp2 = interpolate(zs, BSpline(Quadratic(Flat)), OnGrid)
+sitp2 = scale(itp2, xs, ys)
+
+for x in xs, y in ys
+    @test_approx_eq f(x,y) sitp2[x,y]
+end
+
+# Test gradients of scaled grids
+xs = -pi:.1:pi
+ys = sin(xs)
+itp = interpolate(ys, BSpline(Linear), OnGrid)
+sitp = scale(itp, xs)
+
+for x in -pi:.1:pi
+    g = gradient(sitp, x)[1]
+    @test_approx_eq_eps cos(x) g .05
+end
+
+end

test/scaling/withextrap.jl

@@ -0,0 +1,7 @@
+
+module ScalingWithExtrapTests
+
+using Interpolations
+
+
+end