[Numpy-discussion] interpolating arrays (?)

RayS rays at blue-cove.com
Mon Mar 22 21:39:00 EST 2004


At 08:19 PM 3/22/04 -0800, Warren Focke wrote:

>Would you consider Fourier interpolation?

I'll try it out tomorrow and check the results.
Thanks!
Ray

>interpolated = FFT.inverse_real_fft(FFT.real_fft(data), desiredLength) * \
>         desiredLength / len(data)
>
>This works whether desiredLength is less or greater than len(data).  The
>trick here is that inverse_real_fft discards data or zero-pads at high
>frequencies.  I keep toying with the idea of modifying the complex
>version to do the same, I cannot imagine a use for the current behavior
>(truncate or pad at small negative frequencies).
>
>There can be some pretty wierd edge effects, they can often be defeated
>with padding, or, better yet, overlapping segments.  Geting that right
>might be harder than fixing your approach.
>
>It has the nice feature that it is information-preserving:  if you
>interpolate to a higher number of points, then back to the original, you
>should get the same numbers back, give or take a little bit of rounding
>error.
>
>Warren Focke
>
>On Mon, 22 Mar 2004, Ray Schumacher wrote:
>
> > Howdy,
> >
> > I'd like to resize an array (1-dimensional) a of length n into an array of
> > n+x, and interpolate the values as necessary.
> > x might be+/-200, len(a) is usually ~500
> >
> > I tried it in pure Python but it's not quite right yet. (rubberBand, below)
> > Someone out there must have done this... and linear interpolation for
> > values is probably sufficient (instead of bi-cubic or some such)
> >
> > Cheers,
> > Ray
> >
> >
> >
> >
> > ==================================================================
> >          startSlice = 10000
> >          endOfSlice =
> > 10000+int(round(self.samplesPerRotation))    #self.samplesPerRotation=661
> >          ## self.dataArray is about len=200,000 of Int
> >
> >          for delta in range(-2,3):
> >              ## interpolate the different slices to the original length
> >         ##  0 should return the same array, this is a test
> >              newArray = self.rubberBand(
> > 
> self.dataArray[startSlice:endOfSlice+delta],
> >                                  round(self.samplesPerRotation))
> >              print len(self.dataArray[startSlice:endOfSlice+delta]),
> > len(newArray),
> >              ## show the result
> >              print delta, Numeric.sum(self.dataArray[startSlice:endOfSlice]
> > - newArray)
> >
> >      def rubberBand(self, data, desiredLength):
> >          """ alias/interpolate the data so that it is adjusted to the
> > desired array length"""
> >          currentLength = float(len(data))
> >
> >          ## positive if the new array is to be longer
> >          difference = desiredLength - currentLength
> >          #print  difference, desiredLength, currentLength
> >
> >          ## set up the desired length array
> >          newData = Numeric.zeros(desiredLength, Numeric.Float)
> >
> >          ## accounts for binary rounding errors
> >          smallNum = 2**-14
> >
> >          for index in range(desiredLength):
> >              ## find the ratio of the current index to the desired length
> >              ratio = index / float(desiredLength)
> >
> >              ## find the same (float) position in the old data
> >              currIndex = int(ratio * currentLength)
> >              ## find the decimal part
> >              decimalPart = (ratio * currentLength) - currIndex
> >              print index, currIndex, decimalPart,
> >              if(decimalPart>smallNum):
> >                  ## interpolate
> >                  newData[index] = ((1 - decimalPart) * data[currIndex] +
> > decimalPart * data[currIndex+1])
> >                  print data[currIndex], data[currIndex+1], newData[index]
> >              else:
> >                  newData[index] = data[currIndex]
> >                  print 'else',data[currIndex], newData[index]
> >
> >          return newData
> >
> >
> >
> > -------------------------------------------------------
> > This SF.Net email is sponsored by: IBM Linux Tutorials
> > Free Linux tutorial presented by Daniel Robbins, President and CEO of
> > GenToo technologies. Learn everything from fundamentals to system
> > administration.http://ads.osdn.com/?ad_id=1470&alloc_id=3638&op=click
> > _______________________________________________
> > Numpy-discussion mailing list
> > Numpy-discussion at lists.sourceforge.net
> > https://lists.sourceforge.net/lists/listinfo/numpy-discussion
> >





More information about the NumPy-Discussion mailing list