Hi John, Hi Eric,

In case this has slipped under the radar, I repost my patch for bilinear interp 
in NonUniformImage (including the discussion around it). I think Eric is the 
most concerned by this, but you were both 
on holidays when I raised the issue in matplotlib newsgroup.

Best regards,

Greg.

-------------


Thanks a lot for reviewing my patch!
I have corrected most of the problems (I think ;-) )
I indeed introduced memory leak, I think it is fixed now, I have also
reorganized the code to avoid duplication of the cleanup code. I used an
helper function instead of the goto, because this cleanup is needed for
normal exit too, so the helper function can come in handy for that.
std::vector would have been nice, but I tried to respect the original
coding style, maybe this could be changed but I guess the original
author should have something to say about it....
Only thing I did not change is the allocation of acols/arows even when
they are not used. It is not difficult to do, but the code will be a
little more complex and I think that the extra memory used is not
significant: The memory for those mapping structure is doubled (1 float
and 1 int vector of size N, instead of a single int vector of size N), 
but those mapping structures are an order of magnitude smaller than the
image buffer of size N*N of 4 char anyway...

If one agree that this is to be saved at the (slight) expense of code
conciseness/simplicity, I can add this optimisation...

Thanks for pointing the error at the left/bottom pixel lines, it is
corrected now :-)

And I set interpolation to NEAREST by default, so that the behavior of
NonUniformImage remains the same as before if the user do not specify
otherwise (I forgot to do it in the first patch)...


I include the new patch,

Best regards,

Greg.

