From: <halldor@ve...>  20060330 11:33:37
Attachments:
niceinterp.py

First, I want to thank John Hunter and Jouni Seppanen for their very=20 useful help with my xticklabel problems. I have little experience with=20 python, and am a total newbie wrt. matplotlib ( which is very=20 impressive...). Few days ago I was looking for a method to interpolate 12hourly data,=20 and was not happy with the python alternatives I found. I wanted=20 something in native python, not wrappers to C or Fortran programs (in this case speed is not a big issue). It so happened that few weeks=20 ago a coworker gave me a faint photocopy of a photocopy from an article=20 published 26 years ago with an interpolation method I had not seen before. The main point of the method was that it did not generate spline like=20 excursions. The claim was that the method was "consistantly well behaved". So, I wrote it up in Python (well, I translated my Matlab code). It=20 seems to work as the author claimed. The interpolation method needs x,y=20 and y' (the slope of y(x)) as input, but the article also provided a=20 simple algorithm to estimate the slope. I am including these functions=20 in an attachment. To use them is very simple: Example where y' is known x =3D seq(0,2*pi,10); y =3D sin(x); yp =3D cos(x) xi =3D seq(0,2*pi,40); yi =3D StinemanInterp(x,y,yp,xi); plot(x,y,'o',xi,yi) Example where y' is calculated using slopes() x =3D seq(0,2*pi,10); y =3D sin(x); yp =3D slopes(x,y) xi =3D seq(0,2*pi,40); yi =3D StinemanInterp(x,y,yp,xi); plot(x,y,'o',xi,yi) I make no claim that the method is fast, nor that my coding is optimal.=20  I am sure it can be sped up, but it is fast enough for my uses. Nor do I claim that this is the best interpolation algorithm,  I have=20 no interest in taking part in a religous debate on the merits of=20 different interpolation methods. For me this is simply a consise method that works resonably well. The=20 slope y' can be calculated using other methods than slopes(), using=20 divided differences to get the derivative of a quadratic through any=20 three points is straightforward (but not implemented here). Also if yp=20 is simply the linear slope the method yields linear interpolation. As I said, I have not seen this method before, but it may well be that=20 it is known under a different guise....My coworker got it from a Swiss=20 glaciologist.... Sincerely, Halld=F3r =20 Halldor Bjornsson Bustadavegur 9 IS150 Reykjavik Iceland. e: halldor()vedur.is tel:+354522600 
From: John Hunter <jdhunter@ac...>  20060330 15:22:58

