Code covered by the BSD License  

Highlights from
Surface Fitting using gridfit

4.87879

4.9 | 107 ratings Rate this file 500 Downloads (last 30 days) File Size: 1.38 MB File ID: #8998
image thumbnail

Surface Fitting using gridfit

by

 

11 Nov 2005 (Updated )

Model 2-d surfaces from scattered data

Editor's Notes:

This file was a File Exchange Select file.

Select files are submissions that have been peer-reviewed and approved as meeting a high standard of utility and quality.

| Watch this File

File Information
Description

Those wishing to model a surface from data in the
form of z(x,y) from scattered or semi-scattered
data have had few options in matlab - mainly
griddata.

Griddata is a valuable tool for interpolation of
scattered data. However it fails when there are
replicates or when the data has many collinear
points. Griddata is also unable to extrapolate
beyond the convex hull of the data unless the 'v4'
option is used, which is slow.

Gridfit solves all of these problems, although it
is not an interpolant. It builds a surface over a
complete lattice, extrapolating smoothly into the
corners. You have control of the amount of
smoothing done, as well as interpolation methods,
which solver to use, etc.

This release allows the user to solve much larger problems using a new tiling option. There is essentially no limit on the size of the suface one builds now, as long as you have dense enough data and enough memory to store the final gridded surface.

Example uses are found in the file gridfit_demo.m,
as well as comparisons to griddata on the same
surfaces.

Acknowledgements

This file inspired Bi Dimensional Emperical Mode Decomposition (Bemd) and Regularize Data3 D.

MATLAB release MATLAB 7.0.1 (R14SP1)
Other requirements Gridfit should be accessible to older versions of matlab, as long as symamd and colamd are available. Older releases could use one of the iterative solvers, such as symmlq.
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (178)
12 Nov 2014 Judith

This is a fantastic function and it runs perfectly with my 2D-data. Do you also provide a gridfit3d version? That whould save a lot of time in my post-processing work. Thanks in advance.

10 Nov 2014 John D'Errico

http://blogs.mathworks.com/community/2010/12/13/citing-file-exchange-submissions/

10 Nov 2014 Chaipichit

How to cited your code in research paper?

24 Oct 2014 John D'Errico

Evgheny - If you seriously need to fit such a surface in spherical coordinates, I might be able to help you out with a separate, un-posted toolbox that is capable of fitting a surface in spherical coordinates. You would need to contact me directly.

23 Oct 2014 Evgheny

Thanks, function does do its job great!

Does anyone know alternative of this function for closed surface (sphere-like).
If I try this function with spherical coordinates, it works fine, but there are problems with the seam (line of joint)

30 Sep 2014 John D'Errico

Sorry. I've never heard of a Java implementation. Not that it counts, I've written it 7 times in MATLAB, once each in Fortran and APL, those last were very old though. The difference each time was what I learned from the previous incarnation. Were someone to do a Java implementation perhaps the most important thing is to use sparse linear algebra capabilities for the solve, else it would take forever.

30 Sep 2014 behruz

is there any implementation of this in java? that get some points location and value and give a grdi?

24 Sep 2014 Alex

Just to add my $0.02 to this whole sphere discussion. By math definition of what is function it means for each x we have only one F(x). If there are more than one F(x) relation between x and F(x) is not a function.

20 Sep 2014 Jacob

I wish I could remove my previous question, 2 seconds on google pointed me to Matlab's interp2 which I'm using now on the output of your excellent gridfit. Thanks for sharing!

20 Sep 2014 Jacob

I found that it is very straightforward to fit a nice looking surface to some 55 data-points using gridfit. I was wondering though, is there a way to get values of the surface that are not on the nodes? I wish to interpolate between my 55 samples using the surface.

20 Sep 2014 Jacob  
30 Aug 2014 John D'Errico

For those of you who don't appreciate why Felon's comment is silly, think of it like this. Gridfit fits a surface of the form f(x,y), over a rectangular grid. It does so quite well, as many people have found over the years.

While you may think of the surface of a sphere as a surface, it is not of the form that gridfit can fit. It is multi-valued, so for any single (x,y) pair, there will be zero, one, or two values of z that would apply. As well, that "surface" (better to call it a manifold) has derivative singularities, if we were to look at it as a function of x and y. So even a hemisphere will be problematic for this tool.

You would not use gridfit to fit something that is not representable as a function of two variables over a rectangular grid, any more than you would expect it to do numerical integration, numerical optimization, or compute an FFT. Nor would you expect it to cook dinner for you, do your laundry, etc. Use the right tool to solve your problem, but if you try to force the wrong tool to solve a random problem, expect poor results and don't complain about what you get.

29 Aug 2014 John D'Errico

Felon -

So you are using a tool that builds a single valued function to fit something that is obviously not. What did you expect? Magic?

Software does what it is programmed to do. It does not magically rewrite itself when you give it a problem of a completely different sort. In fact, I fail to understand why you would downrate a tool for not solving a class of problem it is explicitly not designed to solve.

If you have a closed manifold, like a ball or some other multivalued form, then don't use this tool. I have NEVER claimed it would solve that problem. Instead, you might look into tools like convex hulls, alpha shapes, CRUST, etc. Or, you might choose to convert the problem into spherical coordinates, at which point gridfit would be able to build a viable surface.

Or maybe you just wanted to complain with no good reason.

29 Aug 2014 Felon Highman

Still it cannot surface a cluster of points that looks like a ball. The result the interpolation is always a surface of single-value function.

29 Aug 2014 Felon Highman  
30 Jul 2014 Angelo Tafuni  
26 Jul 2014 John D'Errico

No, gridfit does not explicitly allow you to apply derivative constraints. That does not say it is impossible, only that I did not offer it as an option.

The main reason why not, is it would require a set of linear inequality constraints on the unknowns. For a not uncommon grid of size 100x100, there are 100*100=10000 unknowns to solve for. This is not a problem, since the linear system is a sparse one. However, to solve a sparse linear inequality constrained system, one would need to use LSQLIN, or a solver like it. And the last time I checked, LSQLIN was not set up yet to handle sparse large scale inequality constrained problems. (That may have changed with the most recent release, but I have not checked.) If I made all of the matrices full ones, the solve time would probably be incredibly slow and memory intensive.

So I'm sorry, but gridfit will not handle the problem as is.

If you were willing to build a fairly coarse grid, AND add the constraint system, it would probably be doable in a reasonable time. I don't know how small the grid would need to be to make the solve time reasonable. And your definition of reasonable would surely differ from mine, depending on how badly you needed the answer.

26 Jul 2014 Javier

Hey John, would gridfit allow me to constraint the slope of the fitted surface for a given range of data?

I'm needing to model a surface in the form of z(x,y) from scattered data points while keeping the first derivative of z negative in both directions.

25 Feb 2014 Kim

I'm trying to make a surface approximation to describe the required input torque to drive a hydraulic motor operating under certain load conditions. My inputs are: output flow, output pressure and input shaft speed. I have scattered measurement data where the relationship between these variables are to be found. I read in some of the earlier posts that gridfit could be extended to higher-dimension fits. Do you have any version where a 3rd dimension can be included that I could try?

12 Oct 2013 Jamal

Gridfit is very nice. I have written an update that adds bicubic interpolation and improves the internal scaling (makes surface smoothness completely independent of the x/y grid size). If I send you the code, would you update Gridfit?