On Fri, 2008-08-15 at 15:45 -0400, Michael Droettboom wrote:
> Thanks for all the work you put into this patch.
> 
> As you say, it would be nice to have a generic framework for this so 
> that new types of interpolation could be easily added and to be able to 
> support arbitrary (non-axis-aligned) quadmeshes as well.  But that's 
> even more work -- if we keep waiting for everything we want, we'll never 
> get it... ;)  I agree that Agg probably won't be much help with that. 
> 
> There are a couple of comments with the patch as it stands -->
> 
> There seems to be a gap extrapolating over the left and bottom edge (see 
> attached screenshot from pcolor_nonuniform.py).
> 
> Memory management looks problematic, some of which I think you inherited 
> from earlier code.  For example, arows and acols are never freed.  
> Personally, I think these temporary buffers should be std::vector's so 
> they'll be free'd automatically when scope is left.  It might also be 
> nice to move all of the Py_XDECREF's that happen when exceptions are 
> thrown to either a master try/catch block or an "exit" goto label at the 
> bottom.  The amount of duplication and care required to ensure 
> everything will be freed by all of the different exit paths is a little 
> cumbersome.
> 
> Also, acols and arows are only used in BILINEAR interpolation, but they 
> are allocated always.
> 
> Once these issues are addressed, it would be great to have someone who 
> *uses* the nonuniform pcolor functionality (Eric Firing?) have a look at 
> this patch for any regressions etc..  Assuming none, I'll be happy to 
> commit it (but I won't be around for a week or so).
> 
> Cheers,
> Mike
> 
> Grégory Lielens wrote:
> > Hi all,
> >
> > here is a patch which implement bilinear interpolation on irregular grid
> > ( i.e. it allows  NonUniformImage
> > to accept  both 'nearest' and 'bilinear' interpoaltion, instead of only
> > 'nearest'.)
> >
> > It is not perfect, given the current architecture of the image module  I
> > think there is not simple way
> > to specify other interpolations (except for 'nearest', all
> > interpolations are implemented at the AGG level, not in
> > matplotlib itself).
> > For the same reason, it is not possible to interpolate before colormap
> > lookup instead of after (and this 
> > is usually where the biggest effect is, when dealing with coarse grid).
> > However, I think it is already quite useful and the best one ca do
> > without a -very- extensive rewrite
> > of the matrix map modules....
> >
> > BTW, I think it also corrects a small bug in the 'nearest'
> > interpolation: the last intervals was ignored
> > in the CVS version, now it is taken into account.
> >
> > BOTH nearest and bilinear interpolation do the same for values oustside
> > the data grid: they just use boundary values (constant extrapolation).
> > I avoided linear extrapolation for 'bilinear', as extrapolation is
> > usually dangerous or meaningless (we can go outside the colormap very
> > fast)...
> >
> > I have included a small example showing how both interpolation works....
> >
> > Any remarks, could this be added before the next release? ;-)
> >
> >
> > Greg.
> >
> >
> >
> >
> >
> > On Mon, 2008-08-11 at 16:50 +0200, Grégory Lielens wrote:
> >   
> >> On Fri, 2008-08-08 at 16:05 +0200, Grégory Lielens wrote:
> >>     
> >>> Hello everybody,
> >>>
> >>> I have sent this message to the user group, but thinking of it, it may be 
> >>> more
> >>> relevant to the development mailing list...so here it is again.
> >>>
> >>>
> >>>
> >>> We are looking for the best way to plot a waterfall diagram in
> >>> Matplotlib. The 2 functions which could be used 
> >>> to do that are (as far as I have found) imshow and pcolormesh. Here is a
> >>> small script that use both to compare the output:
> >>>
> >>> -----------------
> >>>
> >>> from pylab import *
> >>>
> >>>
> >>> delta = 0.2
> >>> x = arange(-3.0, 3.0, delta)
> >>> y = arange(-2.0, 2.0, delta)
> >>> X, Y = meshgrid(x, y)
> >>> Z1 = bivariate_normal(X, Y, 1.0, 1.0, 0.0, 0.0)
> >>> Z2 = bivariate_normal(X, Y, 1.5, 0.5, 1, 1)
> >>> # difference of Gaussians
> >>> Z = 10.0 * (Z2 - Z1)
> >>> figure(1)
> >>> im = imshow(Z,extent=(-3,3,-2,2))
> >>> CS = contour(X, -Y, Z, 6,
> >>>              colors='k', # negative contours will be dashed by default
> >>>              )
> >>> clabel(CS, fontsize=9, inline=1)
> >>> title('Using imshow')
> >>> figure(2)
> >>> im = pcolormesh(X,-Y,Z)
> >>> CS = contour(X, -Y, Z, 6,
> >>>              colors='k', # negative contours will be dashed by default
> >>>              )
> >>> clabel(CS, fontsize=9, inline=1)
> >>> title('Using pcolormesh')
> >>> show()
> >>>
> >>> ---------------------
> >>>
> >>>
> >>> The problem is that we need some of the flexibility of pcolormesh (which
> >>> is able to map the matrix of value on any deformed mesh), while
> >>> we would like to use the interpolations available in imshow (which
> >>> explain why the imshow version is much "smoother" than the pcolormesh
> >>> one).
> >>>
> >>> In fact, what would be needed is not the full flexibility of pcolormesh
> >>> (which can map the grid to any kind of shape), we "only" have to deal
> >>> with rectangular grids where x- and y- graduations are irregularly spaced.
> >>>
> >>> Is there a drawing function in Matplotlib which would be able to work
> >>> with such a rectangular non-uniform grid?
> >>> And if not (and a quick look at the example and the code make me think 
> >>> that indeed the capability is currently not present),
> >>> what about an extension of imshow which would work as this:
> >>>  
> >>> im = imshow(Z,x_gridpos=x, y_gridpos=y)  #specify the
> >>> position of the grid's nodes, instead of giving the extend and assuming
> >>> uniform spacing.
> >>>
> >>> Longer term, would a pcolormesh accepting interpolation be possible? The
> >>> current behavior, averaging the color of the grids node to get a uniform
> >>> cell color, 
> >>> is quite rough except for a large number of cells...And even then, it
> >>> soon shows when you zoom in...
> >>>
> >>> The best would be to allow the same interpolations as in imshow (or a
> >>> subset of it), and also allows to use interpolation before colormap
> >>> lookup (or after), 
> >>> like in Matlab. Indeed, Matlab allows to finely tune interpolation by
> >>> specifying Gouraud (interpolation after color
> >>> lookup)/Phong(interpolation before color lookup, i.e. for each pixel).
> >>> Phong is usually much better but also more CPU intensive. Phong is
> >>> especially when using discrete colormap, producing  banded colors
> >>> equivalent to countour lines, while Gouraud does not work in those
> >>> cases.
> >>>
> >>> Of course, the performance will be impacted by some of those
> >>> interpolation options, which would degrade performance in animations for
> >>> example.... but I think that having the different options available
> >>> would be very useful, it allows to have the highest map quality, or have
> >>> a "quick and dirty" map depending on situation (grid spacing, type of
> >>> map, animation or not, ...).
> >>>
> >>> Best regards,
> >>>
> >>> Greg.
> >>>       
> >> I have found a method which implement the proposed extension to imshow:
> >> NonUniformImage...
> >>
> >> However, this image instance support only nearest neighbor
> >> interpolation. Trying to set the interpolation (using the
> >> set_interpolation method)
> >> to something akin imshow throw a "NotImplementedError: Only nearest
> >> neighbor supported" exception....
> >>
> >> So basically I am still stuck, it seems that currently there is no way
> >> in matplotlib to plot interpolated
> >> colormap on irregular rectangular grid, and even less on arbitrarily
> >> mapped grid...
> >>
> >> Is there any plans to add support for more interpolation in
> >> NonUniformImage in the future? Or maybe there is another 
> >> drawing function that I did not find yet, with this ability?
> >>
> >>
> >> Best regards,
> >>
> >> Greg.
Index: lib/matplotlib/image.py
===================================================================
--- lib/matplotlib/image.py	(revision 6039)
+++ lib/matplotlib/image.py	(working copy)
@@ -394,6 +394,7 @@
                 ):
         AxesImage.__init__(self, ax,
                            **kwargs)
