トーラス上の格子 三角形

  • こちらで、正方格子座標をトーラス上に描いた
  • 今日は三角格子

http://www.genome.med.kyoto-u.ac.jp/StatGenet/testRY20110208/sankakuzahyou.png

  • まず、平面に三角格子を描き
  • 3方向の平行線をトーラス上に描き
  • 最後に3方向の平行線のすべてをトーラス上に色分けして描く
  • 格子を描くのは、平行線の方向ベクトルとその法線ベクトルの格子を作って、type="l"で描くと描ける
  • 3方向あるときは、その3方向を与えればよい
  • 正三角形のときは、pi/3ずつ回す
  • サムネイル画像を見ると、「吸い込まれ」ている様子が見える
  • さて、平面上の直線で作られた格子をトーラス上に描いた
  • 平面上の直線は「まっすぐな道」
  • では、トーラス上に作った世界で「まっすぐに歩く」とどこを歩くだろう。こんな感じ

http://www.genome.med.kyoto-u.ac.jp/StatGenet/testRY20110208/massugu.png

library(rgl)
# 三角形を描こう
gr<-51 

initangle<-runif(1)*2*pi
v1<-c(sin(initangle),cos(initangle))
u1<-c(sin(initangle+pi/2),cos(initangle+pi/2))
v2<-c(sin(initangle+pi/3),cos(initangle+pi/3))
u2<-c(sin(initangle+pi/3+pi/2),cos(initangle+pi/3+pi/2))
v3<-c(sin(initangle+2*pi/3),cos(initangle+2*pi/3))
u3<-c(sin(initangle+2*pi/3+pi/2),cos(initangle+2*pi/3+pi/2))

N1<-N2<-seq(from=-1,to=1,length=gr)
Ns<-expand.grid(N1,N2)

x1<-v1[1]*Ns[,1]+u1[1]*Ns[,2]
y1<-v1[2]*Ns[,1]+u1[2]*Ns[,2]

tmp<-cbind(x1,y1)

xysq1<-tmp

xyseq<-tmp
x1<-v2[1]*Ns[,1]+u2[1]*Ns[,2]
y1<-v2[2]*Ns[,1]+u2[2]*Ns[,2]

tmp<-cbind(x1,y1)
xysq2<-tmp
xyseq<-rbind(xyseq,tmp)
tmp<-cbind(x1,y1)
#xyseq<-tmp
x1<-v3[1]*Ns[,1]+u3[1]*Ns[,2]
y1<-v3[2]*Ns[,1]+u3[2]*Ns[,2]

tmp<-cbind(x1,y1)
xysq3<-tmp
xyseq<-rbind(xyseq,tmp)

plot(xyseq[,1],xyseq[,2],type="l")

R1<-1/pi
R2<-1/pi

# 第一座標
Xsq<-R1*sin(2*pi*xysq1[,1])+R2*sin(2*pi*xysq1[,2])*sin(2*pi*xysq1[,1])
Ysq<-R1*cos(2*pi*xysq1[,1])+R2*sin(2*pi*xysq1[,2])*cos(2*pi*xysq1[,1])
Zsq<-R2*cos(2*pi*xysq1[,2])
plot3d(Xsq,Ysq,Zsq,type="l")

# 第二座標
Xsq2<-R1*sin(2*pi*xysq2[,1])+R2*sin(2*pi*xysq2[,2])*sin(2*pi*xysq2[,1])
Ysq2<-R1*cos(2*pi*xysq2[,1])+R2*sin(2*pi*xysq2[,2])*cos(2*pi*xysq2[,1])
Zsq2<-R2*cos(2*pi*xysq2[,2])
open3d()
plot3d(Xsq2,Ysq2,Zsq2,type="l")

# 第三座標
Xsq3<-R1*sin(2*pi*xysq3[,1])+R2*sin(2*pi*xysq3[,2])*sin(2*pi*xysq3[,1])
Ysq3<-R1*cos(2*pi*xysq3[,1])+R2*sin(2*pi*xysq3[,2])*cos(2*pi*xysq3[,1])
Zsq3<-R2*cos(2*pi*xysq3[,2])
open3d()
plot3d(Xsq3,Ysq3,Zsq3,type="l")

