Generating a curve around the point on a scatter plot

Question

I am generating a scatter plot containing data from multiple sources, as displayed below.

I would like to be able to generate a curve surrounding an arbitrary query point and passing through points on scatter plot. Final goal is to calculate the area between the lines on the plot.

I have implemented solution using finding points with knnsearch in a circular fashion and then applying hampel filter to eliminate noise. In the example below, I have selected a point right about in the middle of the blue-shaded area. As you can see, the result is far from perfect, and I need more precision.

I am looking for something similar to boundary function, but to work from the inside of the point cloud, not from the outside.

score 1 · Answer 1 · answered Nov 30 '16 at 16:03

1

Final goal is to calculate the area between the lines on the plot.

I would do it differently. Just take any two lines of the plot, calculate the area under the curves with some kind of numerical approximation (for example trapezoidal numerical integration), then subtract the areas and obtain the area between the lines.

answered Nov 30 '16 at 16:03

NoDataDumpNoContribution

10,591
9
64
104

Thanks for the idea to use curves. I had to take much different approach for the rest of the solution, but this got me started. – Nikola Malešević Dec 01 '16 at 16:31

score 1 · Answer 2 · edited Jun 20 '20 at 09:12

Thank to idea in Trilarion's answer, I was able to come up with the better solution.

Note that I use notation for YZ plane instead of XY (to keep consistent with robot coordinate system).

Solution

Generate curves for each set of scatter data

% Scatter data is in iy and iz vectors.
curve = fit(iy, iz, 'smoothingspline', 'SmoothingParam', 0.5);
% Remove outliers.
fdata = feval(curve, iy);
I = abs(fdata - iz) > 0.5 * std(iz);
outliers = excludedata(iy, iz, 'indices', I);
% Final curve without outliers.
curve = fit(iy, iz, 'smoothingspline', 'Exclude', outliers, 'SmoothingParam', 0.5);

Plot curves and scatter data

% Color maps generated by MATLAB's colormap function.
h_curve = plot(curve);
set(h_curve, 'Color', color_map_light(i,:));    
scatter(iy, iz, '.', 'MarkerFaceColor', color_map(i,:))

Let user provide an input by selecting points

User selects one point as a query point and two points for limits along Y axis. This is because some curves come close, but never intersect.

[cs_position.y, cs_position.z] = ginput(1);
[cs_area_limits, ~] = ginput(2);

if cs_area_limits(1) > cs_area_limits(2)
  cs_area_limits = flipud(cs_area_limits);
end

plot_cross_section(cs_position);

Finally calculate and plot surface area

This section uses fantastic answer by Doresoom.

function [ ] = plot_cross_section(query_point)
%PLOT_CROSS_SECTION Calculates and plots cross-section area.
%   query_point  Query point.

  % Find values on query point's Y on each of the curves.
  z_values = cellfun(@(x, y) feval(x, y),...
    curves, num2cell(ones(size(curves)) * query_point.y))

  % Find which curves are right above and below the query point.
  id_top = find(z_values >= query_point.z, 1, 'first')
  id_bottom = find(z_values < query_point.z, 1, 'last')

  if isempty(id_top) || isempty(id_bottom)      
    return
  end

  % Generate points along curves on the range over Y.
  y_range = cs_area_limits(1):0.1:cs_area_limits(2);
  z_top = feval(curves{id_top}, y_range).';
  z_bottom = feval(curves{id_bottom}, y_range).';

  % Plot area.
  Y = [ y_range, fliplr(y_range) ];
  Z = [ z_top, fliplr(z_bottom) ];
  fill(Y, Z, 'b', 'LineStyle', 'none')
  alpha 0.5
  hold on

  % Calculate area and show to user.
  cs_area = polyarea(Y, Z);
  area_string = sprintf('%.2f mm^2', cs_area);
  text(0, -3, area_string, 'HorizontalAlignment', 'center')
end