Consider following series of functions

 Do[f[k][w_] = Exp[-(w-Sin[k])^2], {k, 10}]

depicted as

series of dunctions

Let us tabulate and interpolate them

fdata = Flatten[Table[{w, k, f[k][w]}, {w, -3, 3, 0.1}, {k, 10}], 1];
g = Interpolation[fdata];

Instead of a smooth function some irregularities appear:


as shown here

3d image

Let us zoom in

 Evaluate[{g[w, 1], g[w, 1.5], g[w, 2]}], {w, -2, 3},
 PlotStyle -> {Black,Red,Black},Axes -> False,Frame->True]

The red curve should look similar to 2 black ones and have a peak somewhere in the middle. Instead, due to the 2d interpolation a much larger peak at wrong location arises.

3 curves

My question therefore would be: is there a way to do such interpolation in mathematica properly?

  • 1
    $\begingroup$ But ... you never really use the interpolation in your 3D plot $\endgroup$ – Sektor Mar 17 '16 at 14:07
  • $\begingroup$ @Sektor Indeed, my problem is interpolation, not the 3d plot itself. I am aware of the fact that ListPlot3D performs interpolation as well. If I have a good original function I would be able to get rid of plot-induced features by using a finer mesh. $\endgroup$ – yarchik Mar 17 '16 at 14:11

Assuming we dont know the underlying function, start with a 1D interpolation of each curve:

Do[ intf[k] = Interpolation[Table[{x, f[k][x]}, {x, -3, 3, .1}]] , {k,10}];

find the 'shift' from one curve to the next, assuming they are similar ( this is very slow )

diff[a_?NumericQ, i_, j_] := 
 NIntegrate[(intf[i][x] - intf[j][x + a])^2 , {x, -2, 2}]
offsets = Table[a /. Last@NMinimize[diff[a, i, i + 1]], {i, 9} ]

{0.0678265, -0.768177, -0.897923, -0.202122, 0.679509, 0.936402, 0.332372,-.57724, -0.95614}

now construct a 2d interpolation by linearly interpolating between the 1d interpolation functions:

intf2d[y_?NumericQ, x_] := 
 Module[ {k = Floor[y], ci = FractionalPart[y]},
      (1 - ci) intf[k][x - ci offsets[[k]]] + 
       ci  intf[k + 1][x + (1 - ci) offsets[[k]]]]

now we can let Plot3D sample it as needed to make a reasonably smooth plot

Plot3D[ intf2d[k, x], {k, 1, 10}, {x, -3, 3},PlotPoints->100]

enter image description here

superpose the original curves for validation:

  { Plot3D[intf2d[k, x], {k, 1, 10}, {x, -3, 3}, PlotPoints->100] , 
     Graphics3D[Table[{Thick, Green,
       Line[Table[{k, x, f[k][x]}, {x, -3, 3, .1}]]}, {k, 1, 10}]]}]

enter image description here

You could probably improve on this by fitting a smooth curve through the offsets.

  • $\begingroup$ Aha, so in other words you are suggesting to apply some linear transformation (shift) before doing interpolation. Maybe not perfect, but already something. $\endgroup$ – yarchik Mar 17 '16 at 17:39

What you're seeing is the result of polynomial interpolation. Here is a simple 1D example taken from Wikipedia

Notice specifically that the curve between points $1$ and $2$, as well as between $4$ and $5$ isn't actually between those samples. Still, this seems to be a better fit of the data, even though it's not an optimal local interpolation. The reason to use polynomial interpolation is to improve the smoothness of the overall fit.

If your main concern is that you get a straight, linear interpolation between the individual curves, and don't care about smoothness, then you should tell Interpolate to do so, by setting the InterpolationOrder parameter to 1 (the default is 3):

