%%%建立函数f，包含u0、v0、k1、k2四个未知参数；

function f=myfun(x);

u0=x(1);v0=x(2);k1=x(3);k2=x(4);

%x、y为两个4*4矩阵，存储16个坐标值；

x=[840.224650000000,941.255360000000,1041.80517000000,1141.57622000000;839.663680000000,940.703460000000,1041.51096000000,1141.29222000000;839.365190000000,940.283360000000,1041.03379000000,1140.86218000000;839.423820000000,940.207630000000,1040.62456000000,1140.22984000000];

y=[862.025330000000,862.747530000000,863.674580000000,864.728420000000;962.966060000000,963.522130000000,963.916880000000,964.733440000000;1064.09411000000,1064.46505000000,1064.61563000000,1064.77008000000;1164.86025000000,1164.85888000000,1164.80380000000,1164.52098000000];

%建立两个4*4矩阵，将x、y两个矩阵中每个坐标值更新为含有u0、v0、k1、k2未知参数的值分别存入x1、y1矩阵；

x1=ones(4,4);

y1=ones(4,4);

for i=1:4;

for j=1:4;

r=sqrt((x(i,j)-u0)^2+(y(i,j)-v0)^2);

x1(i,j)=x(i,j)+(x(i,j)-u0)*(k1*r^2+k2*r^4);

y1(i,j)=y(i,j)+(y(i,j)-v0)*(k1*r^2+k2*r^4);

end

end

%对x1、y1矩阵中每行和每列中的四个坐标通过jiaobi函数求一个比例关系存入P矩阵；

for k=1:4;

C_hengxiang=jiaobi(x1(k,1),y1(k,1),x1(k,2),y1(k,2),x1(k,3),y1(k,3),x1(k,4),y1(k,4));

P(k,1)=abs(C_hengxiang-4/3);

C_zongxiang=jiaobi(x1(1,k),y1(1,k),x1(2,k),y1(2,k),x1(3,k),y1(3,k),x1(4,k),y1(4,k));

P(k+4,1)=abs(C_zongxiang-4/3);

end

%对x1、y1矩阵中两对角线中的四个坐标通过jiaobi函数求一个比例关系存入P矩阵；

C_xiexiang1=jiaobi(x1(1,1),y1(1,1),x1(2,2),y1(2,2),x1(3,3),y1(3,3),x1(4,4),y1(4,4));

P(9,1)=abs(C_xiexiang1-4/3);

C_xiexiang2=jiaobi(x1(1,4),y1(1,4),x1(2,3),y1(2,3),x1(3,2),y1(3,2),x1(4,1),y1(4,1));

P(10,1)=abs(C_xiexiang2-4/3);

%待求最小值的函数为矩阵P中各项的和

f=sum(P);