+        AxesImage.set_interpolation(self, 'nearest')
 
     def make_image(self, magnification=1.0):
         if self._A is None:
@@ -407,7 +408,9 @@
         height *= magnification
         im = _image.pcolor(self._Ax, self._Ay, self._A,
                            height, width,
-                           (x0, x0+v_width, y0, y0+v_height))
+                           (x0, x0+v_width, y0, y0+v_height),
+                           self._interpd[self._interpolation])
+
         fc = self.axes.patch.get_facecolor()
         bg = mcolors.colorConverter.to_rgba(fc, 0)
         im.set_bg(*bg)
@@ -450,8 +453,8 @@
         raise NotImplementedError('Method not supported')
 
     def set_interpolation(self, s):
-        if s != None and s != 'nearest':
-            raise NotImplementedError('Only nearest neighbor supported')
+        if s != None and not s in ('nearest','bilinear'):
+            raise NotImplementedError('Only nearest neighbor and bilinear interpolations are supported')
         AxesImage.set_interpolation(self, s)
 
     def get_extent(self):
Index: src/_image.cpp
===================================================================
--- src/_image.cpp	(revision 6039)
+++ src/_image.cpp	(working copy)
@@ -1138,21 +1138,211 @@
   return Py::asObject(imo);
 }
 
+// utilities for irregular grids
+void _bin_indices_middle(unsigned int *irows, int nrows, float *ys1, int ny,float dy, float y_min)
+{ int  i,j, j_last;
+   unsigned  int * rowstart = irows;
+  float *ys2 = ys1+1;
+  float *yl = ys1 + ny ;
+  float yo = y_min + dy/2.0;
+  float ym = 0.5f*(*ys1 + *ys2);
+  // y/rows
+   j = 0;
+   j_last = j;
+  for (i=0;i<nrows;i++,yo+=dy,rowstart++) {
+      while(ys2 != yl && yo > ym) {
+          ys1 = ys2;
+          ys2 = ys1+1;
+          ym = 0.5f*(*ys1 + *ys2);
+          j++;
+      }
+      *rowstart = j - j_last;
+      j_last = j;
+  }
+}
 