25 Jul 2013 Jon Griffiths  
17 Jun 2013 John D'Errico

Alejandro - while I would like to implement that capability one day in the future, that day may be a long time away before I find the time. Until then, all I can offer is to use a finer grid, which will mimic a higher order interpolant.

17 Jun 2013 John D'Errico

Dominik - Sorry, but no. You could use this tool to generate a gridded surface, then use interp2, to interpolate, but that is perhaps too much if you only wanted to interpolate a few scattered points.

17 Jun 2013 dominik

Hey
i wonder, can this code also be used like an interpolant (e.g. gridded interpolant), i.e. such that the output is not the smooth funtion on a grid, but on certain points queried?

24 May 2013 Alejandro Arrizabalaga

Great function!
I was just wondering: would it be difficult to extend it for a cubic interpolation? (it was mentioned in the code as a future enhancement).

11 Feb 2013 Roman Voronov  
11 Feb 2013 Roman Voronov

Exactly what I was looking for, the help could start with a simpler example:
[X,Y] = meshgrid(1:5);
Zsmooth = gridfit(X,Y,A,1:5,1:5,'smooth',10)
Zsmooth =
0.9638 1.9517 2.9512 3.9517 4.9638
1.9517 2.8162 3.7338 4.8162 5.9517
2.9512 3.7338 4.5266 5.7338 6.9512
3.9517 4.8162 5.7338 6.8162 7.9517
4.9638 5.9517 6.9512 7.9517 8.9638

06 Feb 2013 John D'Errico

Andrew - It sounds like your goal is to reproduce the local shape of your surface over a region that lies inside the bounds of your data. Gridfit can do this, if you are willing to expand that grid so that it fully contains the data. Simply add some nodes that extend to the boundaries of your data in both x and y. The nodes used by gridfit need not be equally spaced, so this should not be a problem, but to retain the essential shape of the surface, you may need to add a few such nodes in each direction. Once the fit is obtained, then you can always extract only that portion of the surface that you wish to retain.

If this is not what you are looking for, then feel free to contact me directly, perhaps with some data so that I can better understand your problem.

05 Feb 2013 Andrew

This is an amazing program. I do have one question. I have a scattered grid of latitude(Y) and longitude(X) coordinates corresponding to brightness(Z) values. I've been using GRIDDATA to interpolate to a new grid that is not scattered and within a smaller range of latitude and longitude values, but GRIDDATA smooths the data so much that it becomes useless.

I understand that GRIDFIT is not an interpolant, and from my reading of the help, it does not seem that it was meant to be used in this fashion, but I just want to make sure.

21 Nov 2012 John D'Errico

Collin - What you wish is essentially a 3-d version of gridfit, that could use time as the third dimension. Sorry, but perhaps one day...

21 Nov 2012 John D'Errico

Atul Ingle - Yes, it is possible to extend gridfit to higher dimensions. In practice it tends to turn into a memory & cpu hog before long. 5 or 6 dimensions were the practical extreme when I've done it before, although computers are bigger and faster since those days. 64 bit MATLAB and lots of RAM will help here, but even with that, the curse of dimensionality is huge and arrives quickly. Sadly, the n-d version I wrote many years ago is not mine to give away.

21 Nov 2012 Atul Ingle

Is it possible to extend this idea to fit 3d volume to scattered data, and in general, to higher dimensions? The only limitation I can see is matrix size.

16 Nov 2012 Collin

Excellent bit of code that has made my life a lot easier.

If you have knowledge of the evolution of your surface in time, how would you suggest incorporating knowledge of a previous state to help improve interpolation/extrapolation of the current state?

Thanks

Collin

08 Nov 2012 John D'Errico

Laura, I do have a tool (based on the same philosophy as gridfit) to handle fitting over a non-rectangular domain. That set of tools can also test for points inside a general domain, although inpolygon can handle that problem in 2-dimensions. I've not posted the code although it is fully functional, because I've never been sufficiently happy with the entire toolset. feel free to send me e-mail if you want to try it though.

08 Nov 2012 Laura

I recently (today) submitted a question entitled "availble solutions for 2D interpolation on non monotonic scattered data". Right after I discovered gridfit.
Thanks a lot for such a routine, really solid and with the best fitting I was looking for!.

However I would like to ask one thing that partially has already been addressed. Indeed I am most interested in the good interpolation it performs (although not the objective) than the extrapolation. My x y data points do not define a square domain but more a distorted "trapezoidal" one with curved contours. I would like to get rid of the data out of the original domain.
Ignacio proposed an interesting way (maybe not very clean but working) to know which points fall inside or outside the domain. But I can not apply it to my case since due to these convex borders the Delaunay triangulation produces spurious triangles there that indeed are out of the domain.

Could you suggest me in which way I could identify which of the points of my "square box grid" fall actually inside the trapezoidal domain?

Thanks a lot

31 Oct 2012 John D'Errico

Bruno - There are a couple of issues to watch out for here.

First of all, if you wish to find the area that is strictly above zero and below the surface, be careful if the surface goes below zero. (Some individuals miss this point, but it is crucial.) You could simply compute the area by calling trapz(trapz(max(0,Z))). This would be an approximation of course, missing out on the subtle issue of exactly where the surface passes below zero.

Next, depending on how the surface is intended to be interpolated, if you want the EXACT integral of that volume, then be careful. For the first order tensor product surface, i.e., that which interp2 would call 'linear', or what is called 'bilinear' in gridfit, then it suffices to interpolate at the center of each cell in the lattice of points. Sum up those values and you will have an exact volume. (Be careful to scale this result by the width of the cells in each dimension.) As pointed out, trapz (called twice, so once in each dimension) will also be exact here. It turns out that merely computing the average value over each of these cells is also equivalent to that act, then sum up the volumes in each cell of the lattice.

For the triangular faceted surface, trapz is not exact, but the difference may be subtle. Here we end up with a weighted linear combination of the values of each corner of the square cell in the lattice, but the 4 corners are not weighted equally. Send me an e-mail and I can detail how to do this more clearly.

Really though, the easy answer is to use trapz, but you should not ignore these subtleties I've pointed out.

31 Oct 2012 IainIsm

John - thank you, 'springs' seems to be doing the job!

Bruno - you can use trapz twice: http://blogs.mathworks.com/loren/2011/06/13/calculating-the-area-under-a-surface/

31 Oct 2012 Bruno

Excellent work. I have one question... how can I calculate the volume beneath the interpolated surface and above the plane xy? Maybe it's a simple task but i can't figure out how it's done.

31 Oct 2012 John D'Errico

IainIsm - You cannot (easily) place bound constraints on this fit. That would result in a very large scale bound constrained fit. While doable in theory with a quadratic programming tool or lsqlin, these problems are huge. It would simply take too long to solve.

Having said that, there is a possible solution that might work. Choose the 'springs' regularizer as an option. That option was designed to avoid extrapolation wherever possible, so if your data never goes negative, then the surface will hopefully not do so either.

John

30 Oct 2012 IainIsm

Lovely function - just one question, is it possible to put constraints on the range of values that zgrid can take please? for example, I'm using this to create a surface for a motor rpm/torque/efficiency map (which I am then querying with interp2 to provide an efficiency estimate for a given speed and loading), but every now and then I find that some of the zgrid values are negative. I could set these to zero using find, but I'm hoping to find(!) something a little more elegant!

