How to fit ellipse equation through segment of ellipse?
32 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Dear all,
I have noisy x and y data which should follow more or less an ellipse equation (see plot.png).
x=[1156;1155;1154;1153;1152;1151;1150;1149;1148;1147;1146;1145;1144;1143;1142;1141;1140;1139;1138;1137;1136;1135;1134;1133;1132;1131;1130;1129;1128;1127;1126;1125;1124;1123;1122;1121;1120;1119;1118;1117;1116;1115;1114;1113;1112;1111;1110;1109;1108;1107;1106;1105;1104;1103;1102;1101;1100;1099;1098;1097;1096;1095;1094;1093;1092;1091;1090;1089;1088;1087;1086;1085;1084;1083;1082;1081;1081;1080;1079;1079;1078;1077;1076;1075;1074;1073;1072;1072;1071;1070;1069;1068;1068;1067;1066;1066;1065;1064;1063;1063;1063;1062;1062;1061;1061;1061;1060;1059;1058;1058;1057;1056;1056;1056;1055;1054;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053];
y=[439;439;439;439;439;439;438;438;438;437;437;437;437;437;437;437;437;437;437;436;436;436;435;435;435;435;435;435;435;435;434;434;434;434;433;432;432;432;432;432;431;431;430;429;429;429;428;427;427;427;427;427;426;425;425;425;424;424;424;424;424;424;423;422;422;421;421;420;420;419;419;418;417;416;415;414;413;412;411;410;409;409;408;407;406;405;404;403;402;401;400;399;398;397;396;395;394;393;392;391;390;389;388;387;386;385;384;383;382;381;380;379;378;377;376;375;374;373;372;371;370;369;368;367;366;365;364;363;362;361;360];
The data is only the quarter of an ellipse, but a full or half ellipse should be fitted to it (or rather the ellipse equation should be calculated). Futhermore the distance from the lowest to the highest point (in y direction) of the quarter of the ellipse should be used as a constrain. In this case this it is 80 for the quarter or half of the ellipse and 160 for the full ellipse. So one parameter of the ellipse equation is given through this. I tried several approaches (for example descripes in https://www.mathworks.com/matlabcentral/answers/98522-how-do-i-fit-an-ellipse-to-my-data-in-matlab) but none of them worked.
If the problem is unclear just ask and I try to explain it better, :)
I would be really thankfull for your help.
Best regards
Simon
0 comentarios
Respuestas (2)
Bruno Luong
el 14 de Jul. de 2019
Editada: Bruno Luong
el 27 de Dic. de 2020
Warning: you won't ever get reilable ellipse with such as small portion of data.
I'll still provide you a code (Optimization toolbox is required)
x=[1156;1155;1154;1153;1152;1151;1150;1149;1148;1147;1146;1145;1144;1143;1142;1141;1140;1139;1138;1137;1136;1135;1134;1133;1132;1131;1130;1129;1128;1127;1126;1125;1124;1123;1122;1121;1120;1119;1118;1117;1116;1115;1114;1113;1112;1111;1110;1109;1108;1107;1106;1105;1104;1103;1102;1101;1100;1099;1098;1097;1096;1095;1094;1093;1092;1091;1090;1089;1088;1087;1086;1085;1084;1083;1082;1081;1081;1080;1079;1079;1078;1077;1076;1075;1074;1073;1072;1072;1071;1070;1069;1068;1068;1067;1066;1066;1065;1064;1063;1063;1063;1062;1062;1061;1061;1061;1060;1059;1058;1058;1057;1056;1056;1056;1055;1054;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053;1053];
y=[439;439;439;439;439;439;438;438;438;437;437;437;437;437;437;437;437;437;437;436;436;436;435;435;435;435;435;435;435;435;434;434;434;434;433;432;432;432;432;432;431;431;430;429;429;429;428;427;427;427;427;427;426;425;425;425;424;424;424;424;424;424;423;422;422;421;421;420;420;419;419;418;417;416;415;414;413;412;411;410;409;409;408;407;406;405;404;403;402;401;400;399;398;397;396;395;394;393;392;391;390;389;388;387;386;385;384;383;382;381;380;379;378;377;376;375;374;373;372;371;370;369;368;367;366;365;364;363;362;361;360];
xc=x-mean(x);
yc=y-mean(y);
M=[xc.^2,yc.^2,xc.*yc,xc,yc];
% Fit ellipse through (xc,yc)
P0 = zeros(5,1);
fun = @(P) norm(M*P-1)^2;
nonlcon = @(P) nlcon(P);
opts = optimoptions(@fmincon, 'Algorithm','Active-set'); % EDIT
P = fmincon(fun,P0,[],[],[],[],[],[],nonlcon,opts); % EDIT
% P=M\ones(size(xc)); % for generic conic
xi = linspace(1000,2000);
yi = linspace(0,500);
[XI,YI]=meshgrid(xi-mean(x),yi-mean(y));
M=[XI(:).^2,YI(:).^2,XI(:).*YI(:),XI(:),YI(:)];
z=reshape(M*P-1,size(XI));
close all
h1=plot(x,y,'r.');
hold on
h2=contour(xi,yi,z,[0 0],'b');
axis equal
xlabel('x')
ylabel('y')
legend('data','fit')
%%
function [c,ceq] = nlcon(P)
B = [P(1) P(3)/2;
P(3)/2 P(2)];
smalleig = 1e-4; % empirical value to ensure the bilinear form is strictly positive
c = smalleig-eig(B);
ceq = [];
end
0 comentarios
Image Analyst
el 14 de Jul. de 2019
4 comentarios
Image Analyst
el 14 de Jul. de 2019
If you want to mirror the points then you could do this
mask = poly2mask(x, y, rows, columns);
imshow(mask);
props = regionprops(mask, 'Orientation', 'MajorAxisLength', 'MinorAxisLength');
Rows and columns are the size of the image, which could be anything as long as it's bigger than your max x and max y.
Then, from props.MajorAxisLength, etc. you should be able to get whatever you want about the ellipse.
Ver también
Categorías
Más información sobre Surface and Mesh Plots en Help Center y File Exchange.
Productos
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!