他们只是使用非线性曲线拟合算法对两个方程 但是，拟合的准确性取决于算法的选择，您必须小心确保和始终为正（通过添加约束）。

解决问题的一种方法如下所示（使用 Octave）。

1) 设置点

fig1 = figure()
gca
[x, y, buttons] = ginput(10)
plot(x, y, 'x-');
hold on;

2) 计算凸包

k = convhull(x,y)
plot(x(k), y(k), "r-")

3) 移除船体内部的点

xx = zeros(length(k)-1)
yy = zeros(length(k)-1)
for i=1:length(k)-1
  xx(i) = x(k(i));
  yy(i) = y(k(i));
end

4) 计算船体的质心

x_cen = mean(xx);
y_cen = mean(yy);
plot(x_cen, y_cen, 'mo');

5) 计算角度并绘制 x,y 与角度的关系图

for i=1:length(xx)
  theta(i) = atan2(yy(i)-y_cen, xx(i)-x_cen);
end
angle = theta*180/pi
fig2 = figure();
plot(angle, yy, 'r-'); hold on;
plot(angle, xx, 'b-');

6）加载优化包

pkg load optim

7) 设置曲线拟合

indep_omega = theta(1:length(theta))
obs_x = xx(1:length(theta)) - x_cen
obs_y = yy(1:length(theta)) - y_cen

% Model functions
fx_omega = @ (p, omega) p(1) - p(2) * sin(omega)
fy_omega = @ (p, omega) p(1) + p(2) * cos(omega)

% Initial values
init_p_x = [0; 0.5]
init_p_y = [0; 0.5]

% Constraints (A.' * p + B >= 0)
A = [0; 1] ; B = 0;
settings = optimset ("inequc", {A, B});

8）做曲线拟合

[p_x, model_values, cvg, outp] = nonlin_curvefit(fx_omega, init_p_x, indep_omega, obs_x', settings)
[p_y, model_values, cvg, outp] = nonlin_curvefit(fy_omega, init_p_y, indep_omega, obs_y', settings)