function D=FractalDim(y,cellmax)
%求输入一维信号的计盒分形维数
%y是一维信号
%cellmax:方格子的最大边长,可以取2的偶数次幂次(1,2,4,8...),取大于数据长度的偶数
%D是y的计盒维数（一般情况下D>=1）,D=lim(log(N(e))/log(k/e)),
if cellmax<length(y)
error('cellmax must be larger than input signal!')
end
L=length(y);%输入样点的个数
y_min=min(y);
%移位操作，将y_min移到坐标0点
y_shift=y-y_min;
%重采样，使总点数等于cellmax+1
x_ord=[0:L-1]./(L-1);
xx_ord=[0:cellmax]./(cellmax);
y_interp=interp1(x_ord,y_shift,xx_ord);
%按比例缩放y，使最大值为2^^c
ys_max=max(y_interp);
factory=cellmax/ys_max;
yy=abs(y_interp*factory);
t=log2(cellmax)+1;%叠代次数
for e=1:t
Ne=0;%累积覆盖信号的格子的总数
cellsize=2^(e-1);%每次的格子大小
NumSeg(e)=cellmax/cellsize;%横轴划分成的段数
for j=1:NumSeg(e) %由横轴第一个段起通过计算纵轴跨越的格子数累积N(e)
begin=cellsize*(j-1)+1;%每一段的起始
tail=cellsize*j+1;
seg=[begin:tail];%段坐标
yy_max=max(yy(seg));
yy_min=min(yy(seg));
up=ceil(yy_max/cellsize);
down=floor(yy_min/cellsize);
Ns=up-down;% 本段曲线占有的格子数
Ne=Ne+Ns;%累加每一段覆盖曲线的格子数
end
N(e)=Ne;%记录每e下的N(e)
end
%对log(N(e))和log(k/e)进行最小二乘的一次曲线拟合,斜率就是D
r=-diff(log2(N));%去掉r超过2和小于1的野点数据
id=find(r<=2&r>=1);%保留的数据点
Ne=N(id);
e=NumSeg(id);
P=polyfit(log2(e),log2(Ne),1);%一次曲线拟合返回斜率和截距
D=P(1);
Dec 27, 2009
根据计盒维数原理求一维曲线分形维数的matlab程序
프랙탈의 반복성과 자기유사성 문제
%px,py, r은 각각 큰 원의 원의 중심 좌표와 반지름이며 n는 반복회수이다.
%사용방법:
%cycle6(0,0,1,4);axis equal
t=linspace(-pi,pi);
plot(r*cos(t)+px,r*sin(t)+py)
theta=pi/3*(0:5);
newr=r/3;
newpx=px+2*newr*cos(theta);
newpy=py+2*newr*sin(theta);
if n>1
for i=1:6
hold on;cycle6(newpx(i),newpy(i),newr,n-1);
end
end
Dec 26, 2009
프랙탈개념: Cantor집합(소스코드)
Cantor집합
c=0.2;d=2;
if (bx-ax)>c
x=[ax,bx];y=[ay,by];hold on;
plot(x,y,'LineWidth',5);hold off;
cx=ax+(bx-ax)/3;
cy=ay-d;
dx=bx-(bx-ax)/3;
dy=by-d;
ay=ay-d;
by=by-d;
cantor(ax,ay,cx,cy);
cantor(dx,dy,bx,by);
end
matlab, 프랙탈차원분석(삼각형에서 중간부분 빼기)
clear,clf;
n=10000;
p=rand(n,1);x0=0;y0=0;
x=[x0;zeros(n-1,1)];
y=[y0;zeros(n-1,1)];
for i=2:n
pp=p(i);
if pp<0.333
x(i)=0.5*x(i-1);
y(i)=0.5*y(i-1);
elseif pp<0.666
x(i)=0.5*x(i-1)+0.25;
y(i)=0.5*y(i-1)+3.433;
else
x(i)=0.5*x(i-1)+0.5;
y(i)=0.5*y(i-1);
end
end
plot(x(1:n),y(1:n),'.')
Dec 19, 2009
Matlab동영상제작- 저장하기
MATLAB动画保存只对电影动画有意义，因为其他两种都是实时动画，一眨眼过去了，而电影动画是先将动画一帧一帧的保存下来，在使用movie函数播放。它的好处是，运行一次MATLAB程序就可以播放无数次，只要你的帧数据还在。
但是这还是不方便，由于它没法脱离MATLAB环境，很讨厌。还好MATLAB为我们提供了movie2avi函数，它可以把动画直接转换成avi文件，而avi文件则可以脱离Matalb环境而在其他地方运行了。
2:请教如何保存matlab的动画？
function avimake
warning off;
load data
[M,N,K]=size(data);
data=data/(max(abs(data(:))));
aviobj = avifile('mymovie.avi','fps',10);
for kk=1:10:K
imagesc(data(:,:,kk));
set(gca,'clim',[-1 1]);
colormap(hsv(128))
frame = getframe(gca);
aviobj = addframe(aviobj,frame);
end
aviobj = close(aviobj);
3:直接保存gif动画
m(:,k)=getframe;
%%%构造gif图像的帧,
nn(:,:,:)=getframe;
%%转换为可以直接输出的格式(这会是图像丢失)
%如果要制作彩色的图像,你只能把生成的彩色图像单独制作(使用其他软件)
nn1=nn.cdata;
nn1=rgb2gray(nn1);
imwrite(nn1,'out.gif','gif','WriteMode','append')
MATLAB의 실행과정을 동영상으로 만들 수 있는지?
문제점: MATLAB으로 만든 프로그램을 실행시킬 때, 그 과정을 GIF 포맷으로 저장할 수 있는지?
평면으로 단엽쌍곡선 횡단면을 자르는 GIF 애니메이션 프로그램
clear
clc
filename='mianjiedanyeshuangqu.gif';
a=3;
b=3;
c=4;
z = -1.5*c:c/20:1.5*c;
[r,alpha1] = meshgrid(z,linspace(0,2*pi,30));
X=a*sqrt(1+r.^2/c^2).*cos(alpha1);
Y=b*sqrt(1+r.^2/c^2).*sin(alpha1);
Z=repmat(z,30,1);
backcolor=get(gcf,'Color');
uicontrol(gcf,'style','text','units','normalized','pos',[0.74 0.02 0.25 0.15],...
'string',{'제작자','QQ：56330069','대학교명'},'fontsize',12,'fontweight','bold',...
'fontunits','normalized','Hor','left','ForegroundColor',[0.5,0.5,0.5],...
'bac',backcolor)
set(gcf,'DoubleBuffer','on')
surf(X,Y,Z)
hold on
%colormap gray
[X,Y]=meshgrid(-6:0.5:6,-6:0.5:6);
Z=-6*ones(size(X));
h=surf(X,Y,Z);
shading interp
alpha(0.6);
quiver3(0,0,0,-1,0,0,9,'k','filled','LineWidth',2);
quiver3(0,0,0,0,-1,0,9,'k','filled','LineWidth',2);
quiver3(0,0,0,0,0,1,7,'k','filled','LineWidth',2);
text(0,-0.8,7,'Z')
text(0,-9,0.8,'Y')
text(-9,0,1,'X')
axis equal
%axis([-7,7,-7,7,-7,7])
axis off
title('평면으로 단엽쌍곡선을 자르는 데모','fontsize',14,'fontweight','bold')
theta=0:pi/10:2*pi;
r=sqrt(1+36/c^2);
h1=plot3(a*r*cos(theta),b*r*sin(theta),-6*ones(size(theta)),'-k','LineWidth',1.5);
for z=-6:0.4:6
set(h, 'xdata' ,X, 'ydata' ,Y, 'zdata' ,z*ones(size(X)));
R=sqrt(1+z^2/c^2);
xx=a*R*cos(theta);
yy=b*R*sin(theta);
zz=z*ones(size(xx));
set(h1, 'xdata' ,xx, 'ydata' ,yy, 'zdata' ,zz);
drawnow; % 새로고침
pause(0.05)
f=getframe(gcf);
imind=frame2im(f);
[imind,cm] = rgb2ind(imind,256);
if z==-6
imwrite(imind,cm,filename,'gif', 'Loopcount',inf,'DelayTime',0.3);
else
imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.3);
end
end
delete(h);
[X,Z]=meshgrid(-6:0.5:6,-6:0.5:6);
Y=-5.4083*ones(size(X));
h=surf(X,Y,Z);
shading interp
alpha(0.6);
h1=plot3(0,-5.4083,-6,'-k','LineWidth',1.5);
h2=plot3(0,-5.4083,6,'-k','LineWidth',1.5);
xx=-5.4083:0.1:5.4083;
for y=-5.4083:0.25:-3
set(h, 'xdata' ,X, 'ydata' ,y*ones(size(X)), 'zdata' ,Z);
yy=y*ones(size(xx));
zz1=-c*sqrt(xx.^2/a^2+yy.^2/b^2-1);
zz2=c*sqrt(xx.^2/a^2+yy.^2/b^2-1);
set(h1, 'xdata' ,xx(zz1>-6), 'ydata' ,yy(zz1>-6), 'zdata' ,zz1(zz1>-6));
set(h2, 'xdata' ,xx(zz2<6), 'ydata' ,yy(zz2<6), 'zdata' ,zz2(zz2<6));
drawnow; % 새로고침
pause(0.05)
f=getframe(gcf);
imind=frame2im(f);
[imind,cm] = rgb2ind(imind,256);
imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.3);
end
zz=-6:0.5:6;
for y=-3:0.25:3
set(h, 'xdata' ,X, 'ydata' ,y*ones(size(X)), 'zdata' ,Z);
yy=y*ones(size(zz));
xx1=-a*sqrt(1+zz.^2/c^2-yy.^2/b^2);
xx2=a*sqrt(1+zz.^2/c^2-yy.^2/b^2);
set(h1, 'xdata' ,xx1, 'ydata' ,yy, 'zdata' ,zz);
set(h2, 'xdata' ,xx2, 'ydata' ,yy, 'zdata' ,zz);
drawnow; % 새로고침
pause(0.05)
f=getframe(gcf);
imind=frame2im(f);
[imind,cm] = rgb2ind(imind,256);
imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.3);
end
for y=3:0.25:5.4083+0.25
set(h, 'xdata' ,X, 'ydata' ,y*ones(size(X)), 'zdata' ,Z);
yy=y*ones(size(xx));
zz1=-c*sqrt(xx.^2/a^2+yy.^2/b^2-1);
zz2=c*sqrt(xx.^2/a^2+yy.^2/b^2-1);
set(h1, 'xdata' ,xx(zz1>-6), 'ydata' ,yy(zz1>-6), 'zdata' ,zz1(zz1>-6));
set(h2, 'xdata' ,xx(zz2<6), 'ydata' ,yy(zz2<6), 'zdata' ,zz2(zz2<6));
drawnow; % 새로고침
pause(0.05)
f=getframe(gcf);
imind=frame2im(f);
[imind,cm] = rgb2ind(imind,256);
imwrite(imind,cm,filename,'gif','WriteMode','append','DelayTime',0.3);
end
delete(h);
hold off
imshow()와 image()의 구별점
matlab에서 연속적으로 폴더의 파일 읽어드리는 방법(1)
dList =dir('*.jpg');
%만약 JPG파일포맷이 아닌 다른 파일 포맷이 있을 경우,반복적으로 [dir 함수]를 사용하기 바람.
% 만약 이미지가 다른 디렉토리에 있을 경우, ["경로" ".확장명"] 으로 처리가능
k =length(dList);
for i=1:1:k
image_data{i}=imread(dList(i).name);
end
%이러한 방식으로 이미지들을 모두 셀형 변수(cell) image_data 에 저장함.
%다음은 사용방식
%1초에 이미지 하나씩 보여줌.
for j=1:1:k
image(image_data{j});
pause(1);
end[/code]
Dec 14, 2009
Rossler 카오스시스템 Lyapunov지수(Matlab 코드)
function dX = Rossler_ly(t,X)
% Rossler 吸引子，用来计算Lyapunov指数
% a=0.15,b=0.20,c=10.0
% dx/dt = -y-z,
% dy/dt = x+ay,
% dz/dt = b+z(x-c),
a = 0.15;
b = 0.20;
c = 10.0;
x=X(1); y=X(2); z=X(3);
% Y的三个列向量为相互正交的单位向量
Y = [X(4), X(7), X(10);
X(5), X(8), X(11);
X(6), X(9), X(12)];
% 输出向量的初始化，必不可少
dX = zeros(12,1);
% Rossler吸引子
dX(1) = -y-z;
dX(2) = x+a*y;
dX(3) = b+z*(x-c);
% Rossler吸引子的Jacobi矩阵
Jaco = [0 -1 -1;
1 a 0;
z 0 x-c];
dX(4:12) = Jaco*Y;
求解LE代码：
% 计算Rossler吸引子的Lyapunov指数
clear;
yinit = [1,1,1];
orthmatrix = [1 0 0;
0 1 0;
0 0 1];
a = 0.15;
b = 0.20;
c = 10.0;
y = zeros(12,1);
% 初始化输入
y(1:3) = yinit;
y(4:12) = orthmatrix;
tstart = 0; % 时间初始值
tstep = 1e-3; % 时间步长
wholetimes = 1e5; % 总的循环次数
steps = 10; % 每次演化的步数
iteratetimes = wholetimes/steps; % 演化的次数
mod = zeros(3,1);
lp = zeros(3,1);
% 初始化三个Lyapunov指数
Lyapunov1 = zeros(iteratetimes,1);
Lyapunov2 = zeros(iteratetimes,1);
Lyapunov3 = zeros(iteratetimes,1);
for i=1:iteratetimes
tspan = tstart:tstep:(tstart + tstep*steps);
[T,Y] = ode45('Rossler_ly', tspan, y);
% 取积分得到的最后一个时刻的值
y = Y(size(Y,1),:);
% 重新定义起始时刻
tstart = tstart + tstep*steps;
y0 = [y(4) y(7) y(10);
y(5) y(8) y(11);
y(6) y(9) y(12)];
%正交化
y0 = ThreeGS(y0);
% 取三个向量的模
mod(1) = sqrt(y0(:,1)'*y0(:,1));
mod(2) = sqrt(y0(:,2)'*y0(:,2));
mod(3) = sqrt(y0(:,3)'*y0(:,3));
y0(:,1) = y0(:,1)/mod(1);
y0(:,2) = y0(:,2)/mod(2);
y0(:,3) = y0(:,3)/mod(3);
lp = lp+log(abs(mod));
%三个Lyapunov指数
Lyapunov1(i) = lp(1)/(tstart);
Lyapunov2(i) = lp(2)/(tstart);
Lyapunov3(i) = lp(3)/(tstart);
y(4:12) = y0';
end
% 作Lyapunov指数谱图
i = 1:iteratetimes;
plot(i,Lyapunov1,i,Lyapunov2,i,Lyapunov3)
程序中用到的ThreeGS程序如下：
%G-S正交化
function A = ThreeGS(V) % V 为3*3向量
v1 = V(:,1);
v2 = V(:,2);
v3 = V(:,3);
a1 = zeros(3,1);
a2 = zeros(3,1);
a3 = zeros(3,1);
a1 = v1;
a2 = v2-((a1'*v2)/(a1'*a1))*a1;
a3 = v3-((a1'*v3)/(a1'*a1))*a1-((a2'*v3)/(a2'*a2))*a2;
A = [a1,a2,a3];
%计算得到的Rossler系统的LE为———— 0.063231 0.092635 -9.8924
Dec 10, 2009
국제특허출원(3차원 메쉬모델 워터마킹)
WIPO(World Intellectual Property Organization), URL:
http://www.wipo.it/ptcdb/ja/ia.jsp
국제공개번호: WO/2009/107982 국제출원번호: PCT/KR2009/000906
국제공개일::03.09.2009
국제출원일:25.02.2009
국제특허분류: H04N 5/913 (2006.01)
출원인1:MARKANY INC.
[KR/KR]; 10F Ssangnim Bldg. 151-11, Ssangnim-dong Jung-gu, Seoul
100-400 (KR) (All Except US).
출원인: CHOI, Jong-Uk [KR/KR]; (KR) (US Only).
출원인:CUI, Jizhe [KR/KR]; (KR) (US Only).
출원인:KIM, Jong-Weon [KR/KR]; (KR) (US Only).
발명자: CHOI, Jong-Uk; (KR). CUI, Jizhe; (KR). KIM, Jong-Weon; (KR).
대리인:KOREANA PATENT FIRM; Dong-Kyong Bldg. 824-19, Yoksam-dong
Gangnam-gu Seoul 135-080 (KR).
특허 우선권정보:10-2008-0016919 25.02.2008 KR
발명의 명칭: (EN) WATERMARKING METHOD FOR THREE-DIMENSIONAL MESH MODEL AND
APPARATUS THEREOF
영문요약:
(EN) Provided are a watermarking method and a watermarking apparatus
for a mesh model, applicable to a system which requires high
precision, such as a rapid prototyping system. A reference coordinate
system is set by using 1-ring values of vertexes of a
three-dimensional mesh model. The vertexes of the model are aligned on
the basis of the set reference coordinate system, and bit information
of a watermark in a bit string is inserted into each of the polygonal
faces of the three-dimensional mesh model in accordance with the
sequence of the alignment. The present invention is advantageous in
that the shape of the model is not deformed even when watermarked, and
therefore the present invention is useful in a rapid prototyping
system which requires high precision for the purpose of authenticating
data integrity. The present invention can be also used for the purpose
of marking contents without the need for storing, and for the purpose
of information hiding.
지정국가:
AE, AG, AL, AM, AO, AT, AU, AZ, BA, BB, BG, BH, BR, BW, BY, BZ, CA,
CH, CN, CO, CR, CU, CZ, DE, DK, DM, DO, DZ, EC, EE, EG, ES, FI, GB,
GD, GE, GH, GM, GT, HN, HR, HU, ID, IL, IN, IS, JP, KE, KG, KM, KN,
KP, KZ, LA, LC, LK, LR, LS, LT, LU, LY, MA, MD, ME, MG, MK, MN, MW,
MX, MY, MZ, NA, NG, NI, NO, NZ, OM, PG, PH, PL, PT, RO, RS, RU, SC,
SD, SE, SG, SK, SL, SM, ST, SV, SY, TJ, TM, TN, TR, TT, TZ, UA, UG,
US, UZ, VC, VN, ZA, ZM, ZW.
African Regional Intellectual Property Org. (ARIPO) (BW, GH, GM, KE,
LS, MW, MZ, NA, SD, SL, SZ, TZ, UG, ZM, ZW)
Eurasian Patent Organization (EAPO) (AM, AZ, BY, KG, KZ, MD, RU, TJ, TM)
European Patent Office (EPO) (AT, BE, BG, CH, CY, CZ, DE, DK, EE, ES,
FI, FR, GB, GR, HR, HU, IE, IS, IT, LT, LU, LV, MC, MK, MT, NL, NO,
PL, PT, RO, SE, SI, SK, TR)
African Intellectual Property Organization (OAPI) (BF, BJ, CF, CG, CI,
CM, GA, GN, GQ, GW, ML, MR, NE, SN, TD, TG).
국제공개에 사용한 언어: Korean (KO)
국제출원에 사용한 언어: Korean (KO)
참고 링크:
http://www.wipo.int/pctdb/ja/ia.jsp?ia=KR2009%2F000906&IA=KR2009000906&DISPLAY=STATUS
Dec 8, 2009
통계
Posts