25 Oct 2012 Chenjie  
27 Aug 2012 Rijk  
12 Jul 2012 Jesse  
25 Jun 2012 Jordan Malof

Great stuff! Additionally, and refreshingly, it worked "right out of the box". Thanks John!

07 Jun 2012 bouncing  
18 May 2012 Josh

John,

Thanks a lot. Not only did I enjoy an excellent bit of code, but also for your variety of stimulating and entertaining responses. Thank you for breaking up my day with at least 5 minutes of enjoyable rhetoric and exceptional references (Mark Twain / god and little green apples) ... brilliant! Have a few extra stars.

Josh

17 May 2012 K E

Lovely, thanks

06 May 2012 Rodrigo R. Oliveira

It is wonderfull. Easy solve my problem!

17 Apr 2012 Gonzalo  
16 Apr 2012 Peter

Oops. In my explanation above P = 1/sqrt(2), not sqrt(2).

13 Apr 2012 Peter

Lanis, John, et al.

This is how to get consistent behavior with different grid sizes.

--------

Save a copy of gridfit.m to gridfit2.m.

In gridfit2.m make the following changes at the locations of the commented-out code (near Lines 550-600):

% Sorry, John D.
% dy1 = dy(j(:)-1)/params.yscale;
% dy2 = dy(j(:))/params.yscale;
dy1 = dy(j(:)-1);
dy2 = dy(j(:));


% Sorry, John D.
% Append the regularizer to the interpolation equations,
% scaling the problem first. Use the 1-norm for speed.
%NA = norm(A,1);
%NR = norm(Areg,1);
%A = [A;Areg*(smoothparam*NA/NR)];
FidelityEquationCount = size(A, 1);
RegularizerEquationCount = nx * (ny - 2) + ny * (nx - 2);
NewSmoothParam = smoothparam * sqrt(FidelityEquationCount / RegularizerEquationCount);
A = [A; Areg * NewSmoothParam];

--------

Run the following script to demonstrate:

x = [0, 1.333333333, 2.666666667, 4, 0, 0, 4, 4, 2];
y = [2, 2, 2.1, 1.8114, 0, 4, 0, 4, 0];
z = [0.5, 0.82, 0.4, 0.7, 0.5, 0.7, 0, 0.7, 0];

% Set the smoothness for the original GridFit. This smoothness produces inconsistent results if the grid size changes.
Smoothness = 10;

% Set up the grids with different spacings...
x_coarse = linspace(0, 4, 11);
y_coarse = x_coarse;
x_fine = linspace(0, 4, 41);
y_fine = x_fine;

% We're taking a slice through each surface at the same coordinate.
% They have different indexes but the same coordinate of 1.2.
Index_coarse = find(y_coarse == 1.2);
Index_fine = find(y_fine == 1.2);

% Run GridFit.
z_coarse = gridfit(x, y, z, x_coarse, y_coarse, 'Smoothness', Smoothness);
z_fine = gridfit(x, y, z, x_fine, y_fine, 'Smoothness', Smoothness);

% Try the same thing with the updated GridFit. Note that the smoothness parameter is on a different scale compared to
% the original GridFit. However, all choices of smoothness value should produce consistent results for any grid size.
UpdatedGridFitSmoothness = 1;
% Run the updated version of GridFit, GridFit2.
z_coarse2 = gridfit2(x, y, z, x_coarse, y_coarse, 'interp', 'bilinear', 'Smoothness', UpdatedGridFitSmoothness);
z_fine2 = gridfit2(x, y, z, x_fine, y_fine, 'interp', 'bilinear', 'Smoothness', UpdatedGridFitSmoothness);

% Show the results for GridFit...
% These lines only line up when you use a very large or very small smoothness value.
subplot(1, 2, 1)
hold on
plot(x_coarse, z_coarse(Index_coarse, :), 'color', 'blue');
plot(x_fine, z_fine(Index_fine, :), 'color', 'red')
title('Profiles from the Original GridFit');
axis([1, 4, 0.3, 0.6]);
hold off

% Show the results for GridFit2...
% These lines have the same profile for all smoothness values.
subplot(1, 2, 2)
hold on
plot(x_coarse, z_coarse2(Index_coarse, :), 'color', 'blue');
plot(x_fine, z_fine2(Index_fine, :), 'color', 'red')
title('Profiles from the Updated GridFit');
axis([1, 4, 0.3, 0.6]);
hold off

13 Apr 2012 Peter

Lanis, John, et al.

How it works:

