Skip to content
This repository has been archived by the owner on Feb 9, 2020. It is now read-only.

Support unequally spaced grids #1

Open
stevengj opened this issue Mar 13, 2013 · 8 comments
Open

Support unequally spaced grids #1

stevengj opened this issue Mar 13, 2013 · 8 comments

Comments

@stevengj
Copy link

Similar to the interp functions in Matlab, it would be good to specify unequally spaced (Cartesian) grids.

@timholy
Copy link
Owner

timholy commented Mar 13, 2013

Agreed, it would be useful functionality, although I don't personally have any need for it (at least not in the immediate future).

I haven't thought it through in detail, but I suspect the computations have extra steps. If so, it would be best (for reasons of performance) to separate this functionality from the evenly-spaced algorithms.

@timholy
Copy link
Owner

timholy commented Jul 10, 2013

This now has limited support via the InterpIrregular type. However, it's a pretty incomplete implementation (it's enough to do what I need right now...), so I am leaving this open.

@rsrock
Copy link
Contributor

rsrock commented Mar 19, 2014

I've run into a problem with InterpIrregular, it doesn't seem to handle the combination of BCnearest and InterpLinear.

# First, the tests
julia> using Grid

julia> Eps = sqrt(eps())
1.4901161193847656e-8

julia> x = [100.0,110.0,150.0]
3-element Array{Float64,1}:
 100.0
 110.0
 150.0

julia> y = rand(3)
3-element Array{Float64,1}:
 0.083211
 0.944004
 0.358623

julia> iu = InterpIrregular(x, y, -200, InterpNearest)
3-element InterpIrregular{Float64,1,BCfill,InterpNearest}:
 -200.0
 -200.0
 -200.0

julia> @assert iu[99] == -200

julia> @assert iu[101] == y[1]

julia> @assert iu[106] == y[2]

julia> @assert iu[149] == y[3]

julia> @assert iu[150.1] == -200

julia> iu = InterpIrregular(x, y, BCna, InterpLinear)
3-element InterpIrregular{Float64,1,BCna,InterpLinear}:
 NaN
 NaN
 NaN

julia> @assert isnan(iu[99])

julia> @assert abs(iu[101] - (0.9*y[1] + 0.1*y[2])) < Eps

julia> @assert abs(iu[106] - (0.4*y[1] + 0.6*y[2])) < Eps

julia> @assert abs(iu[149] - (y[2]/40 + (39/40)*y[3])) < Eps

julia> @assert isnan(iu[150.1])

# Everything above is fine, the problem is here
julia> iu = InterpIrregular(x, y, BCnearest, InterpLinear)
3-element InterpIrregular{Float64,1,BCnearest,InterpLinear}:
 #undef
 #undef
 #undef

julia> iu[99]
ERROR: no method _getindexii(InterpIrregular{Float64,1,BCnearest,InterpLinear}, Int64)
 in getindex at /Users/rrock/.julia/v0.3/Grid/src/interp.jl:248

julia> iu[100]
ERROR: no method _getindexii(InterpIrregular{Float64,1,BCnearest,InterpLinear}, Int64)
 in getindex at /Users/rrock/.julia/v0.3/Grid/src/interp.jl:248

julia> iu[101]
ERROR: no method _getindexii(InterpIrregular{Float64,1,BCnearest,InterpLinear}, Int64)
 in getindex at /Users/rrock/.julia/v0.3/Grid/src/interp.jl:248

@rsrock
Copy link
Contributor

rsrock commented Mar 19, 2014

Ok, I think I've fixed this. Pull request coming shortly...

#11

@CorySimon
Copy link

+1, unevenly spaced grid interpolation would be great.

@CorySimon
Copy link

I added a type that is similar to MATLAB's interp1d function. It only supports linear interpolation.

https://github.com/CorySimon/Grid.jl/blob/interp1d/src/interp1d.jl

How can I overload the interp1d.interpolate function so that it can take an array as well? Without the commented out portion, it just sees the interp1d.interpolate(x::Array) function and not the Float64 one.

What do you think?

You can use this by:

x = [0, 5.0, 6.0]
y = 2.0 * x
interp1d = Interp1d(x, y)
@printf("y=2x, x= %f, y(x) = %f\n" , .001, interp1d.interpolate(0.001))
@printf("y=2x, x= %f, y(x) = %f\n" , 5.5, interp1d.interpolate(5.5))

@mbauman
Copy link
Contributor

mbauman commented Mar 9, 2015

Have you seen InterpIrregular? I haven't looked at your interp1d, but I think it should do exactly what you want.

Create a grid like this:

x = cumsum(rand(100)/10)
y = sin(x)
g = Grid.InterpIrregular(x, y, Grid.BCnan, Grid.InterpLinear)

And then you can index into g with respect to your previous x-axis to get the interpolation result:

x2 = .1:.1:maximum(x)
y2 = g[x2]

I think this issue could probably be closed. Note, though, that this package is slated to be replaced by Interpolations.jl, which still needs some irregular interpolation love.

@CorySimon
Copy link

@mbauman Thanks, we should add this example to the README.md; I wasn't aware of this...

Sign up for free to subscribe to this conversation on GitHub. Already have an account? Sign in.
Labels
None yet
Projects
None yet
Development

No branches or pull requests

5 participants