3D curve fitting a periodic surface

Suppose I have a 3D concave surface defined in an MxN mesh r,tht,z where:

tht = repmat(linspace(0, 2*pi, N), M, 1) is MxN;

r is monotonically increasing MxN

z is monotonically increasing MxN

I need to fit a smooth periodic spline to this data where the innermost r for

each tht is at the fixed specified value r(1, :) > 0

Sample data can be generated by:

N = 10;

M = 5;

tht = repmat(linspace(0, 2*pi, N), M, 1); % tht is uniform

r = rand(M, N);

r(1, :) = 1; % Fix inner circle at r==1

r(:, end) = r(:, 1); % Ensure r at 0 deg == r at 360 deg

r = cumsum(r, 1); % Ensure r is monotonically increasing

[x,y] = pol2cart(tht, r);

rmax = max(r(:)) * 1.1;

z = sqrt(rmax^2 - x.^2 - y.^2);

figurewsd

surf(x,y,z)

axis equal

grid on

I suspect that I need to use csape in the curve fitting toolbox to fit this

periodic surface with this end condition, however I can't quite understand the

documentation and examples.

From the examples provided, I suspect that I need to do something in the form

sph = csape({x,y},v,{'clamped','periodic'});

However, the form of x, y, and z in the provided example confused me.

Does anyone have a suggestion how I could use this (or another) function

to fit my defined surface?

Added complications:

1. Would this change if tht was slightly different than uniform for each r?

2. What if I wanted to clamp the fit at r(end, :), as well?

Thanks so much for any advice.

Answer 1

[Nithin]

Hi @Scott Sycamore,

To achieve a smooth periodic spline fit for a 3D surface data, we can use either "NURBS" (Non-Uniform Rational B-Splines) or "scatteredInterpolant" function. These methods are particularly effective for advanced surface fitting with periodic constraints.

For "NURBS" fitting, you can use the "NURBS" toolbox available on MATLAB Central File Exchange: https://www.mathworks.com/matlabcentral/fileexchange/26390-nurbs-toolbox-by-d-m-spink

Preprocess the data to define control points for the surface and then create periodic knot vectors for both "r" and "tht". Use the control points and knot vectors to construct a "NURBS" surface and fit in the data, ensuring that it is periodic and respect the constraints at "r(1, :)".

You can refer to the below sample code snippet for better understanding:

% NURBS requires control points to be in a 4D homogeneous coordinate format
coefs = zeros(4, M, N);
coefs(1, :, :) = x;
coefs(2, :, :) = y;
coefs(3, :, :) = z;
coefs(4, :, :) = 1; % Homogeneous coordinate

% Define knot vectors
knotU = [zeros(1, 3), linspace(0, 1, M-2), ones(1, 3)];
knotV = [zeros(1, 3), linspace(0, 1, N-2), ones(1, 3)];

surface = nrbmak(coefs, {knotU, knotV});

The image below shows the fitting after implementing the "NURBS" fitting:

For spline fitting using the "scatteredInterpolant" function, first flatten the data and create a grid for interpolation using "meshgrid". Then, use the "scatteredInterpolant" function for fitting the flattened data. Finally, evaluate the spline fit on the grid using "spline_fit" function.

You can refer to the below sample code snippet for better understanding:

% Flatten the data for fitting
x_flat = x(:);
y_flat = y(:);
z_flat = z(:);

% Create a grid for interpolation
[xq, yq] = meshgrid(linspace(min(x_flat), max(x_flat), 100), ...
                    linspace(min(y_flat), max(y_flat), 100));

spline_fit = scatteredInterpolant(x_flat, y_flat, z_flat, 'linear', 'linear');
zq = spline_fit(xq, yq);

The image below shows the fitting after implementing the "scatteredinterpolant" funciton:

If "tht" is non-uniform, you will need to adjust the knot vector for the angular direction accordingly. This may involve using non-uniform spacing in the knot vector to reflect the spacing in "tht"

To clamp the fit at "r(end, :)", ensure that your control points and knot vectors are correctly defined to enforce these boundary conditions. This might involve setting additional constraints or adjusting the weights of the "NURBS" surface.

Kindly refer to the following MathWorks documentations to know more about the function used:

scatteredInterpolant: https://www.mathworks.com/help/matlab/ref/scatteredinterpolant.html

meshgrid: https://www.mathworks.com/help/matlab/ref/meshgrid.html

spline_fit: https://www.mathworks.com/help/matlab/ref/spline.html