Do[f[k][w_] = Exp[-(w - Sin[k])^2], {k, 10}]
fdata = Flatten[Table[{w, k, f[k][w]}, {w, -3, 3, 0.1}, {k, 10}], 1];
g = Interpolation[fdata, InterpolationOrder -> 1];

Plot[Evaluate[{g[w, 1], g[w, 1.5], g[w, 2]}], {w, -2, 3}, 

PlotStyle -> {Black, Red, Black}, Axes -> False, Frame -> True]

enter image description here

Also note that the 3rd-order interpolation of Interpolate is actually smoother than the samples displayed by ListPlot3D:

g = Interpolation[fdata];
Plot3D[Evaluate[g[w, k]], {w, -3, 3}, {k, 1, 10}]

enter image description here

So it really depends on what you're looking for. If you want a smooth interpolation, there are bound to be adjacent samples where not every interpolated curve lies exactly between. If you want to ensure that, you'll have to use linear interpolation, but then your surface won't be smooth.

  • $\begingroup$ I actually want both, a smooth interpolation with interpolated curve in between based on a priory knowledge that the shape of all functions in the series is similar, but the curves may be shifted. $\endgroup$ – yarchik Mar 17 '16 at 14:31
  • $\begingroup$ @yarchik I don't think that is going to happen with Interpolation. (And in fact, if you fit the correct solution, i.e. a sine wave through your samples you will notice that even the correct solution, at its extrema, will not always lie between the two adjacent samples.) I don't know what other what other interpolation functions Mathematica provides, but if that's not what you're looking for you might have to roll your own. $\endgroup$ – Martin Ender Mar 17 '16 at 14:47
  • $\begingroup$ Something like a smoothstep interpolation might work (it'll be smooth and always between two adjacent samples), although I think the resulting fit would look a lot worse than the polynomial interpolation you already have. $\endgroup$ – Martin Ender Mar 17 '16 at 14:47

I think you need to put more points in between to get what you want:

i = 1; Do[(f[i++][w_] = Exp[-(w - Sin[k])^2]), {k, 1, 10, .1}]
fdata = Table[{w, k, f[k][w]}, {w, -3, 3, 0.1}, {k, i}];
ListPlot3D[Flatten[fdata, 1], InterpolationOrder -> 3, 
 PlotRange -> All]

enter image description here

  • $\begingroup$ This actually defeats the purpose of my example. For real world case it is overshoot and inefficient to use so many points. $\endgroup$ – yarchik Mar 17 '16 at 16:22

Not an answer but a comment (with code) that hopefully illuminates the problem.

This demonstrates that ad-hoc connection into polygons of the points generated regularly over the curves brings a result similar to that of ListPlot3D:

Block[{fdata = fdata},
 fdata = SortBy[#, #[[1]] &] & /@ 
   SortBy[GatherBy[fdata, #[[2]] &], #[[All, 2]] &];
 ps = Map[
    MapThread[Polygon[Join[#1, #2]] &, 
     Partition[#, 2, 1] & /@ rpair]], Partition[fdata, 2, 1]];
 Graphics3D[{EdgeForm[Black], ps}, BoxRatios -> {1, 1, 1/2}]

enter image description here

  • $\begingroup$ Would you then say (based on your plot) that linear interpolation is superior for such functions? $\endgroup$ – yarchik Mar 18 '16 at 7:26
  • $\begingroup$ @yarchik Hmm... No. I would say linear interpolation is also inadequate. I was considering giving an answer similar to that of george2029. $\endgroup$ – Anton Antonov Mar 18 '16 at 12:50
  • $\begingroup$ @yarchik I see the code and plot I put as a way to explain the results from ListPlot3D. If we compare the plots there are obvious similarities. $\endgroup$ – Anton Antonov Mar 18 '16 at 12:54
  • $\begingroup$ yes, indeed, there are similarities. I think a lot can be done along the lines suggested by @george2079, but it is probably not easy to find a general solution. $\endgroup$ – yarchik Mar 18 '16 at 17:09

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.