GridFit is a function that regularizes a dataset. It constructs a series of linear equations that are the "fidelity equations" (how closely the output matches your input data) and the "smoothness equations" (how close the surface's second derivatives are to zero). It assembles these equations and solves them in a way that minimizes the sum of the squared error (SSE) of the output surface.

The key is to recognize that the SSE is just a sum of squared residuals, the fidelity squared residuals (SSE_F) and the smoothness squared residuals (SSE_Sm).
SSE = SSE_F + SSE_Sm

If we have a certain number of input points, that tells us how many terms are in SSE_F. The number of grid points (more precisely, the number of second derivatives within the grid) tells us how many terms are in SSE_Sm.

Imagine that we are trying to fit this parabola: y = x^2 to our data. The fidelity terms are some finite value (which doesn't matter for this example), such as 100. Because we're using this parabola, the second derivative is always 2. If we have a grid such that there are 20 smoothness terms, then SSE_Sm = 20 * (2^2) = 80.

The ratio SSE_F / SSE_Sm = 100 / 80 = 1.25. This ratio determines how to weigh the smoothness against the goodness of fit.

Now suppose that we suddenly double the number of grid points. Now SSE_Sm = 2 * 20 * (2^2) = 160. That means SSE_F / SSE_Sm = 100 / 160 = 0.625. Uh-Oh! The ratio has changed, which means our fit will favor smoothness more than it did before!

The way to fix this is to reduce the values of the second derivatives by multiplying them by a constant, which we'll call P. In this example, if P = sqrt(2), then the individual terms in SSE_Sm will have half of the value they had before. That's ok because there are twice as many of them. In other words, 0.5 * 2 = 1.

Now using P, SSE_Sm = 2 * 20 * ((P * 2) ^ 2) = 80. This provides the same balance of fidelity vs. smoothness that we originally had.

What I did in GridFit2 is a little more complicated because it also considers the smoothness parameter and the number of fidelity equations. I hope this helps!

12 Apr 2012 John D'Errico

Peter - if you can do better, feel free to write your own code. (It is NOT the use of the 1-norm that impacts your question.)

12 Apr 2012 Peter

I'm not convinced that the smoothness parameter produces consistent results for different grid sizes. I think the problem comes from the use of the 1-norm instead of just using the numbers of fidelity equations vs. derivatives.

Here is an example.

x = [0, 1.333333333, 2.666666667, 4, 0, 0, 4, 4, 2];
y = [2, 2, 2.1, 1.8114, 0, 4, 0, 4, 0];
z = [0.5, 0.82, 0.4, 0.7, 0.5, 0.7, 0, 0.7, 0];

Smoothness = 10;

x_coarse = linspace(0, 4, 11);
y_coarse = x_coarse;
Index_coarse = find(y_coarse == 1.2);

x_fine = linspace(0, 4, 41);
y_fine = x_fine;
Index_fine = find(y_fine == 1.2);

z_coarse = gridfit(x, y, z, x_coarse, y_coarse, 'Smoothness', Smoothness);
z_fine = gridfit(x, y, z, x_fine, y_fine, 'Smoothness', Smoothness);

hold on
plot(x_coarse, z_coarse(Index_coarse, :), 'color', 'blue');
plot(x_fine, z_fine(Index_fine, :), 'color', 'red')
hold off

04 Apr 2012 Daniel

yes, Hendawi Mohamed I removed points below the edge of my data. My data is roughly bounded by a triangle in X-Y so I used:

gf=gridfit(....);
Z=griddata(... same limits);
N=find(isnan(Z));

gf(N)=NaN;

to removed the complete estimated results. When I plotted my data, the points on the graph corresponded to area were I had data.

02 Apr 2012 John D'Errico

You would normally first decide what size result you wish to see, thus the number of nodes in x and y. Then you would choose some amount of smoothing.

I've set it up so that as well as I could so that changing the grid size for a fixed smoothness will have a minimal impact on the overall look of the surface. For example, try this example:

xy = rand(100,2);
z = rand(100,1);
zgrid50 = gridfit(xy(:,1),xy(:,2),z,50,50);
zgrid200 = gridfit(xy(:,1),xy(:,2),z,200,200);
subplot(1,2,1)
surf(zgrid50)
shading interp
subplot(1,2,2)
surf(zgrid200)
shading interp

The two surfaces will look qualitatively similar, although the second one will look less jagged because it is 4 times finer in each dimension. As expected, this is what you gain by using a finer grid, at some cost in time for the solve. In this sense, smoothness is fully decoupled (as well as I could do that) from the size of the grid itself, for any given set of data.

Given a choice of grid, a smoothness of 1 means that there is an equal amount of importance associated with reducing the residuals compared to the overall smoothness of the surface. Be careful however, in making the smoothness too large or too small, as then you may find numerical artifacts generated from the solution. Floating point arithmetic only goes so far.

By the way, a good rule of thumb is if you cannot get the result you wish with a smoothness somewhere between 0.001 and 1000, then I would look to see if there may be a fundamental problem in the data.

As well, if one chose to fit a 2x2 surface through 1000 data points, you will never get a perfect fit regardless of how small a smoothing parameter you choose, UNLESS the data happens to perfectly fit the implicit model used to interpolate inside a unit square. (For example, tensor product linear interpolation uses an implicit local model of a0 + a1*x + a2*y + a3*x*y within any cell.) There simply are not sufficient degrees of freedom in the model for a very coarse grid if you have many data points.

02 Apr 2012 Ianis

Hi! Thank you for this useful function ! I needed to fit scattered points with noise added, and Gridfit handles that nicely ! Now, I would like to test with different values of smoothness and number of x-y- nodes in order to find the best couple of parameters for my problem. So, I have some questions about the 'smoothness' parameter. You explain in the help that the smoothness parameter is ''normalized'' so a value of '1' may generate reasonable results. I would like to know with respect to which parameter do you normalize the ''smoothness'' value ? If I use either 'gradient ' or ' laplacian' as regularizer, is the smoothness defined in the same way ?
Thanks.

02 Apr 2012 Ianis  
31 Mar 2012 John D'Errico

Gridfit is not a tool to do 2-d interpolation. Use interp2 instead.

30 Mar 2012 MAN CHIA

How could I use gridfit to interpolate a set of 2D data on a image?
Just let all z equals 0 ?

30 Mar 2012 MAN CHIA  
19 Mar 2012 R

John,

All your submissions I have tried are great quality! Thank you very much. I have a directory for your functions, they could easily be part of Matlab's distribution in my opinion.
Re gridfit, do you still have in mind submitting a 3D/ND version of it?

Rodrigo.

25 Feb 2012 Katie  
19 Dec 2011 Tomas

Great work!

17 Nov 2011 John D'Errico

What can I say? You can have a good result, or a faster one. Pick only one.

The fact is, gridfit is actually blazingly fast given what it does, which in this case requires the solution of a linear system with 90,000 unknowns. A quick test shows that took only about 2 seconds on my system. Perhaps you want me to do magic? Even a wave of a magic wand takes about 2 seconds.

17 Nov 2011 Yan

Hi, this is an awesome program! But it seems to be very memory consuming. I can solve for Zgrid points with a rather small dimension only (e.g. 300x300). Griddata seems to be able solver for a larger matrix faster, but the interpolation is not nearly as good. Is there anyways to obtain higher resolution?

12 Nov 2011 John D'Errico

By default, it uses triangular interpolation. If you read the help, you will see that a bilinear (quad based) interpolation method is an option. Is one better than another? In general, I'd claim that if you can see the difference, then your grid is far too coarse! In any event, HAD you actually read the help before posting, in your example you might have tried adding one more option to the call.

[zGrid, xGrid, yGrid] = gridfit(SourceData(:, 1), SourceData(:, 2), SourceData(:, 3), xNodes, yNodes, 'smoothness', 0.1,'interp','bi');

Your second comment is about the smoothing parameter. Note that gridfit allows a common scheme for the smoothness penalty function - that the Laplacian is biased towards zero. This is the sum of the second partials, and it is arguably the logical choice for SOME physical surfaces. HOWEVER, note that the default is NOT that method. In fact, the default is a method that DOES uncouple those second partials!!!!!!!!!!!!!!!! Again, read the help, rather than assuming you know the answer and then asking a question based on that wrong assumption.

As for the optimal smoothness parameter changing based on the grid spacing, I have found that for virtually any method of choosing a smoothness parameter, I can also come up with a case where it will not be the best. Gridfit uses a default smoothing parameter that is "reasonable" for many problems out of the box. Is it optimal? Probably not. Any degree of true optimality might involve algorithms that would be slower to run, and algorithmic speed is desirable here. Should there be an adaptive option in gridfit, that would allow the user to set it and walk away, knowing that it will be slower to run, but it will give an always optimal result, for all users? Truly good adaptive methods that will never fail are also terribly difficult to write. Look at the numerical integration tools in MATLAB (quad, quadl, quadgk, etc.) I can easily make those tools fail by a careful choice of integrand and interval.

All that I can suggest is to use the facilities built into gridfit as a computational tool. If it produces by default not exactly what you want, then I have provided a few knobs you can adjust. Only you know if you like the surface you create with a tool like gridfit. But at least the knobs are there to turn to generally produce something that will make you happy.

10 Nov 2011 Peter

GridFit is pretty good, but I have a few questions.

1) The interpolation seems to be triangular, which causes some symmetry problems. Below is an example program to illustrate this. Does anybody know if there is an advantage to using triangular interpolation (three points) instead of rectangular interpolation (four points)?

2) If I change the grid spacing, I also have to adjust the smoothness parameter or GridFit gives different results. This appears to be related to the way GridFit sets up its derivative equations. It mixes horizontal and vertical second derivatives in the same equation. For example, if the residual in the horizontal direction is +0.1 and the residual in the vertical direction is -0.1, the regularizer wrongly deems that point a "perfect fit" because the derivatives are added together. Would it be better to separate the derivatives into separate equations (separate rows in Areg)? If so, is it possible to produce the same results (same smoothness or curve profile) for a single "smoothness" parameter, no matter what grid spacing is used? I think that would be really helpful and convenient to GridFit users so that we wouldn't have to fiddle with the smoothness number manually.