>>>>> "Halldor" =3D=3D Halldor Bj=F6rnsson <halldor@...> writes: Halldor> First, I want to thank John Hunter and Jouni Seppanen for Halldor> their very useful help with my xticklabel problems. I Halldor> have little experience with python, and am a total newbie Halldor> wrt. matplotlib ( which is very impressive...). Thanks! Halldor> Few days ago I was looking for a method to interpolate Halldor> 12hourly data, and was not happy with the python Halldor> alternatives I found. I wanted something in native Halldor> python, not wrappers to C or Fortran programs (in this I think the best reason for writing it up in python is that it is fun and instructive, but I don't agree with the "not wrappers of C or Fortran". Note that pylab is built on top of Numeric, which is written in C, the font handling is built on top of freetype, another C library, and the antigrain renderer engine is a C++ library. matplotlib ships with a fair amount of C/C++ code, and scipy does the same for C, C++ and FORTRAN. The strength of python is its fluid integration with other languages. Also, the agg engine provides many 2d interpolation algorithms: nearest, bilinear, bicubic, spline16, spline36, hanning, hamming, hermite, kaiser, quadric, catrom, gaussian, bessel, mitchell, sinc, lanczos, blackman It would be nice to provide a more direct interface to these, and other agg image functionality, in matplotlib. Thanks for your code  you may also want to look at scipy.interpolate NAME scipy.interpolate FILE /usr/lib/python2.4/sitepackages/scipy/interpolate/__init__.py DESCRIPTION Interpolation Tools =3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D Wrappers around FITPACK functions: splrep  find smoothing spline given (x,y) points on curve. splprep  find smoothing spline given parametrically defined curve. splev  evaluate the spline or its derivatives. splint  compute definite integral of a spline. sproot  find the roots of a cubic spline. spalde  compute all derivatives of a spline at given points. bisplrep  find bivariate smoothing spline representation. bisplev  evaluate bivariate smoothing spline. Interpolation Class interp1d  Create a class whose instances can linearly interpolate to compute unknown values. 
From: John Hunter <jdhunter@ac...>  20060331 05:22:08

>>>>> "Halldor" =3D=3D Halldor Bj=F6rnsson <halldor@...> writes: Halldor> The main point of the method was that it did not generate Halldor> spline like excursions. The claim was that the method was Halldor> "consistantly well behaved". I took a second look at this, and am curious about the notion of "consistently well behaved". Do you have some example use cases where this algorithm "does the right thing" but typical interpolation algorithms like splines fail? Could you look over the scipy.interpolate algorithms and see if this one does something that they do not? test cases in which your algorithm produced good results but scipy.interpolate algorithms fail would be particularly useful. Are the excursions you refer to the case of wild extrapolations beyond the end points of the data being interpolated, or wild fluctuations between interpolated points? If you can demonstrate some usecases, a cleaned up version might be a useful contribution to scipy.interpolate. Cheers, JDH 
From: Norbert Nemec <Norbert.N<emec.list@gm...>  20060331 13:03:36
Attachments:
stinemaninterp.py

Nice work, Halldor! I've spent a bit of time on data interpolation recently, but this Stineman interpolation method beats everything I came up with in quality and simplicity. I took the freedom of going over your code and putting in all the experience I gathered before working with Python on data interpolation and related issues. See attached what I came up with. Compared with your version this is heavily modified: * The core change was to get rid of the explicit loop in the interpolation routine. The method now fully exploits the power of numpy. * The interface of the interpolation routine is changed. Intention was to move yp to the back and make it optional. If it is not given, the "slopes" routine is called automatically. * The "slopes" routine is changed in its core. Instead of a circular interpolation (which is problematic if the aspect ratio between the x and the yaxis is not known) it now uses a parabola interpolation to estimate the slope at inner points. For most curves where the original gave reasonable results, this new version should give something similar. It should, however, also cover many cases where the original "slopes" routine would have produced garbage. * The slopes at the endpoints are now also extrapolated in a much simpler manner. The intent of the original version was not clear to me, but it would definitely caused problems in several corner cases. If you could privately send my a scan of the original paper, I would be very grateful. I believe, it would be quite some gain for the matplotlib library to have this algorithm incorporated. As far as I can see, it should indeed be very robust for producing nice looking interpolations. One should, however, warn against its use in scientific context: Interpolation is always a way of making data look nicer than it actually was measured... Greetings, Norbert 
From: John Hunter <jdhunter@ac...>  20060331 18:20:24

>>>>> "Norbert" == Norbert Nemec <Norbert.Nemec.list@...> writes: Norbert> Nice work, Halldor! I've spent a bit of time on data Norbert> interpolation recently, but this Stineman interpolation Norbert> method beats everything I came up with in quality and Norbert> simplicity. Wow, you two make a potent team! Norbert, I took your revision and made some minor modifications to support numeric and numarray and added it to matplotlib.mlab. It was a crime to take your nice code utilizing numpy/numarray array indexing and backport it to the Numeric.take ugly equivalents, but that's the price we pay for supporting Numeric. See examples/interp_demo.py in svn revision 2244 JDH 
From: Norbert Nemec <Norbert.N<emec.list@gm...>  20060401 08:32:57

I just wanted to start this reply with the words "Very nice!", when I realized that "nice" probably is the last thing one would say about the "uglyfication" necessary for numerix compatibility. :) So lets instead say: "Good work!" and thank you for taking in the code. Greetings, Norbert John Hunter wrote: > > Norbert> Nice work, Halldor! I've spent a bit of time on data > Norbert> interpolation recently, but this Stineman interpolation > Norbert> method beats everything I came up with in quality and > Norbert> simplicity. > >Wow, you two make a potent team! > >Norbert, I took your revision and made some minor modifications to >support numeric and numarray and added it to matplotlib.mlab. It was >a crime to take your nice code utilizing numpy/numarray array indexing >and backport it to the Numeric.take ugly equivalents, but that's the >price we pay for supporting Numeric. > >See examples/interp_demo.py in svn revision 2244 > >JDH > > > >This SF.Net email is sponsored by xPML, a groundbreaking scripting language >that extends applications into web and mobile media. Attend the live webcast >and join the prime developer group breaking into this new coding territory! >http://sel.asus.falkag.net/sel?cmd=lnk&kid=110944&bid=241720&dat=121642 >_______________________________________________ >Matplotlibusers mailing list >Matplotlibusers@... >https://lists.sourceforge.net/lists/listinfo/matplotlibusers > > > > 