# 連結する
XsqA<-c(Xsq,Xsq2,Xsq3)
YsqA<-c(Ysq,Ysq2,Ysq3)
ZsqA<-c(Zsq,Zsq2,Zsq3)
# 経線と緯線の色を塗り分ける
open3d()
plot3d(XsqA,YsqA,ZsqA,type="l",col=c(rep(2,length(xysq[,1])),rep(3,length(xysq[,1])),rep(4,length(xysq[,1]))))
  • まっすぐな道のソース
library(rgl)
# 三角形を描こう
gr<-101 

initangle<-runif(1)*2*pi
#initangle<-pi/2
v1<-c(sin(initangle),cos(initangle))
u1<-c(sin(initangle+pi/2),cos(initangle+pi/2))
v2<-c(sin(initangle+pi/3),cos(initangle+pi/3))
u2<-c(sin(initangle+pi/3+pi/2),cos(initangle+pi/3+pi/2))
v3<-c(sin(initangle+2*pi/3),cos(initangle+2*pi/3))
u3<-c(sin(initangle+2*pi/3+pi/2),cos(initangle+2*pi/3+pi/2))

N1<-N2<-seq(from=-3,to=3,length=gr)
Ns<-expand.grid(N1,N2)

x1<-v1[1]*Ns[,1]+u1[1]*Ns[,2]
y1<-v1[2]*Ns[,1]+u1[2]*Ns[,2]

tmp<-cbind(x1,y1)

xysq1<-tmp

xyseq<-tmp
x1<-v2[1]*Ns[,1]+u2[1]*Ns[,2]
y1<-v2[2]*Ns[,1]+u2[2]*Ns[,2]

tmp<-cbind(x1,y1)
xysq2<-tmp
xyseq<-rbind(xyseq,tmp)
tmp<-cbind(x1,y1)
#xyseq<-tmp
x1<-v3[1]*Ns[,1]+u3[1]*Ns[,2]
y1<-v3[2]*Ns[,1]+u3[2]*Ns[,2]

tmp<-cbind(x1,y1)
xysq3<-tmp
xyseq<-rbind(xyseq,tmp)

#plot(xyseq[,1],xyseq[,2],type="l")

R1<-1/pi
R2<-1/pi

# 第一座標
Xsq<-R1*sin(2*pi*xysq1[,1])+R2*sin(2*pi*xysq1[,2])*sin(2*pi*xysq1[,1])
Ysq<-R1*cos(2*pi*xysq1[,1])+R2*sin(2*pi*xysq1[,2])*cos(2*pi*xysq1[,1])
Zsq<-R2*cos(2*pi*xysq1[,2])
#plot3d(Xsq,Ysq,Zsq,type="l")

# 第二座標
Xsq2<-R1*sin(2*pi*xysq2[,1])+R2*sin(2*pi*xysq2[,2])*sin(2*pi*xysq2[,1])
Ysq2<-R1*cos(2*pi*xysq2[,1])+R2*sin(2*pi*xysq2[,2])*cos(2*pi*xysq2[,1])
Zsq2<-R2*cos(2*pi*xysq2[,2])
#open3d()
#plot3d(Xsq2,Ysq2,Zsq2,type="l")

# 第三座標
Xsq3<-R1*sin(2*pi*xysq3[,1])+R2*sin(2*pi*xysq3[,2])*sin(2*pi*xysq3[,1])
Ysq3<-R1*cos(2*pi*xysq3[,1])+R2*sin(2*pi*xysq3[,2])*cos(2*pi*xysq3[,1])
Zsq3<-R2*cos(2*pi*xysq3[,2])
#open3d()
#plot3d(Xsq3,Ysq3,Zsq3,type="l")

# 連結する
XsqA<-c(Xsq,Xsq2,Xsq3)
YsqA<-c(Ysq,Ysq2,Ysq3)
ZsqA<-c(Zsq,Zsq2,Zsq3)
# 経線と緯線の色を塗り分ける
#open3d()
#plot3d(XsqA,YsqA,ZsqA,type="l",col=c(rep(2,length(Xsq)),rep(3,length(Xsq)),rep(4,length(Xsq)-gr),rep(1,gr)))
#open3d()
plot3d(XsqA[(length(Xsq)*3-gr+1):(length(Xsq)*3)],YsqA[(length(Xsq)*3-gr+1):(length(Xsq)*3)],ZsqA[(length(Xsq)*3-gr+1):(length(Xsq)*3)],type="l")
open3d()
plot3d(XsqA[1:gr],YsqA[1:gr],ZsqA[1:gr],type="l")