% Choose four points at each of four corners and one point in the center.
% Make sure everything is symmetric around the center.
SourceData = [
0, 0, 0
0, 1, 0
1, 0, 0
1, 1, 0
0.5, 0.5, 1
];

% Create a 4x4 grid.
xNodes = linspace(0, 1, 4);
yNodes = xNodes;

% Run GridFit. The smoothness parameter has little effect on the result
% for this demonstration.
[zGrid, xGrid, yGrid] = gridfit(SourceData(:, 1), SourceData(:, 2), SourceData(:, 3), xNodes, yNodes, 'smoothness', 0.1);

% Display the surface so that we can see the asymmetry.
surf(xGrid, yGrid, zGrid);
xlabel('x');
ylabel('y');
zlabel('z');

% Calculate the asymmetry.
disp(['GridFit asymmetry is ', num2str(zGrid(2, 2) - zGrid(3, 2)), '.']);

12 Oct 2011 John D'Errico

This will definitely have much better performance around the edges compared to a delaunay based tool. They often show serious interpolation artifacts around the edges.

In terms of how good is the fit, I'm surprised that I did not return an array of residuals or predicted values, but this is easy enough to build. Use interp2 to compute predicted values, since what comes out of gridfit is now a nice regular one that interp2 will handle. Once you have the predictions, residuals are easy enough to get by subtraction, or you could compute a standard error, etc.

12 Oct 2011 Mark

Very nice, currently using this to perform a calibration between an imager and a gimbaled laser pointer. It's not quite optimized to get the best expected performance over the full FOV, but it's close. Toward the edges, it works much better than other routines I've tried. Is there a way to see how good the fit is?

12 Oct 2011 Mark  
02 Sep 2011 Darin

Thanks for the great tool. Trying to find a way around "pull ins" when the regularizer has a significant low spatial frequency bias over extended areas... for example, around the top of a higly resolved gaussian peak with 'gradient'. Using a higher order regularizer doesn't help (I've implemented a 4th order gradient): that just gives a surface that's highly "puckered" around indiviual sample errors. Any ideas?

27 Aug 2011 Christian

I'm sorry - that was a typo: I meant gridfit when I typed griddata. I since figured it out it was a MATLAB environment issue, unrelated to gridfit. Contour plots just weren't showing up in the display. Restart solved the problem. I would edit or remove my original post, but this forum doesn't allow for it. Thanks for the quick reply.

27 Aug 2011 John D'Errico

How can anybody guess what you are doing wrong? How exactly does it fail? What are you doing exactly? What is your data? And why are you saying that griddata output is a problem, in a question about gridfit?

26 Aug 2011 Christian

I'm trying to use gridfit to produce a contour map (functions: contour and contourf). MATLAB seems to interpret the griddata output okay for the surf() function, but fails with contour(). Anyone else run into this problem? I think it has something to do with grid orientation.

01 Aug 2011 M. C Ertem

Thank you!!! Saved the day. (Actually, night.)

27 Jul 2011 Yang

Hi John,
I tried some smaller parameters (say 0.1 for smoothness), but then I get wave-like artifact pattern in my image.
What I am looking for is a combination: griddata in the inner part (so no smearing at all), and gridfit at the edge (where griddata produces NaN). Hence smearing is minimized. Can I achieve that by tuning parameters? Thank you very much!

22 Jul 2011 John D'Errico

Smaller values of the smoothing parameter give less smoothing, although zero as a value may result in a singularity.

22 Jul 2011 Yang

Hi John, I have a set of 2d scattering data, and I want to have interpolation values at a regular (1:N,1:N) lattice. What parameters should I use to minimize the smoothing? (I need to iterate the interpolation for 50+ times, setting "smoothness"=1 for 50 times makes the resulting image very blur).
Thanks!

21 Jul 2011 Qiang

hi John,
I got error information now, but not last year. Seemingly it is because the code is not consistent with the new Matlab version (2011a). Maybe I am wrong. Do you have ideas about this? The error is as following:

??? Error using ==> mldivide
MLDIVIDE is not supported for one sparse input and one single input.

Error in ==> gridfit at 616
zgrid = reshape(A\rhs,ny,nx);

Thanks

03 Jul 2011 stefano pasquali

hi John.
i'm pretty new in matlab. so maybe i'm going to do a stupid question.

i have my data already in a matrix form.
tipically 200x200.

i cannot use your function because it expect 3 vectors and nodes ...

can you help me ?

basically i have this surface, results of observation (along x and y there are categorical variables) and i'm expecting a monotonic decreasing surface ... i got it but with some outliers ... i think your function can help me.

any suggestion ?

26 May 2011 John D'Errico

Melissa - gridfit merely generates a function defined at a set of discrete node points, on a regular lattice. It does not deliver an interpolated prediction at any point. But you can easily enough get that prediction from the surface produced from gridfit.

There are three interpolation methods that are allowed in gridfit. Those methods define how gridfit treated your data when it built the surface itself. Nothing requires that you use the same method to interpolate that gridfit used when it built the surface though.

The nice thing is, we already have interp2 to do 2-dimensional interpolation once the surface is defined on a regular lattice, with several good methods already there. The case of 'linear' for interp2 is equivalent to that employed in gridfit for the 'bilinear' case. It turns out that this is also a method used by tools like photoshop for image interpolation.

You can also use the other methods of interp2 though - splines for example. Or, if you prefer a simplicial interpolant, where each square of the lattice is split into a pair of triangles, then a truly linear interpolation is done across those triangles, you can find my interpns on the file exchange.

http://www.mathworks.com/matlabcentral/fileexchange/30932

So yes, you can easily do interpolation at any arbitrary point that lies inside the lattice as defined by gridfit.

26 May 2011 Melissa

This is an amazing piece of work John. I bow to your superior programming skills and knowledge. I have a question though, is it possible for me to evaluate a point or a series of points on the surface?

16 May 2011 John Reinert Nash

Excellent tool that much improved my visualization of some photographic exposure data. Thanks, John, great to use one of your tools (again).

16 May 2011 John Reinert Nash  
12 Apr 2011 Januka

Great piece of work John. No hassle and easy to use. Thanks!

06 Apr 2011 Ignacio

Really good

One way of using gridfit just for interpolation is :
- Use gridfit in the whole grid zifit=gridfit(X,Y,Z,xi,yi)
- Use griddata in the same grid
zidata=griddata(X,Y,Z,xi,yi)
- Use isnan to know which elements are NaN in zidata, and apply it to zifit

At least for me it worked
- Use a isnan(griddata) to know which points are inside the original data, and apply it to the points created with gridfit

31 Mar 2011 jose tapia

Thank's for your answer John D'Errico.

26 Mar 2011 John D'Errico

You have asked this question in multiple places. I have no idea what simulink does, so I cannot answer your question, nor can I even be sure what it is you wish to do.

Matt already answered (as an answer) that if you just desire to interpolate the function at a specific location, then interpne (from the file exchange) will solve that problem. I believe that is the correct answer from what little I know about your problem.

If you wish to create a new surface that extends further out, gridfit can do that. Just supply the new coordinates that include the old x and y, but will create a new grid that extends out as far as you wish.

Remember that any extrapolation is a risky business, trying to predict the unknown. It is RARELY accurate or even intelligent, and is only as good as the data it is built on. If you are trying to extrapolate a long distance out, don't be upset at the result.

