# Lookup table interpolation of an array

**URL:** <https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787>\
**Category:** New to Julia\
**Tags:** question\
**Created:** [May 3, 2024, 10:58am UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787 "2024-05-03T10:58:56Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![ahrensma](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahrensma/32/209243_2.png) [@ahrensma](https://discourse.julialang.org/u/ahrensma)\
**Post date:** [May 3, 2024, 10:58am UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/1 "2024-05-03T10:58:56Z")

</div>

Hello,

I am not yet fully familiar with Julia and am currently trying to rewrite a Matlab code. I am stumbling with the (Matlab) function interp1. Summarized the code looks like this:

A → Type: double(2048x1)  
B → Type: complex-single(2048,512)  
C → Type: double(1x2048)

D = interp1(A,B,C,“spline”)

Can anyone give me a tip on how I can write this in Julia, previous attempts have only led to wild error messages. Many thanks in advance.

Kind regards,

Martin

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [May 3, 2024, 11:41am UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/2 "2024-05-03T11:41:21Z")

</div>

Have you tried [Home · Interpolations.jl](http://juliamath.github.io/Interpolations.jl/latest/)?

---

<div class="post-metadata">

**Author:** ![mkitti](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkitti/32/12459_2.png) [@mkitti](https://discourse.julialang.org/u/mkitti)\
**Post date:** [May 3, 2024, 11:55am UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/3 "2024-05-03T11:55:48Z")

</div>

I would also look at [GitHub - jipolanco/BSplineKit.jl: A collection of B-spline tools in Julia](https://github.com/jipolanco/BSplineKit.jl) for irregular data.

I also list several other packages here:  
[https://juliamath.github.io/Interpolations.jl/latest/other\_packages/#Other-Interpolation-Packages](https://juliamath.github.io/Interpolations.jl/latest/other_packages/#Other-Interpolation-Packages)

Also post your wild error messages.

---

<div class="post-metadata">

**Author:** ![ahrensma](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahrensma/32/209243_2.png) [@ahrensma](https://discourse.julialang.org/u/ahrensma)\
**Post date:** [May 3, 2024, 6:43pm UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/4 "2024-05-03T18:43:06Z")

</div>

Thanks for the answer, I tried "interpolate\* from Interpolations.jl as well as from BSplineKit.jl. The code I use is

`using BSplineKit;`

with

A → typeof : Matrix{ComplexF64}, size : (2048, 1)  
B → typeof : Matrix{ComplexF64}, size : (2048, 512)  
C → typeof : Vector{Int64}, size : (2048,)

and in case of

`itp = interpolate(A[:,1], B[:,1], BSplineOrder(1))`

I received

`MethodError: no method matching make_knots(::Vector{ComplexF64}, ::BSplineOrder{1}, ::Nothing)`

In case of

`itp = interpolate(A, B, BSplineOrder(1))`

I received

`MethodError: no method matching interpolate(::Matrix{ComplexF64}, ::Matrix{ComplexF64}, ::BSplineOrder{1})`

I also tried

`spl = Spline1D(A, B , C , k=3, bc="nearest")`

which results in

`MethodError: no method matching Spline1D(::Matrix{ComplexF64}, ::Matrix{ComplexF64}, ::Vector{Int64}; k::Int64, bc::String)`

So, I’m not sure what is the best approach to replace the interp1 function from MatLab.

Kind regards,

Martin

---

<div class="post-metadata">

**Author:** ![bertschi](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bertschi/32/33462_2.png) [@bertschi](https://discourse.julialang.org/u/bertschi)\
**Post date:** [May 3, 2024, 10:18pm UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/5 "2024-05-03T22:18:51Z")

</div>

Cannot test Matlab directly, but from what I understand 1D interpolation needs real inputs, i.e., would not work on complex inputs `A`.  
Here is an example – assuming real inputs – using `BSplineKit`:

```julia
using Distributions, BSplineKit
N = 2048
A = sort(rand(Uniform(0, 10), N)) # Note: Real input vector
B1 = rand(ComplexF64, N) # Single complex output vector
sp1 = BSplineKit.interpolate(A, B1, BSplineKit.BSplineOrder(3)) # Get interpolation function
C = rand(Uniform(0, 10), N)
sp1.(C) # Note: Broadcast spline function over desired new inputs

# In contrast to Matlab, to handle multiple outputs, we need to be more explicit and broadcast manually
B = rand(ComplexF64, N, 512) # Matrix of multiple outputs
sps = BSplineKit.interpolate.(Ref(A), eachcol(B), Ref(BSplineKit.BSplineOrder(3))) # Get vector of splines interpolating each column of B
# Finally, apply each of these splines to all input values and combine results
sps .|> (f -> f.(C)) |> stack

```

---

<div class="post-metadata">

**Author:** ![ahrensma](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahrensma/32/209243_2.png) [@ahrensma](https://discourse.julialang.org/u/ahrensma)\
**Post date:** [May 4, 2024, 7:53pm UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/6 "2024-05-04T19:53:38Z")

</div>

Thank you for the example, it works well for real data; nevertheless I need to interpolate complex data …

---

<div class="post-metadata">

**Author:** ![mkitti](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkitti/32/12459_2.png) [@mkitti](https://discourse.julialang.org/u/mkitti)\
**Post date:** [May 4, 2024, 8:25pm UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/7 "2024-05-04T20:25:24Z")

</div>

He demonstrated with complex data. I’m unclear what is missing.

By the way, the above methods work more similarly to MATLAB’s gridded interpolant:

> **[Gridded data interpolation - MATLAB](https://www.mathworks.com/help/matlab/ref/griddedinterpolant.html)**
>
> Use griddedInterpolant to perform interpolation on a 1-D, 2-D, 3-D, or N-D gridded data set.

That is first we input the sample points and the values at those sample points and get an object that we can call to obtain other points.

Let’s breakdown what @bertschi demonstrated above. We first create a real vector of sample points in `A` and a complex 32-bit floating point in `B1`. These are both column vectors.

```julia-repl
julia> using Distributions, BSplineKit

julia> N = 2048
2048

julia> A = sort(rand(Uniform(0, 10), N));

julia> B1 = rand(ComplexF32, N)
2048-element Vector{ComplexF32}:
   0.8409336f0 + 0.19596744f0im
   0.9809899f0 + 0.5659446f0im
  0.79642314f0 + 0.0299052f0im
   0.5722098f0 + 0.89307314f0im
  0.96031356f0 + 0.1327033f0im
  0.71136796f0 + 0.862921f0im
   0.5338061f0 + 0.76945704f0im
   0.5737985f0 + 0.1366781f0im
  0.27081704f0 + 0.2231226f0im
  0.25618458f0 + 0.18294215f0im
               ⋮
  0.17688131f0 + 0.44664085f0im
   0.7323388f0 + 0.94128615f0im
  0.35286158f0 + 0.29547954f0im
   0.8608282f0 + 0.7629858f0im
 0.047004223f0 + 0.70044565f0im
   0.8486275f0 + 0.44268453f0im
   0.8413368f0 + 0.25111085f0im
  0.23773956f0 + 0.5144081f0im
   0.9474749f0 + 0.6850947f0im

julia> typeof(A)
Vector{Float64} (alias for Array{Float64, 1})

julia> typeof(B1)
Vector{ComplexF32} (alias for Array{Complex{Float32}, 1})

```

Then we create the `SplineInterpolation object`. Here I use `BSplineOrder(4)` for a cubic spline, which is what MATLAB uses.

```julia-repl
julia> sp1 = BSplineKit.interpolate(A, B1, BSplineKit.BSplineOrder(4))
SplineInterpolation containing the 2048-element Spline{ComplexF32}:
 basis: 2048-element BSplineBasis of order 4, domain [0.000857447419672086, 9.992905397874818]
 order: 4
 knots: [0.000857447, 0.000857447, 0.000857447, 0.000857447, 0.00779442, 0.0123642, 0.0176373, 0.0214538, 0.0264004, 0.0297231 … 9.9547, 9.95925, 9.96464, 9.96861, 9.97387, 9.98188, 9.99291, 9.99291, 9.99291, 9.99291]
 coefficients: ComplexF32[0.840934+0.195967im, 1.13049+1.21022im, 0.780093-1.07272im, 0.235346+1.87558im, 1.19961-0.603868im, 0.938954+1.03391im, 0.0414462+0.989143im, 2.25266-1.20595im, -0.204628+0.862867im, 0.804739-1.07816im … 2.00787+2.41364im, -2.50337-2.01189im, 5.14302+4.93931im, -1.16833-0.828984im, 1.90401+0.996272im, -0.889589+0.651825im, 2.48836+0.234553im, 1.07138+0.456499im, -4.38962-0.590771im, 0.947475+0.685095im]
 interpolation points: [0.000857447, 0.00305412, 0.00675145, 0.0135777, 0.0167635, 0.0225707, 0.0250273, 0.0316031, 0.032539, 0.0327446 … 9.94789, 9.95161, 9.95218, 9.9603, 9.96527, 9.96834, 9.97223, 9.98104, 9.99238, 9.99291]

```

We then query the `SplineInterpolation` object at various points via column vector `C` which is also randomly generated.

```julia-repl
julia> C = rand(Uniform(0,10), N)
2048-element Vector{Float64}:
 9.483733019514059
 5.43085430896525
 8.318840782645099
 5.288133871920023
 4.0813642928007
 6.78805933501573
 6.291783894803595
 5.170340091879966
 5.295940125355911
 6.382796375189246
 ⋮
 8.079857717607496
 4.940743863230217
 0.8397228433210491
 8.198655558986061
 5.212823437408405
 9.677186608905854
 1.7191592119233634
 9.625444117783816
 5.039112932203945

```

We can then query the spline as follows using the dot syntax.

```julia
julia> sp1.(C)
2048-element Vector{ComplexF64}:
   0.6628213223809246 + 0.13845306052661216im
    7.266404546511955 + 1.4918834901186397im
    2.529257327289871 - 3.0934900229509243im
    2.082829065843516 + 0.35194745428180896im
  0.42051821379475274 + 0.11991317115395697im
 -0.06217035130524348 + 0.42779375900348354im
    -1.29880373722419 + 0.410014463695485im
   0.7360561045406451 + 1.4965432408578285im
  0.29754144965019985 + 0.6998459950323489im
   0.8115314971604752 + 0.7237799271783811im
                      ⋮
   3.3251618089535615 - 5.701583525945178im
  -21.331010022778777 + 73.7176687071869im
   1.2751402913844252 + 1.268245785368036im
 -0.06294013906448183 + 1.631632769479881im
    0.766038339089081 + 1.10628834593872im
  0.20352877087166413 + 0.6675186029274851im
  -0.2737631672819028 + 0.2230333372408741im
  0.41110098498380204 + 0.1973548701882316im
    18.76530445342889 + 22.771700462031767im

julia> typeof(ans)
Vector{ComplexF64} (alias for Array{Complex{Float64}, 1})

```

We get a complex valued answer back.

For your original question, you provided multiple columns, 512 of them. This corresponds to 512 independent splines. bertschi’s second example demonstrated how to do this, but I will spell it out in a simple to use function for you that I will save in a file called `interp1.jl`.

```julia
import BSplineKit

function interp1(
    x::AbstractVector,
    v::AbstractVector,
    xq::AbstractVector,
    method::String = "spline"
)
    spline = interp1(x, v, method)
    return spline.(xq)
end

function interp1(
    x::AbstractVector,
    v::AbstractArray,
    xq::AbstractArray,
    method::String = "spline"
)
    splines = interp1(x, v, method)
    results = map(splines) do spline
        spline.(xq)
    end
    return stack(results)
end

function interp1(
    x::AbstractVector,
    v::AbstractVector,
    method::String = "spline"
)
    order = if method == "spline" || method == "cubic"
        BSplineKit.BSplineOrder(4)
    elseif method == "quadratic"
        BSplineKit.BSplineOrder(3)
    elseif method == "linear"
        BSplineKit.BSplineOrder(2)
    elseif method == "constant"
        BSplineKit.BSplineOrder(1)
    end

    return BSplineKit.interpolate(x, v, order)
end

function interp1(
    x::AbstractVector,
    v::AbstractArray,
    method::String = "spline"
)
   v_matrix = reshape(v, size(v, 1), prod(size(v)[2:end]))
   splines = map(eachcol(v_matrix)) do col
       interp1(x, col, method)
   end
   return splines
end

function interp1(
    x::AbstractVector,
    v::AbstractArray,
    method::String = "spline"
)
   v_matrix = reshape(v, size(v, 1), prod(size(v)[2:end]))
   splines = map(eachcol(v_matrix)) do col
       interp1(x, col, method)
   end
   return reshape(splines, size(v)[2:end])
end

function interp1(
    x::AbstractArray,
    v::AbstractArray,
    method::String = "spline"
)
   if !all(size(x) .== size(v))
      error("The size of x and v must the be same. Got $(size(x)) and $(size(v))")
   end
   x_matrix = reshape(v, size(x, 1), prod(size(x)[2:end]))
   v_matrix = reshape(v, size(v, 1), prod(size(v)[2:end]))
   splines = map(eachcol(x_matrix), eachcol(v_matrix)) do x_col, v_col
       interp1(x_col, v_col, method)
   end
   return reshape(splines, size(v)[2:end])
end

```

Now we can use this as follows.

```julia-repl
julia> using BSplineKit, Distributions

julia> include("interp1.jl")
interp1 (generic function with 10 methods)

julia> N = 2048; M = 512;

julia> A = sort(rand(Uniform(0,10), N)); size(A)
(2048,)

julia> B = rand(ComplexF32, N, M); size(B)
(2048, 512)

julia> C = rand(Uniform(0,10), N); size(C)
(2048,)

julia> D = interp1(A, B, C, "spline")
2048×512 Matrix{ComplexF64}:
  -3.56976+0.601276im -0.219194-0.447994im … -1.13067+0.473447im
  -6.33133-0.066161im 0.965402+2.17368im 0.11428+7.66711im
  0.265947+0.564525im 0.0940509+0.639697im 0.0618999+0.717675im
 -0.525911-0.268237im 0.030246+0.119171im 0.433661+0.519061im
    0.5511+0.237039im 0.589046+0.306925im -0.0207684+1.00121im
  0.601093+0.613032im 0.681022+0.923576im … 0.622523+0.701978im
 -0.310867+2.16622im 0.057356-0.517554im 2.94041+1.27067im
  0.513746+0.668517im -0.145421+0.247885im 0.966718-0.471837im
   1.90147+0.734572im 0.1871+1.62462im 3.60835-1.48335im
  0.747279+0.0694702im 1.03762+0.918896im 0.798624+0.0975493im
          ⋮ ⋱  
  -1.75224-0.431601im 1.33743-2.43919im 0.67119+2.0219im
   0.26693+0.842739im 0.0213721+0.180948im … 0.147373+0.755189im
   2.36164+0.665174im 0.710751+0.182283im 1.62779+1.1795im
   4.38767+3.99549im -2.68346-0.922896im -1.24872+0.693219im
  0.691161+0.844551im 0.664908+0.00167152im 0.635339+0.898001im
  0.552734+0.512182im 0.31708+0.480791im 0.807158+0.48275im
  0.570123+0.512312im 0.215266+0.285696im … 0.640996+0.0682268im
  0.322452-0.137867im -0.147389+0.277156im 0.00488787+0.61134im
 -0.204593+0.740612im 0.89808-0.270704im 1.12213-0.0284527im

```

Also you gave `C` as a row vector for some reason. This will also work. What I’m not entirely sure about is the expected size of the output `D` for your input.

```julia-repl
julia> C2 = transpose(C)
1×2048 transpose(::Vector{Float64}) with eltype Float64:
 9.15085 0.0196092 5.6302 3.97698 … 9.05284 8.40113 5.87534 7.91804

julia> D2 = interp1(A, B, C2, "spline")
1×2048×512 Array{ComplexF64, 3}:
...

```

---

<div class="post-metadata">

**Author:** ![ahrensma](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahrensma/32/209243_2.png) [@ahrensma](https://discourse.julialang.org/u/ahrensma)\
**Post date:** [May 5, 2024, 6:50pm UTC](https://discourse.julialang.org/t/lookup-table-interpolation-of-an-array/113787/8 "2024-05-05T18:50:10Z")

</div>

Hello,

I think now works. Thank you very much.

Kind regards,  
/Martin
