Fitting Low-Curvature Volume to Scattered Data

Hi there, I am trying to figure out a way to form something like a low-curvature volume fit in spherical coordinates.

So, I first have a grid of azimuth (-pi,pi) and elevation (-pi/2, pi/2) values and a corresponding array for radii, that is R = R(az,el) of dimension NxN. I then determine all maximum values in R using imregionalmax and write teh corrsponding (az,el,R) values into three vectors.

My goal would be now to fit a volume to these maximum values in R on the original az/el grid (currently using griddata,cubic). That is, if I had a flat set of original data in R with only a few equal peaks in it, that surface should be a sphere.

Right now I have tried to incorporate a periodic boundary solution by [R R R; R R R; R R R] and mapping the grid from (-3pi,3pi) and (-3pi/2, 3pi/2) and deleting duplicate boundary entries in R.

The best I could do so far leads to a solution like in the Figure (same data in the plots, once with transparency to see the original data in the center of the volume fit. The original data is a sphere of radius 1 and some peaks distributed on the surface of height 1 and 0.5 (the method works for all equal peaks already, that seemed trivial). But: Obviously, this would not be the surface with lowest curvature, and I guess I am missing something (or maybe even a lot).

Has anyone ever dealt with this kind of problem? I would appreciate any support or ideas to find a solution.

Thank you in advance, Regards from Germany

Answer 1

[John D'Errico]

The problem with converting to az/elev form is it forces you to add the BCs to the problem.

Anyway, it is not clear what you mean by "equal peaks". Do you have scattered data, but with some outliers that you need to discard? Is there noise in the data?

Personally, I'd probably use some sort of a triangulation in 3-d, then relax the shape to fit the surface, so use a penalty for missing the data, but also for high curvature between adjacent triangles. Adjust the relative weights between the two terms, and you will get a smooth surface that tries to fit your data. Since it is already in 3-d, there is no need to worry about the boundary conditions.

Answer 2

[Dennis]

Thank you for your response! So later, the real data will indeed be noisy and might contain outliers. Still, I want to fit a volume to the data points farthest from the center of that point cloud (be it outliers or not). As the data will be given in spherical coordinates, this would imply finding max®. So far, I am using conversion to cartesian coordinates just for reasons of visualization. The whole thing is part of a smoothing procedure I am trying to implement.

Anyway, to test my function I am synthezing data:

[AZ, EL] = meshgrid(linspace(-pi,pi,100),linspace(-pi/2,pi/2,100));
R = ones(100);
R([25,75],[25,75]) = 2;

This will give me a sphere with four distinct but "equal" peaks on its surface. The resulting maxima smoothing volume gives a sphere as well.

If I now added four other peaks with R = 1.5 that volume should transform into a bean- or ellipsoid-shaped volume.

If have not yet worked my way into triangulation methods - can you recommend any examples or the Matlab help to start with?