25 Mar 2011 jose tapia

Hi,
I'm having the following problem. I have a 2d table that contains data from some measurements. How can I do the same kind of Extrapolation that is possible in SIMULINK 2-D table lookup using interpolation-extrapolation lookup method, but in Matlab. As I figured out 'griddata' and 'interp2' can not do the job for me.
This is the dimension of my data:
x <1x37 double>
y <1x28 double>
z <28x37 double>

Thanks in advance.

15 Mar 2011 Mikkel Lund

Very nice.. Just what I have been looking for

13 Mar 2011 John D'Errico

Jason - This is a problem of extrapolation, something often difficult to do intelligently. My standard response is to quote Mark Twain:

“In the space of one hundred and seventy six years the Lower Mississippi has shortened itself two hundred and forty-two miles. That is an average of a trifle over a mile and a third per year. Therefore, any calm person, who is not blind or idiotic, can see that in the Old Oölitic Silurian Period, just a million years ago next November, the Lower Mississippi was upwards of one million three hundred thousand miles long, and stuck out over the Gulf of Mexico like a fishing-pole. And by the same token any person can see that seven hundred and forty-two years from now the Lower Mississippi will be only a mile and three-quarters long, and Cairo [Illinois] and New Orleans will have joined their streets together and be plodding comfortably along under a single mayor and a mutual board of aldermen. There is something fascinating about science. One gets such wholesale returns of conjecture out of such a trifling investment of fact.”

"Life on the Mississippi", Mark Twain, 1884

The point is, extrapolation is difficult, especially if all you use is a tool that tries to take data that it sees and predict into the unknown. How should matlab know that zero is a special number? How should gridfit know that it is ok to predict smoothly beyond the extent of your data, yet stop at zero?

Having said all of that, there are several options open to you, both of which might work nicely. First, you can use the 'springs' method in gridfit. It tries to prevent extrapolation beyond your data. This is sometimes the proper solution. Try it, and you might be happy, or not. This I cannot say.

A second option is a transformation. Very often when a system has a property that it cannot go less than zero, you are working in the wrong "space". The trick here is to use a transformation to fix it. Here, I might log your data. Now, use gridfit to model the surface, allowing it to extrapolate smoothly as it wishes with the default options. Now, exponentiate the result. In effect, you have done an interpolation in the log domain, but then transformed back. Do you care if the logs went negative? Of course not. exp is a function that is positive for all real inputs. This trick often works very nicely. In effect, it is just recognizing that the interpolation is best done in the proper domain, one where your system is truly more additive.

02 Mar 2011 Jason

I have scattered 3d data of a membrane subject to tensile loads at the corners resulting in transverse displacement. I'm trying to use gridfit to look at a the 3d surface representing the data. The problem is around the edges the representation gridfit and surf produce is very inaccurate. In fact, at the corners the resulting plot shows negative values when there are no negative z displacements in the data! Any ideas? Thanks.

02 Mar 2011 John D'Errico

There is no functional form for the surface, at least unless you are willing to accept a tool like interp2 to interpolate it at any point. The best that you can do is to view that surface as a low order spline in two dimensions, but as a spline, all you have are a large number of tiny linear segments all neatly connected together. (This applies to the default method, which uses a triangulation of the regular lattice. The alternative is the bilinear interpolant, which in effect is not even truly linear, but an oddly and mildly piecewise quadratic form.)

The virtue of this approach is you can fit any set of data that follows the general form of a function z = f(x,y). You need not formulate a model as other modeling tools force you to do. But you can't get a model out of this either.

If you really want to do empirical modeling, you need to make the effort of posing a model, of choosing some form that will represent your surface. That model may be polynomial (as my polyfitn would allow you to fit) or it may be nonlinear. Then you need to use a tool that can fit that model to your data. There are many such tools on the file exchange to serve that purpose, or you can use the curve fitting toolbox.

02 Mar 2011 Lennert

Hi there! I have a question in relation to this tool. Suppose you have already fit a surface through the data. Is there a command to extract the formula/function of that surface? Thanks!

21 Feb 2011 Michael Brown

Just wanted to extend my gratitude that you made this available. The documentation is excellent, it works as advertised and saved me several days of setting up and debugging a regularization scheme for fitting an appropriately smoothed surface through some unequally spaced thermodynamic data that have error and do not extend to all corners of the needed P-T grid. It was a relief to be able to focus on the science rather than on debugging the scheme to process the data.

29 Dec 2010 Toni

Sorry!
Still stumbling about these issues.
Thanks for the right directions!

29 Dec 2010 John D'Errico

Argh.

Gridfit is set up to ignore any nan data. So your question is meaningless and irrelevant, because gridfit in fact ignores the nan data. See the following code fragment, taken directly from the code...

% also drop any NaN data
x=x(:);
y=y(:);
z=z(:);
k = isnan(x) | isnan(y) | isnan(z);
if any(k)
x(k)=[];
y(k)=[];
z(k)=[];
end

So giving this code a rating of only 4 stars because a different code (griddata) did not work for you seems a bit silly.

If your goal is to fill in the nan elements, while leaving the actual data alone, then use inpaint_nans to interpolate. Use the right tool to solve your problem.

If your goal is to do smoothing, while also interpolating the empty (nan) data, then use gridfit. It will solve your problem, but only if you call gridfit and not griddata.

Finally, if your problem is that griddata is not working as you wish or that you don't understand griddata, then why in the name of god and little green apples, why are you asking it here, in a place intended for comments on gridfit?

29 Dec 2010 Toni

Hello!
I have an 512x512 array z where ~11000 points are non zero and the rest must be nans! Why is the following code giving me an array zgd filled with only nans? How does this program deal with nans?

xi = linspace(0,1,512);
[xg,yg]=meshgrid(xi,xi);
x = 1:512;
y = 1:512;
zgd = griddata(x,y,z,xg,yg);

09 Dec 2010 Lars Robert

excellent job on meteorological balloon data!

05 Oct 2010 Alexis

+1 vote for gridfitn ;)

16 Sep 2010 Ben

Thanks John for your code. I enjoy it very much!

19 Aug 2010 xfactor y  
30 Jul 2010 Adam Chapman

the new setting "'smoothness',[xsmooth ysmooth]" is a Godsend. Many Thanks

22 Jul 2010 John D'Errico

I like Adam's idea, and this is easily enough put into the code.

22 Jul 2010 Adam Chapman

This is brilliant, but a capability to vary the degree of smoothing in x and y exclusively would be very useful

13 Jun 2010 Jeff  
02 Jun 2010 John D'Errico

I've responded via direct e-mail, but the gist of it is that I have plans to introduce a gridfitn one day, but it will take some time to write, and I have many things to write on that same list.

02 Jun 2010 Janine Bijsterbosch

Hi there! Is there a gridfit3d available at all? I have a 3d array of electric field data (369x551x325) on a irregular but highly collinear grid and griddata3 therefore does not work. As I understand, the gridfit code would be perfect for this if it could deal with 3d data. Thanks!

18 Apr 2010 John D'Errico

Alternatively, a simple citation style for a FEX submission is found here, if all you wish is a citation:

http://matlabwiki.mathworks.com/Citing_Files_from_the_File_Exchange

18 Apr 2010 John D'Errico