+void _bin_indices_middle_linear(float *arows, unsigned int *irows, int nrows, float *y, int ny,float dy, float y_min)
+{ int i;
+    int ii = 0;
+    int iilast = ny-1;
+    float sc = 1/dy;
+    int iy0 = (int)floor(sc * (y[ii]  - y_min) );
+    int iy1 = (int)floor(sc * (y[ii+1]  - y_min) );
+    float invgap=1.0f/(iy1-iy0);
+    for (i=0; i<nrows && i<=iy0; i++) {
+        irows[i] = 0;
+        arows[i] = 1.0;
+        //std::cerr<<"i="<<i<<"  ii="<<0<<" a="<< arows[i]<< std::endl;
+    }
+    for (; i<nrows; i++) {
+        while (i > iy1 && ii < iilast) {
+            ii++;
+            iy0 = iy1;
+            iy1 = (int)floor(sc * (y[ii+1] - y_min));
+            invgap=1.0f/(iy1-iy0);
+        }
+        if (i >= iy0 && i <= iy1) {
+            irows[i] = ii;
+            arows[i]=(iy1-i)*invgap;
+            //std::cerr<<"i="<<i<<"  ii="<<ii<<" a="<< arows[i]<< std::endl;
+        }
+        else break;
+    }
+    for (; i<nrows; i++) {
+        irows[i] =iilast-1;
+        arows[i] = 0.0;
+        //std::cerr<<"i="<<i<<"  ii="<<iilast-1<<" a="<< arows[i]<< std::endl;
+    }
+}
+
+void _bin_indices(int *irows, int nrows, double *y, int ny,
+                                            double sc, double offs)
+{
+    int i;
+    if (sc*(y[ny-1] - y[0]) > 0)
+    {
+        int ii = 0;
+        int iilast = ny-1;
+        int iy0 = (int)floor(sc * (y[ii]  - offs));
+        int iy1 = (int)floor(sc * (y[ii+1]  - offs));
+        for (i=0; i<nrows && i<iy0; i++) {
+            irows[i] = -1;
+        }
+        for (; i<nrows; i++) {
+            while (i > iy1 && ii < iilast) {
+                ii++;
+                iy0 = iy1;
+                iy1 = (int)floor(sc * (y[ii+1] - offs));
+            }
+            if (i >= iy0 && i <= iy1) irows[i] = ii;
+            else break;
+        }
+        for (; i<nrows; i++) {
+            irows[i] = -1;
+        }
+    }
+    else
+    {
+        int iilast = ny-1;
+        int ii = iilast;
+        int iy0 = (int)floor(sc * (y[ii]  - offs));
+        int iy1 = (int)floor(sc * (y[ii-1]  - offs));
+        for (i=0; i<nrows && i<iy0; i++) {
+            irows[i] = -1;
+        }
+        for (; i<nrows; i++) {
+            while (i > iy1 && ii > 1) {
+                ii--;
+                iy0 = iy1;
+                iy1 = (int)floor(sc * (y[ii-1] - offs));
+            }
+            if (i >= iy0 && i <= iy1) irows[i] = ii-1;
+            else break;
+        }
+        for (; i<nrows; i++) {
+            irows[i] = -1;
+        }
+    }
+}
+
+void _bin_indices_linear(float *arows, int *irows, int nrows, double *y, int ny,
+                                            double sc, double offs)
+{
+    int i;
+    if (sc*(y[ny-1] - y[0]) > 0)
+    {
+        int ii = 0;
+        int iilast = ny-1;
+        int iy0 = (int)floor(sc * (y[ii]  - offs));
+        int iy1 = (int)floor(sc * (y[ii+1]  - offs));
+        float invgap=1.0/(iy1-iy0);
+        for (i=0; i<nrows && i<iy0; i++) {
+            irows[i] = -1;
+        }
+        for (; i<nrows; i++) {
+            while (i > iy1 && ii < iilast) {
+                ii++;
+                iy0 = iy1;
+                iy1 = (int)floor(sc * (y[ii+1] - offs));
+                invgap=1.0/(iy1-iy0);
+            }
+            if (i >= iy0 && i <= iy1) {
+                irows[i] = ii;
+                arows[i]=(iy1-i)*invgap;
+            }
+            else break;
+        }
+        for (; i<nrows; i++) {
+            irows[i] = -1;
+        }
+    }
+    else
+    {
+        int iilast = ny-1;
+        int ii = iilast;
+        int iy0 = (int)floor(sc * (y[ii]  - offs));
+        int iy1 = (int)floor(sc * (y[ii-1]  - offs));
+        float invgap=1.0/(iy1-iy0);
+        for (i=0; i<nrows && i<iy0; i++) {
+            irows[i] = -1;
+        }
+        for (; i<nrows; i++) {
+            while (i > iy1 && ii > 1) {
+                ii--;
+                iy0 = iy1;
+                iy1 = (int)floor(sc * (y[ii-1] - offs));
+                invgap=1.0/(iy1-iy0);
+            }
+            if (i >= iy0 && i <= iy1) {
+                irows[i] = ii-1;
+                arows[i]=(i-iy0)*invgap;
+            }
+            else break;
+        }
+        for (; i<nrows; i++) {
+            irows[i] = -1;
+        }
+    }
+}
+
+
+
 char __image_module_pcolor__doc__[] =
 "pcolor(x, y, data, rows, cols, bounds)\n"
 "\n"
