Since you have X, Y, Z data,you do not need to depend on surf function to find out what the average height (z-coordinate) of the surface is. However, since the surface is non-linear, you cannot simply use MEAN.
There are two methods to accomplish this, depending on the whether or not there are singularities in the surface fit.
If there are discontinuities, the best way to estimate the average height of the surface is to sample a large number of points on the surface and divide by the total number of points sampled. This is similar to a Monte Carlo estimation. The following code shows how do to this.
num_samples = 100;
x_samples = linspace(min(X), max(X), num_samples);
y_samples = linspace(min(Y), max(Y), num_samples);
[xx, yy] = meshgrid(x_samples, y_samples);
z_samples = fitted_surface(xx(:), yy(:));
average_z = mean(z_samples);
In the above code "fitted_surface" is the surface returned by Curve Fitting Toolbox and "X" and "Y" are the X and Y data.
If the surface fit does not have singularities, the INTEGRAL2 function can be used to numerically evaluate the volume under the surface. So you can have a measure of the average height as the volume of the surface divided by the area projected onto the XY plane. The following code demonstrates how to do this
my_fun = @(x, y) fitted_surface(x, y);
volume = integral2(my_fun, min(X), max(X), min(Y), max(Y));
average_height = volume / ((max(X) - min(X)) * (max(Y) - min(Y)));
In the above code "fitted_surface" is the surface fit returned by Curve Fitting Toolbox. The anonymous function "my_fun" is necessary to pass the surface fit into the INTEGRAL2 function. The division divides the volume by the area the surface projects onto the XY plane.