Perhaps I should add this...

If you really desire something written, many years ago I was granted a pair of US patents on a similar idea, there applied to higher dimensional problems in color modeling. Gridfit only works in 2-dimensions of course.

US Patent # 4992861, 4941039

At least for me, actually reading a patent is about as miserable a pastime as I can imagine.

18 Apr 2010 Matteo Niccoli

Excellent work! I needed a sound extrapolation algorithm to pad some gravity data to use with both vertical continuation and source depth inversion. It worked very well for me. Might I add this has one of the best and most easily understood helps I've seen for Matlab code, even for someone like me whose first language is not English.
SOmeone asked before about journal reference and I did not see your reply. I'd like to know as well if there's anything published we can reference. Thank you.

08 Apr 2010 Jeff Evans  
17 Mar 2010 Panagiotis

John thank you for you prompt reply.

I really need a tool that only interpolates in a given arbitrary domain, as I need to create contour lines in it, which wouldn't make sense out of the domain boundaries.

Could you please post a message here when you post your new tool?
many thanks!

16 Mar 2010 John D'Errico

In answer to the question of extrapolation, there are several issues here.

Gridfit tries to extrapolate gracefully, that is, it will extrapolate as smoothly and linearly as possible. This is the default mode. It is also possible to specify the 'springs' method for gridfit. This method tries to extrapolate minimally, while still generating solutions over the entire domain of interest. Think of that method as extrapolating as a constant, where it can do so, at least as a continuous, smooth function.

Regardless, these methods still do extrapolation. If you absolutely must avoid any extrapolation, leaving the entire domain outside of the convex hull empty, then an option is to use griddata. Griddata is of course an interpolant. It will do no smoothing as gridfit allows you to control. But griddata will interpolate inside the convex hull, and leave points outside as NaN.

The last possibility is if you need to do smoothing, but only inside a given, arbitrary domain, convex or not. This requires a new tool that is like gridfit, one that I have written but not posted on the FEX yet, but will do so.

16 Mar 2010 Panagiotis

Is it possible to totally avoid any extrapolation?

11 Feb 2010 Félix Martin  
19 Dec 2009 Christine A.  
02 Oct 2009 Chamane

excellent

02 Oct 2009 Ulrich  
29 Sep 2009 Michael Adsetts Edberg Hansen  
09 Sep 2009 John D'Errico

Use of an interpolant followed by a smooth is a poor second choice, for several reasons.

Gridfit finds the surface that is as smooth as possible, that is consistent with the data. Smoothing a interpolated surface after the fact does not ensure that the result is consistent with the data. When you do a posterior smooth of the surface, the act of smoothing is now disconnected from the data.

Next, if you use griddata to interpolate a surface in advance, you will only get a result that lies within the convex hull of the data. Griddata will not extrapolate unless you use the v4 option, and that option is VERY slow for any significant number of points. Gridfit can extrapolate using several methods, depending upon your goals.

Extrapolation is an important capability of gridfit. But, extrapolation can come in many different forms. For example, consider a data set which is just slightly concave along one edge of the data. The published demo has a good example of this. See the second example, fitting a trigonometric surface. Along the edges of that data, see that the griddata interpolant generates long, thin triangles. Long thin triangles are terrible for interpolation, so what happens is you see strange interpolation artifacts along the edges.

Gridfit allows you to have replicates in your data, treating them properly in a least squares sense to generate the surface. Try out griddata with replicate data points. Even near replicate points can introduce nasty artifacts in the interpolant. Worse, the delaunay triangulation used in griddata will often have problems if you have sets of collinear data points. This is no problem at all for gridfit.

Finally, you can easily control the extent of smoothing done by gridfit.

In short, griddata has its purposes. There are general circumstances when I recommend griddata. But I would never recommend the use of griddata to be then followed by a smoothing operator.

09 Sep 2009 Yaroslav Bulatov

What is the advantage of using gridfit instead of griddata followed by smoothing?

25 Aug 2009 Daniel Arevalo

Hello John, the function is really nice but I'm wondering if it would be possible to use it with non-rectangular domains where the nodes are just specified with (xn,yn) pairs n=1,..., N
Thanks for any help !

03 Jul 2009 Angel Atanasov

Thank you mr. D'Errico, for making my life easier.You piece of code is really good.

12 Jun 2009 Chris Sherwood

This routine is almost magic...it is the only routine I have found that finds surfaces I know are there amidst very noisy elevation data. Thanks for contributing it.

06 Jun 2009 Karl

This is spectacular. I hope the Mathworks is paying you some sort of royalty for your efforts!

11 May 2009 Michael Dupin  
04 Mar 2009 Bryan Raines  
13 Jan 2009 J Ú  
25 Dec 2008 Michael Jordan  
25 Dec 2008 Jveer

wow! looks awesome! was wondering if there was any way to fit a surface around a 3D scatter? maybe this could be used on the data of for e.g x,y data corresponding to max(z) and then repeating for max(x) and max(y)?

05 Dec 2008 Alexis

Good job John, beautiful tool! Very parameterizable and helps a lot in giving you a feel of how your data is behaving.

Does anyone know how I can anchor the surface to any given point? I'd like to force the surface to pass through a given set of (x,y,z) anchor points, something like establishing an "infinite" weight to these keypoints, while not altering the weight of the other "normal" points...

25 Oct 2008 Cristina Velasco  
20 Oct 2008 Pauline Wong

Versatile smoothing function.

14 Oct 2008 Cedric Testaz  
20 Sep 2008 John D'Errico

The review by Andrey is a bit of self serving braggadocio, since he does not actually offer Matlab code on the file exchange that replicates in any form what gridfit does. By the way, the link he gives is for a .exe file - BEWARE.

His time comparisons are also meaningless of course, since the times referred to are for a now wildly obsolete computer, and comparing an interpolation code to this class of modeling code makes no sense anyway.

Finally, converting an interpolation code to a noise reducing surface modeling code is not at all trivial. I welcome Andrey to do so and provide MATLAB code for the purpose, if that is truly so easy.

20 Sep 2008 Andrey Matseevskiy

You wrote "For example, my computer took 1020 seconds to solve a 500x500 problem." What an old comp do you use? I need 5 second to solve 500x500 problem. Download my prog from http://www.smartfills.com/Html/2D.zip . Compare with your own method. My own is in fact interpolation, but could be modified for approximation without much problems.
Andrey, Kamchatka

11 Sep 2008 gai cs

it is very good toolbox ?thanks

29 Jul 2008 Anmol Agrawal

Excellent utility!! A default setting did everything I needed with minimal error.

08 Mar 2008 Carlos Adrián Vargas Aguilera

Hi John:

I've been using your incredible program for some years and in a wide range of applications. My question is, do you have written any kind of documentation in a scientific journal, for example, for any citation and also to know a little bit more about this super "surface fitting". Any reference can help too.

BTW: VP uses the more silly test to "compare" GRIDDATA and GRIDFIT, and, as John says, he doesn't bother to read that by default the former is 'linear' and the latter uses 'smooth' to 1. But, any way...

21 Feb 2008 Alex Gardner  
13 Feb 2008 V P

