貝の関数を作る

• 昨日はWikipediaやWolframの貝のパラメタ表示のメモをした
• 自分なりのそれも作りたい
• ムール貝風

phis <- seq(from=0,to=pi/8,length=100)
rs <- sin(seq(from=0,to=pi,length=length(phis)))
t <- seq(from=-50,to=-1,length=300)
k <- 0.3
theta <- 0.1
X <- matrix(0,0,3)

for(i in 1:length(phis)){
	rot <- matrix(c(cos(phis[i]),sin(phis[i]),-sin(phis[i]),cos(phis[i])),2,2)
	tmp <- cbind(rs[i]*exp(k*t)*cos(theta*t),rs[i]*exp(k*t)*sin(theta*t),rep(0,length(t)))
	tmp[,c(1,3)] <- t(rot %*% t(tmp[,c(1,3)]))
	X <- rbind(X,tmp)
}
X. <- rbind(X,rep(min(X),3))
X. <- rbind(X.,rep(max(X),3))
plot3d(X.)

• ホタテガイ風

# ホタテガイ風
phis <- seq(from=0,to=pi/1.1,length=100)
rs <- sin(seq(from=0,to=pi,length=length(phis)))
t <- seq(from=-50,to=-1,length=300)
k <- 0.1
theta <- 0.05
X <- matrix(0,0,3)

for(i in 1:length(phis)){
	rot <- matrix(c(cos(phis[i]),sin(phis[i]),-sin(phis[i]),cos(phis[i])),2,2)
	tmp <- cbind(rs[i]*exp(k*t)*cos(theta*t),rs[i]*exp(k*t)*sin(theta*t),rep(0,length(t)))
	tmp[,c(1,3)] <- t(rot %*% t(tmp[,c(1,3)]))
	X <- rbind(X,tmp)
}
X. <- rbind(X,rep(min(X),3))
X. <- rbind(X.,rep(max(X),3))
plot3d(X.)

• オウムガイ風

# オウムガイ風
phis <- seq(from=-pi/1,to=pi/1,length=30)
rs <- sin(seq(from=pi/3,to=pi-pi/3,length=length(phis)))
t <- seq(from=-30,to=2,length=400)
k <- 0.1
kz <- 0.00001

theta <- 1
X <- matrix(0,0,3)

for(i in 1:length(phis)){
	#rot <- matrix(c(cos(phis[i]),sin(phis[i]),-sin(phis[i]),cos(phis[i])),2,2)
	tmp.r <- exp(k*t)
	tmp.R <- exp(k*t)*(1-exp(-k*2*pi))
	tmp <- cbind((tmp.r+cos(phis[i])*tmp.R)*cos(theta*t),(tmp.r+cos(phis[i])*tmp.R)*sin(theta*t),rep(sin(phis[i])*tmp.R,length(t)))
	#tmp[,c(1,3)] <- t(rot %*% t(tmp[,c(1,3)]))
	X <- rbind(X,tmp)
}
X. <- rbind(X,rep(min(X),3))
X. <- rbind(X.,rep(max(X),3))
plot3d(X.)