-"Generate a psudo-color image from data on a non-univorm grid using\n"
-"nearest neighbour interpolation.\n"
+"Generate a pseudo-color image from data on a non-uniform grid using\n"
+"nearest neighbour or linear interpolation.\n"
 "bounds = (x_min, x_max, y_min, y_max)\n"
+"interpolation = NEAREST or BILINEAR \n"
 ;
+
+void _pcolor_cleanup(PyArrayObject* x, PyArrayObject* y,  PyArrayObject *d,
+                                        unsigned int * rowstarts ,unsigned int*colstarts ,
+                                        float *acols , float *arows) {
+    if (x)
+      Py_XDECREF(x);
+    if (y)
+      Py_XDECREF(y);
+    if(d)
+      Py_XDECREF(d);
+    if(rowstarts)
+      PyMem_Free(rowstarts);
+    if(colstarts)
+      PyMem_Free(colstarts);
+    if(acols)
+      PyMem_Free(acols);
+    if(arows)
+      PyMem_Free(arows);
+    return;
+}
+    
 Py::Object
 _image_module::pcolor(const Py::Tuple& args) {
   _VERBOSE("_image_module::pcolor");
 
 
-  if (args.length() != 6)
-      throw Py::TypeError("Incorrect number of arguments (6 expected)");
+  if (args.length() != 7)
+      throw Py::TypeError("Incorrect number of arguments (7 expected)");
 
   Py::Object xp = args[0];
   Py::Object yp = args[1];
@@ -1160,6 +1350,7 @@
   unsigned int rows = Py::Int(args[3]);
   unsigned int cols = Py::Int(args[4]);
   Py::Tuple bounds = args[5];
+  unsigned int interpolation = Py::Int(args[6]);
 
   if (rows > 1 << 15 || cols > 1 << 15) {
     throw Py::ValueError("rows and cols must both be less than 32768");
@@ -1179,26 +1370,29 @@
   // Check we have something to output to
   if (rows == 0 || cols ==0)
       throw Py::ValueError("Cannot scale to zero size");
-
+  
+  PyArrayObject *x = NULL; PyArrayObject *y = NULL; PyArrayObject *d = NULL;
+  unsigned int * rowstarts = NULL; unsigned int*colstarts = NULL;
+  float *acols = NULL; float *arows = NULL;
+  
   // Get numpy arrays
-  PyArrayObject *x = (PyArrayObject *) PyArray_ContiguousFromObject(xp.ptr(), PyArray_FLOAT, 1, 1);
-  if (x == NULL)
+  x = (PyArrayObject *) PyArray_ContiguousFromObject(xp.ptr(), PyArray_FLOAT, 1, 1);
+  if (x == NULL) {
+      _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::ValueError("x is of incorrect type (wanted 1D float)");
-  PyArrayObject *y = (PyArrayObject *) PyArray_ContiguousFromObject(yp.ptr(), PyArray_FLOAT, 1, 1);
+  }
+  y = (PyArrayObject *) PyArray_ContiguousFromObject(yp.ptr(), PyArray_FLOAT, 1, 1);
   if (y == NULL) {
-      Py_XDECREF(x);
-      throw Py::ValueError("y is of incorrect type (wanted 1D float)");
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
+     throw Py::ValueError("y is of incorrect type (wanted 1D float)");
   }
-  PyArrayObject *d = (PyArrayObject *) PyArray_ContiguousFromObject(dp.ptr(), PyArray_UBYTE, 3, 3);
+  d = (PyArrayObject *) PyArray_ContiguousFromObject(dp.ptr(), PyArray_UBYTE, 3, 3);
   if (d == NULL) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::ValueError("data is of incorrect type (wanted 3D UInt8)");
   }
   if (d->dimensions[2] != 4) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
-      Py_XDECREF(d);
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::ValueError("data must be in RGBA format");
   }
 