1. In strict accordance with the help, xi and yi WERE x-nodes and y-nodes of interest, generated by MESHGRID. If this kind of input is not accepted by GRIDFIT, this means that GRIDFIT must be improved.
2. This is also necessary for comparison with GRIDDATA, which has no problems with my example.
3. If these arguments are still insufficient, and you insist that I am incorrect, then you have to improve the description of xi and yi in your help and define them as VALID INPUT arguments for MESHGRID. However, this will again mean an improvement.
4. If this is not convincing for you, look at the example below. If it is not a strong indication of necessary improvements, then I have to follow the principle of "haende khokh" - this is much better than the cold war around an empty egg.

Z=rand(100);
a=(1:100)';
[X,Y]=meshgrid(a,a);
z=gridfit(X(:),Y(:),Z(:),a,a);
norm(z(:)-Z(:))
ans =
4.9225

12 Feb 2008 John D'Errico

The problem that VP has is he failed to read the help, or apparently bother to understand what gridfit does. I'd recommend reading the help, and perhaps looking at the examples when you don't understand what you are using. I do understand that real men don't need no steenkin help. Perhaps this is the approach that VP has followed.

Gridfit does not interpolate a surface at some random scattered list of points, as would griddata.

Gridfit generates a surface from scattered data on a complete, regular lattice of points. The node arguments are what he might have passed into meshgrid. Again, READ THE HELP.

12 Feb 2008 V P

What is wrong with calling GRIDFIT? Needs improvement.

Z=rand(100);
a=1:100;
[X,Y]=meshgrid(a,a);
good=(X-50).^2+(Y-50).^2>100;
x=X(good);
y=Y(good);
z=Z(good);
bad=~good;
xi=X(bad);
yi=Y(bad);
zi=gridfit(x,y,z,xi,yi);
??? Error using ==> gridfit at 404
xnodes and ynodes must be monotone increasing

22 Jan 2008 Tim Marvin  
10 Dec 2007 Gassilloud rémy

How to use gridfit with polar coordinates ?
[th,r] = meshgrid((0:5:360)*pi/180,0:5:300);
[X,Y] = pol2cart(th,r);
Z=griddata(x,y,z,X,Y);

Thank you,
Rémy

28 Nov 2007 John D'Errico

Henrique - My guess is that you have called gridfit with no arguments. This code is not a gui or command. See the demos for examples of the use of this code. It truly does work, at least when called with data. John

28 Nov 2007 Henrique Mendes

It is not working! Matlab says:
"??? Input argument "x" is undefined.
Error in ==> gridfit at 325
x=x(:);"
I tried to change x(:) to x(:,1), but the problem kepps the same. My Matlab is version 7. It had to be more recent?

05 Nov 2007 Darren Paget

Thanks heaps for the function, it saved me a lot of work. Surface fitting has come a long way! absolutely awesome fit to my data. and really easy to use.

02 Aug 2007 Hendawi Mohamed

I just write to say..."Thank you"!

PS. but I add: It would be perfect if there was a way to sidestep the extrapolation of the the gridnodes beyond the datapoint boundary. It should be added in matlab next release!

24 May 2007 Juan Perez  
17 Apr 2007 gabriele flauti  
01 Feb 2007 David Sprinkle

Great tool and insightful docs! Gridfit succeeds where griddata fails; that is, with noisy & ill-spaced data such as exist in the real world. Thanks!

19 Jan 2007 Ken Lin

Thanks for your programs. It solves my problem now:)
by the way, would you like to work out some way to interpolate the data just like the surfer does

07 Jan 2007 Zhijun Wang

This is very very good work!

14 Dec 2006 E Farhi

This is simply great, with all we can expect for a replacement to griddata (which fails too often).

09 Nov 2006 Carlos Adrián Vargas Aguilera

Very very good John! Take good care of Amy!

02 Nov 2006 Thomas Clark

Excellent function, John - this has been extremely useful to me.

The ability to turn off extrapolation beyond the bounds of the data would be very useful, though.

A note for the inexperienced user:
If you have repetitions (multiple points with the same (x,y) but different (z)), this algorithm will handle where griddata might fall. In such a case, think hard about why you have these repetitions; it may be that volumetric plotting is more appropriate than surfaces. If so, you could get a good-looking surface out of this algorithm; which isn't totally meaningful.

That's not a fault with this algorithm; it just means that you have to think hard about why you are getting repeated points in Z: If it's acceptable to average these out to a surface, this algorithm is for you!

26 Oct 2006 Phani Ivatury

Excellent!!!

13 Oct 2006 Wolfgang Schwanghart

Excellent! This is exactly what I was searching for. Thank you very much for this great contribution.

21 Aug 2006 Christopher Schwalm

Excellent program and documentation.

A question: the help states "Gridfit is not an interpolant." Why is that? It seems to be perfectly able to be used as such, given the outputs and how they were generated. One use I see in my work is to use this function to interpolate where repeated values are present. I've yet to find another Matlab solution that works as well as this does.

19 Jun 2006 John D'Errico

I'll invite Håvard Torpe to contact me via e-mail to discuss this idea. But while I COULD allow the user to turn off extrapolation in gridfit, this will likely cause singularities - serious numerical problems, for most users who might do so.

16 Jun 2006 Håvard Torpe

Works like a charm right from the moment you start using it.

However: I'd like gridfit() NOT to extrapolate values for gridnodes that has no datapoint in the vicinity. Any ideas?

29 May 2006 Andy Lulham  
27 May 2006 thank you

Marvelous!

05 May 2006 Aslak Grinsted

Should be included into the default distribution of matlab! perhaps somehow integrated into griddata?

25 Apr 2006 Wen Soong

Does a great job!

08 Mar 2006 T O

Excellent Work!

06 Mar 2006 Kevin Denis  
21 Feb 2006 David Barker

This is an excellent program. It makes it easy to fit surfaces to 3D data.

11 Feb 2006 Julio Oliveira

Great ! Thank you!

10 Jan 2006 CASTRY Alen

Thank you...

16 Dec 2005 Kong Zuor

Well done John, I have been looking for this while ago and I think it is a must have code. Thank you very much for being so helpful to mathlab community

13 Dec 2005 Alberto Baraldi

I was struggling against griddata since days, beacause of its problems when extrapolating...this is absolutely great!
Thank you for your great help to the Matlab community. We love you

13 Dec 2005 Chaowen Liang

This is exactly what I am looking for... Thanks you!!! Save me a lot of time to write one...

16 Nov 2005 Urs (us) Schwarz

simply said: just another extremely useful member of the d'erriconian-draconian family of must-haves if you're professionally serious about nd-data investigation / reconstruction / interpolation / visualization, which smoothly unites with its now-famous siblings (consolidator, inpaint-nans) into a very nice toolbox for those everyday cheating-with-data endeavors...
the code is very clean and the profiler reveals no unnecessary hot-spot...
HOWEVER, in contrast, the help bit is just unbelievably clumsy and almost indigestible without sever anxiolytic medication - JOHN, JOHN! this part requires some serious overhaul and sprucin'-up...
altogether, very well done - and thank you for this snippet
us

11 Nov 2005 Michael Robbins

Neat picture!

Updates
13 Dec 2005

Allows the user to specify a number of nodes from
min to max of the data. Also fixed it so if the user specifies a list of nodes which does not bound the data, it can either recover smoothly or issue an error, all at the discretion of the user.

23 May 2006

Release 2.0: Adds the option of tiling to solve much larger problems

29 Jul 2010

Allow unequal smoothing parameters in x and y directions

Contact us