大贤者
精华
|
战斗力 鹅
|
回帖 0
注册时间 2008-5-24
|
本帖最后由 玉米黍 于 2014-5-26 11:20 编辑
起因是这一帖： [讨论] 以数学方法验证条河麻耶的脚线美的形式原理
http://bbs.saraba1st.com/2b/thread-1026428-1-2.html
个人认为日本人的验证方法过于粗糙。比如无法说明两对比图之间的比例问题。于是本着严谨的态度，利用简单的工程数据分析法对小腿曲线函数为： y=1/3*sinθ 的说法进行论证。
这里我采用南宫那月的腿部作为本次研究的对象。
首先用PS裁剪所需的腿部区域，并去掉地面的纹理方便图像处理。接着上传图片至外链网站。
用Matlab导入图片 我们得到一个319*117 的3D矩阵。 由于此图黑丝的颜色对比大，所以单选其中一个（红色）矩阵即可计算出腿部的曲线数据。然后通过图像处理得到一个 长度为 x=300的点数据。
利用Matlab的命令。我尝试用多项式和一些非多项式（比如三角函数）来符合数据（此数据为整个腿部曲线）。发现在没有经验公式可用的情况下效果都不理想。于是只能退而求其次研究简单的小腿弧线。
传说中的经验公式y=1/3*sinθ 是三角函数，因此猜测的公式为： y=a*sin(b*x+c).
结果如图：
红色为腿部曲线，蓝色为符合的小腿曲线函数。
具体数据 a=15.9943 b=0.0191 c= 5.8051
原图中南宫那月的身高约为1000像素，其实际身高约为1.4米（1400mm).所以压缩比为1:1.4
通过还原比例最后的结果表明其小腿曲线函数为： y=22.4sin(0.0136x+4.15） (单位cm）
对比图：
通过比例换算我们可以轻易的得出（15.9943*0.0191=0.305） y=0.305sin(th+0.184)
可见其小腿曲线是接近1/3sin(th) 这个经验公式的。
代码：
% load an image from a web site,
% save the scaled image to a local file
clear all; close all; clc;
url = 'http://image16.poco.cn/mypoco/myphoto/20140525/03/4528375020140525035544094.jpg';
localfile1 = 'tui.jpg';
webfile = 'tui.htm';
img = imread(url); % load the desired image
imwrite(img,localfile1); % save it to a local file
% make a web page to display the images
fp = fopen(webfile,'wt');
fprintf(fp,'<center>\n');
fprintf(fp,'<img src=%s><br>Original Image<p>\n',localfile1);
fprintf(fp,'</center>\n');
fclose(fp);
web(webfile); % display the web page
% use red noly
red1=img(:,:,1);
x1=length(red1)-19;
x=1:x1;
y = [];
for i=1:x1
yy= red1(:,i);
h=length(yy)
hh = 0;
for n=1:h
if yy(n)>=252
hh=hh+1;
else
break
end
end
y(i)=54-hh;%max(hh)=54
end
plot(x,y,'r-');
axis equal % x and y axis units are set equal
title( 'tui bu qu xian');
y1=y(58:162);%use the region 58-162 for calf
x1=x(58:162);
% nlinfit
funky=@(b,x1)b(1)*sin(b(2)*x1+b(3));
beta = nlinfit(x1,y1,funky,[0.1 0.1 0.1 ])
hold on;
plot(x1,funky(beta,x1),'b-')
A=beta(1);
w=beta(2);
p=beta(3);
txt=sprintf('Curve fit: y= %.2fsin(%.2f x + %.2f),(58<=x<=62)',A,w,p);
legend('Data',txt);
xlabel('x');
ylabel('y');
-------5.26更新-------------------------------------------------------------------
分段函数这么简单的方法都忘了。。。。。
最后结果如图：
| 0.0016x^2 + 0.0728x - 0.2705 , 0<=x<=58
y(x) = | 15.9943sin( 0.0191x + 5.8051) , 58<x<160
| 0.0034x^2 - 1.1308x + 103.0777 , 160<=x<=210
| -27.9595sin( -0.0111x - 17.0257) , 210<x<=300
找了一个垂直视角的图进一步分析：
简单的通过椭圆公式模拟曲面，最后结果如图：
更新最终的3D效果图：
所用腿部图片：
代码：（这个乱点，懒得改备注了）
% load an image from a web site,
% save the scaled image to a local file
clear all; close all; clc;
url2 = 'http://image16.poco.cn/mypoco/myphoto/20140526/04/4528375020140526045514062.jpg';
localfile1 = 'tui2.jpg';
webfile2 = 'tui2.htm';
img2 = imread(url2); % load the desired image
imwrite(img2,localfile1); % save it to a local file
% make a web page to display the images
fp = fopen(webfile2,'wt');
fprintf(fp,'<center>\n');
fprintf(fp,'<img src=%s><br>Original Image<p>\n',localfile1);
fprintf(fp,'</center>\n');
fclose(fp);
web(webfile2); % display the web page
% use red noly
red1=img2(:,:,1);
x1=length(red1);
x=1:x1;
y2= [];
for i=1:x1
yy= red1(:,i);
h=length(yy);
hh = 0;
for n=1:h
if yy(n)>=198
hh=hh+1;
else
break
end
end
y(i)=63-hh;%max(hh)=54
end
plot(x,y,'r-');
axis equal % x and y axis units are set equal
title( 'tui bu qu xian');
for i=1:x1
yy2= red1(:,i);
h=length(yy2);
hh = 0;
for n=1:h
if yy2(n)<=248
hh=hh+1;
else
break
end
end
y2(i)=64-hh;%max(hh)=54
end
hold on
% nlinfit
plot(x,y2,'b-')
xlabel('x');
ylabel('y');
red2=flipud(red1);
x1=length(red2);
x=1:x1;
y3= [];
for i=1:x1
yy= red2(:,i);
h=length(yy);
hh = 0;
for n=1:h
if yy(n)>=198
hh=hh+1;
else
break
end
end
y3(i)=hh-60;%max(hh)=54
end
hold on
% nlinfit
plot(x,y3,'b-')
ya=[];
yb=[];
for i=1:x1
ya(i)=y(i)+y2(i);
yb(i)=abs(y3(i)-y2(i));
end
url = 'http://image16.poco.cn/mypoco/myphoto/20140525/03/4528375020140525035544094.jpg';
localfile1 = 'tui.jpg';
webfile = 'tui.htm';
img = imread(url); % load the desired image
imwrite(img,localfile1); % save it to a local file
% use red noly
red1=img(:,:,1);
xs=length(red1)-19;
ys = [];
for i=1:xs
yy= red1(:,i);
h=length(yy)
hh = 0;
for n=1:h
if yy(n)>=252
hh=hh+1;
else
break
end
end
ys(i)=62-hh;%max(hh)=54
end
yx=ys;
%plot(x,yx(1:167),'b-')
A = [];
B = [];
for n=1:x1
for i=1:50
if i<=ya(n)
yc(i)=((1-(i-1)^2/(ya(n))^2)*yx(n)^2)^0.5;
else
yc(i)=0;
end
end
for ii=1:50
if ii<=yb(n)
yc2(ii)=(1-(ii-1)^2/(yb(n))^2)^0.5*yx(n);
else
yc2(ii)=0;
end
end
A = [A; yc];
B = [B; yc2];
end
figure(2)
yd=A';
ye=flipud(B');
z=[ye ;yd];
[x,y] = meshgrid([1:1:167],[1:1:100]);
% draw the surface
h = surf(x,y,z);
axis equal % x and y axis units are set equal
%img = imread(['http://image16.poco.cn/mypoco/myphoto/20140525/07/452837502014052507454102.jpg?690x951_120']);
%set(h,'CData',img,'FaceColor','texturemap')
axis off
------------------------------------------------------
把原图载入一下直接变成3D黑丝。。。。。
|
|
沪公网安备 31010702007642号 )
发表于 2014-5-25 08:42
发表于 2014-5-25 08:55
发表于 2014-5-25 10:51