@@ -1206,26 +1400,21 @@
   int nx = x->dimensions[0];
   int ny = y->dimensions[0];
   if (nx != d->dimensions[1] || ny != d->dimensions[0]) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
-      Py_XDECREF(d);
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::ValueError("data and axis dimensions do not match");
   }
 
   // Allocate memory for pointer arrays
-  unsigned int * rowstarts = reinterpret_cast<unsigned int*>(PyMem_Malloc(sizeof(unsigned int)*rows));
-  if (rowstarts == NULL) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
-      Py_XDECREF(d);
+  rowstarts = reinterpret_cast<unsigned int*>(PyMem_Malloc(sizeof(unsigned int)*rows));
+  arows = reinterpret_cast<float *>(PyMem_Malloc(sizeof(float)*rows));
+  if (rowstarts == NULL || arows == NULL ) {
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::MemoryError("Cannot allocate memory for lookup table");
   }
-  unsigned int * colstarts = reinterpret_cast<unsigned int*>(PyMem_Malloc(sizeof(unsigned int*)*cols));
-  if (colstarts == NULL) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
-      Py_XDECREF(d);
-      PyMem_Free(rowstarts);
+  colstarts = reinterpret_cast<unsigned int*>(PyMem_Malloc(sizeof(unsigned int)*cols));
+  acols = reinterpret_cast<float*>(PyMem_Malloc(sizeof(float)*cols));
+  if (colstarts == NULL || acols == NULL) {
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::MemoryError("Cannot allocate memory for lookup table");
   }
 
@@ -1238,54 +1427,17 @@
   size_t NUMBYTES(rows * cols * 4);
   agg::int8u *buffer = new agg::int8u[NUMBYTES];
   if (buffer == NULL) {
-      Py_XDECREF(x);
-      Py_XDECREF(y);
-      Py_XDECREF(d);
-      PyMem_Free(rowstarts);
-      PyMem_Free(colstarts);
+     _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
       throw Py::MemoryError("Could not allocate memory for image");
   }
+  
 
   // Calculate the pointer arrays to map input x to output x
-  unsigned int i, j, j_last;
+  unsigned int i, j;
   unsigned int * colstart = colstarts;
   unsigned int * rowstart = rowstarts;
   float *xs1 = reinterpret_cast<float*>(x->data);
   float *ys1 = reinterpret_cast<float*>(y->data);
-  float *xs2 = xs1+1;
-  float *ys2 = ys1+1;
-  float *xl = xs1 + nx - 1;
-  float *yl = ys1 + ny - 1;
-  float xo = x_min + dx/2.0;
-  float yo = y_min + dy/2.0;
-  float xm = 0.5*(*xs1 + *xs2);
-  float ym = 0.5*(*ys1 + *ys2);
-  // x/cols
-  j = 0;
-  j_last = j;
-  for (i=0;i<cols;i++,xo+=dx,colstart++) {
-      while(xs2 != xl && xo > xm) {
-          xs1 = xs2;
-          xs2 = xs1+1;
-          xm = 0.5f*(*xs1 + *xs2);
-          j++;
-      }
-      *colstart = j - j_last;
-      j_last = j;
-  }
-  // y/rows
-  j = 0;
-  j_last = j;
-  for (i=0;i<rows;i++,yo+=dy,rowstart++) {
-      while(ys2 != yl && yo > ym) {
-          ys1 = ys2;
-          ys2 = ys1+1;
-          ym = 0.5f*(*ys1 + *ys2);
-          j++;
-      }
-      *rowstart = j - j_last;
-      j_last = j;
-  }
 
 
   // Copy data to output buffer
@@ -1297,87 +1449,82 @@
   agg::int8u * position = buffer;
   agg::int8u * oldposition = NULL;
   start = reinterpret_cast<unsigned char*>(d->data);
-  for(i=0;i<rows;i++,rowstart++)
-  {
-      if (i > 0 && *rowstart == 0) {
-          memcpy(position, oldposition, rowsize*sizeof(agg::int8u));
-          oldposition = position;
-          position += rowsize;
-      } else {
-          oldposition = position;
-          start += *rowstart * inrowsize;
-          inposition = start;
-          for(j=0,colstart=colstarts;j<cols;j++,position+=4,colstart++) {
-              inposition += *colstart * 4;
-              memcpy(position, inposition, 4*sizeof(agg::int8u));
+    int s0 = d->strides[0];
+    int s1 = d->strides[1];
+  
+  if(interpolation == Image::NEAREST) {
+      _bin_indices_middle(colstart, cols, xs1,  nx,dx,x_min);
+      _bin_indices_middle(rowstart, rows, ys1,  ny, dy,y_min);
+      for(i=0;i<rows;i++,rowstart++)
+      {
+          if (i > 0 && *rowstart == 0) {
+              memcpy(position, oldposition, rowsize*sizeof(agg::int8u));
+              oldposition = position;
+              position += rowsize;
+          } else {
+              oldposition = position;
+              start += *rowstart * inrowsize;
+              inposition = start;
+              for(j=0,colstart=colstarts;j<cols;j++,position+=4,colstart++) {
+                  inposition += *colstart * 4;
+                  memcpy(position, inposition, 4*sizeof(agg::int8u));
+              }
           }
       }
   }
+  else if(interpolation == Image::BILINEAR) {
+      _bin_indices_middle_linear(acols, colstart, cols, xs1,  nx,dx,x_min);
+      _bin_indices_middle_linear(arows, rowstart, rows, ys1,  ny, dy,y_min);
+      double a00,a01,a10,a11,alpha,beta;
+      
+      
+        agg::int8u * start00;
+        agg::int8u * start01;
+        agg::int8u * start10;
+        agg::int8u * start11;
+      // Copy data to output buffer
+      for (i=0; i<rows; i++)
+        {
+            for (j=0; j<cols; j++)
+            {
+                alpha=arows[i];
+                beta=acols[j];
+                
+                a00=alpha*beta;
+                a01=alpha*(1.0-beta);
+                a10=(1.0-alpha)*beta;
+                a11=1.0-a00-a01-a10;
+                
+                start00=(agg::int8u *)(start + s0*rowstart[i] + s1*colstart[j]);
+                start01=start00+s1;
+                start10=start00+s0;
+                start11=start10+s1;
+                position[0] =(agg::int8u)(start00[0]*a00+start01[0]*a01+start10[0]*a10+start11[0]*a11);
+                position[1] =(agg::int8u)(start00[1]*a00+start01[1]*a01+start10[1]*a10+start11[1]*a11);
+                position[2] =(agg::int8u)(start00[2]*a00+start01[2]*a01+start10[2]*a10+start11[2]*a11);
+                position[3] =(agg::int8u)(start00[3]*a00+start01[3]*a01+start10[3]*a10+start11[3]*a11);
+                position += 4;
+            }
+        }
+        
+    }
 
+
   // Attatch output buffer to output buffer
   imo->rbufOut = new agg::rendering_buffer;
   imo->bufferOut = buffer;
   imo->rbufOut->attach(imo->bufferOut, imo->colsOut, imo->rowsOut, imo->colsOut * imo->BPP);
 
-  Py_XDECREF(x);
-  Py_XDECREF(y);
-  Py_XDECREF(d);
-  PyMem_Free(rowstarts);
-  PyMem_Free(colstarts);
+  
+   _pcolor_cleanup(x,y,d,rowstarts,colstarts,acols,arows);
 
   return Py::asObject(imo);
-}
+    
 
-void _bin_indices(int *irows, int nrows, double *y, int ny,
-                                            double sc, double offs)
-{
-    int i;
-    if (sc*(y[ny-1] - y[0]) > 0)
-    {
-        int ii = 0;
-        int iilast = ny-1;
-        int iy0 = (int)floor(sc * (y[ii]  - offs));
-        int iy1 = (int)floor(sc * (y[ii+1]  - offs));
-        for (i=0; i<nrows && i<iy0; i++) {
-            irows[i] = -1;
-        }
-        for (; i<nrows; i++) {
-            while (i > iy1 && ii < iilast) {
-                ii++;
-                iy0 = iy1;
-                iy1 = (int)floor(sc * (y[ii+1] - offs));
-            }
-            if (i >= iy0 && i <= iy1) irows[i] = ii;
-            else break;
-        }
-        for (; i<nrows; i++) {
-            irows[i] = -1;
-        }
-    }
-    else
-    {
-        int iilast = ny-1;
-        int ii = iilast;
-        int iy0 = (int)floor(sc * (y[ii]  - offs));
-        int iy1 = (int)floor(sc * (y[ii-1]  - offs));
-        for (i=0; i<nrows && i<iy0; i++) {
-            irows[i] = -1;
-        }
-        for (; i<nrows; i++) {
-            while (i > iy1 && ii > 1) {
-                ii--;
-                iy0 = iy1;
-                iy1 = (int)floor(sc * (y[ii-1] - offs));
-            }
-            if (i >= iy0 && i <= iy1) irows[i] = ii-1;
-            else break;
-        }
-        for (; i<nrows; i++) {
-            irows[i] = -1;
-        }
-    }
+  
 }
 
+
 char __image_module_pcolor2__doc__[] =
 "pcolor2(x, y, data, rows, cols, bounds, bg)\n"
 "\n"
@@ -1477,7 +1624,7 @@
         Py_XDECREF(bg);
         throw Py::MemoryError("Cannot allocate memory for lookup table");
     }
-    int * jcols = reinterpret_cast<int*>(PyMem_Malloc(sizeof(int*)*cols));
+    int * jcols = reinterpret_cast<int*>(PyMem_Malloc(sizeof(int)*cols));
     if (jcols == NULL) {
         Py_XDECREF(x);
         Py_XDECREF(y);
-------------------------------------------------------------------------
This SF.Net email is sponsored by the Moblin Your Move Developer's challenge
Build the coolest Linux based applications with Moblin SDK & win great prizes
Grand prize is a trip for two to an Open Source event anywhere in the world
http://moblin-contest.org/redirect.php?banner_id=100&url=/
_______________________________________________
Matplotlib-devel mailing list
Matplotlib-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-devel

Reply via email to