$[0,1]^4$的投影

hbghlyj

Rotating Tesseract
codetrain

1. // Daniel Shiffman
2. // http://youtube.com/thecodingtrain
3. // http://codingtra.in
4. // JavaScript transcription: Chuck England
6. // Coding Challenge #113: 4D Hypercube
7. // https://youtu.be/XE3YDVdQSPo
9. // Matrix Multiplication
10. // https://youtu.be/tzsgS19RRc8
12. let angle = 0;
14. let points = [];
16. function setup() {
17. let size = min(windowWidth, windowHeight);
18. createCanvas(size, size, WEBGL);
19. points[0] = new P4Vector(-1, -1, -1, 1);
20. points[1] = new P4Vector(1, -1, -1, 1);
21. points[2] = new P4Vector(1, 1, -1, 1);
22. points[3] = new P4Vector(-1, 1, -1, 1);
23. points[4] = new P4Vector(-1, -1, 1, 1);
24. points[5] = new P4Vector(1, -1, 1, 1);
25. points[6] = new P4Vector(1, 1, 1, 1);
26. points[7] = new P4Vector(-1, 1, 1, 1);
27. points[8] = new P4Vector(-1, -1, -1, -1);
28. points[9] = new P4Vector(1, -1, -1, -1);
29. points[10] = new P4Vector(1, 1, -1, -1);
30. points[11] = new P4Vector(-1, 1, -1, -1);
31. points[12] = new P4Vector(-1, -1, 1, -1);
32. points[13] = new P4Vector(1, -1, 1, -1);
33. points[14] = new P4Vector(1, 1, 1, -1);
34. points[15] = new P4Vector(-1, 1, 1, -1);
35. }
37. function draw() {
38. background(0);
39. rotateX(-PI / 2);
40. let projected3d = [];
42. for (let i = 0; i < points.length; i++) {
43. const v = points[i];
45. const rotationXY = [
46. [cos(angle), -sin(angle), 0, 0],
47. [sin(angle), cos(angle), 0, 0],
48. [0, 0, 1, 0],
49. [0, 0, 0, 1],
50. ];
52. const rotationZW = [
53. [1, 0, 0, 0],
54. [0, 1, 0, 0],
55. [0, 0, cos(angle), -sin(angle)],
56. [0, 0, sin(angle), cos(angle)]
57. ];
59. let rotated = matmul(rotationXY, v);
60. rotated = matmul(rotationZW, rotated);
62. let distance = 2;
63. let w = 1 / (distance - rotated.w);
65. const projection = [
66. [w, 0, 0, 0],
67. [0, w, 0, 0],
68. [0, 0, w, 0],
69. ];
71. let projected = matmul(projection, rotated);
72. projected.mult(width / 8);
73. projected3d[i] = projected;
75. stroke(255, 200);
76. strokeWeight(32);
77. noFill();
79. point(projected.x, projected.y, projected.z);
80. }
82. // Connecting
83. for (let i = 0; i < 4; i++) {
84. connect(0, i, (i + 1) % 4, projected3d);
85. connect(0, i + 4, ((i + 1) % 4) + 4, projected3d);
86. connect(0, i, i + 4, projected3d);
87. }
89. for (let i = 0; i < 4; i++) {
90. connect(8, i, (i + 1) % 4, projected3d);
91. connect(8, i + 4, ((i + 1) % 4) + 4, projected3d);
92. connect(8, i, i + 4, projected3d);
93. }
95. for (let i = 0; i < 8; i++) {
96. connect(0, i, i + 8, projected3d);
97. }
99. //angle = map(mouseX, 0, width, 0, TWO_PI);
100. angle += 0.02;
101. }
103. function connect(offset, i, j, points) {
104. strokeWeight(4);
105. stroke(255);
106. const a = points[i + offset];
107. const b = points[j + offset];
108. line(a.x, a.y, a.z, b.x, b.y, b.z);
109. }

Copy the Code

解释let w = 1 / (distance - rotated.w);这行：perspective projection distance是成像平面的w值，rotated.w是变换后的点的w值。把点到投影平面（这里是3D平面）的距离缩放成1，即1=攝像機到成像平面的w距离。

GeoGebra

Tesseract by KyeongYong Kim.ggb (21.11 KB, Downloads: 2) 2024-3-12 06:01 Upload Click on the file to download the attachment	Resource
Tesseract rotations by Juan Carlos Ponce Campuzano.ggb (50.92 KB, Downloads: 3) 2024-3-12 06:00 Upload Click on the file to download the attachment	Resource

Video.webm (311.02 KB, Downloads: 222)

2024-3-12 05:31 Upload

Click on the file to download the attachment

英壬画板

bilibili.com/video/BV1t7411B7hV（视频存档）

hbghlyj

英壬画板的作图方法非常简短：

在Oxy平面作$AC$和过A的AC垂线
$AC$、过A的AC垂线交椭圆于$E,F,E_1,F_1$
把$E,F,E_1,F_1$旋转45°得4个点$E',F',E_1',F_1'$， 再绕Ox转90°得4个点， 再绕Ox转90°得4个点， 再绕Ox转90°得4个点， 就得到了全部16个点。

当两垂直线绕$A$旋转时，$E,F,E_1,F_1$在椭圆上运动，那16个点就跟着运动。

问题：为什么可以这样做呢，为什么出来的16个点就是$[0,1]^4$在3维空间的投影？

另外可以发现1楼的Tesseract rotations by Juan Carlos Ponce Campuzano.ggb的作法是使用矩阵算出点的坐标，非常清晰；而 Tesseract by KyeongYong Kim.ggb 的作法是类似于上面的英壬画板的作图方法，画一个椭圆和两条垂直直线，它还多做了一些仿射。

hbghlyj

hbghlyj 发表于 2024-3-11 21:43
英壬画板的作图方法非常简短：
在Oxy平面作垂直线$AC\perp BC$
$AC$、$BC$交椭圆于$E,F,E_1,F_1$

感觉是錯的吧？$xw$-旋转后，那投影四边形的对角线不会保持垂直啊

hbghlyj

1楼解释let w = 1 / (distance - rotated.w);这行提到的distance值是成像平面的w值，而且它把攝像機到成像平面的w距离固定为1了，所以可以说distance值是成像平面的w值与攝像機到成像平面的w距离之比。
画个简图：

以点$(0,0,0,\text{distance})$为中心投影到（3维仿射子空间）$w=\text{distance}-1$上
例如：点$(0,1,-1,0)$投影到$(0,\frac1{\text{distance}},-\frac1{\text{distance}},\text{distance}-1)$.
因为$$(0,0,0,\text{distance})=k(0,1,-1,0)+(1-k)(0,\frac1{\text{distance}},-\frac1{\text{distance}},\text{distance}-1)$$其中$k=\frac1{1-\text{distance}}$.

hbghlyj

hbghlyj 发表于 2024-3-11 22:12
感觉是錯的吧？$xw$-旋转后，那投影四边形的对角线不会保持垂直啊

仔细计算发现真的是正确的！把1楼提到的distance值取成sqrt(2)就一直会保持垂直。

1. P2 = {{1}, {1}, {-1}, {1}}
2. P11 = {{-1}, {1}, {-1}, {-1}}
3. P10 = {{1}, {1}, {-1}, {-1}}
4. P3 = {{-1}, {1}, {-1}, {1}}
6. rotXW = {{Cos[angleXW], 0, 0, Sin[angleXW]}, {0, 1, 0, 0}, {0, 0, 1,
7.    0}, {-Sin[angleXW], 0, 0, Cos[angleXW]}}
9. Q2 = rotXW . P2
10. Q11 = rotXW . P11
11. Q10 = rotXW . P10
12. Q3 = rotXW . P3
14. W2 = 1/(distance - Q2[[4, 1]])
15. M2 = {{W2, 0, 0, 0}, {0, W2, 0, 0}, {0, 0, W2, 0}}
17. W11=1/(distance-Q11[[4,1]])
18. M11={{W11,0,0,0},{0,W11,0,0},{0,0,W11,0}}
20. W10 = 1/(distance - Q10[[4, 1]])
21. M10 = {{W10, 0, 0, 0}, {0, W10, 0, 0}, {0, 0, W10, 0}}
23. W3=1/(distance-Q3[[4,1]])
24. M3={{W3,0,0,0},{0,W3,0,0},{0,0,W3,0}}
26. MP2=M2.Q2
27. MP11=M11.Q11
28. MP10=M10.Q10
29. MP3=M3.Q3
31. A2 = Transpose[MP2][[1]]
32. A11=Transpose[MP11][[1]]
33. A10=Transpose[MP10][[1]]
34. A3=Transpose[MP3][[1]]

Copy the Code

然后验证对角线垂直：$A_2A_{11}\perp A_{10}A_3$，即内积为0

1. Dot[A2-A11,A10-A3]//FullSimplify

Copy the Code

输出\[\frac{8 \left(\text{distance}^2-2\right) \cos (2\ \text{angleXW})}{\cos (4\ \text{angleXW})+2\ \text{distance}^4-4\ \text{distance}^2+1}\]
所以，设$\text{distance}^2-2=0$，则当$\text{angleXW}$变动时，恒成立$A_2A_{11}\perp A_{10}A_3$.

hbghlyj

在1#的Tesseract rotations by Juan Carlos Ponce Campuzano.ggb 中，设置常数distance为sqrt(2)，可以验证当$\text{angleXW}$变动时，恒成立$A_2A_{11}\perp A_{10}A_3$，且$A_2A_{11}$与$A_{10}A_3$的交点$\left(0,\frac{1}{\text{distance}},-\frac{1}{\text{distance}}\right)$是不变的（它是$(0,1,-1,0)$的投影，而$(0,1,-1,0)$在xw旋转下不变），且$A_2,A_{11},A_{10},A_3$所在的2D平面是不变的.

hbghlyj

是不是与求射影变换将二次曲线映为圆且保持顶点在P的角度不变有关？
把外接圆投影成了一条二次曲线（看上去像抛物线），且保持顶点在正方形中心的角度不变？

hbghlyj

hbghlyj 发表于 2024-3-12 11:00
把外接圆投影成了一条二次曲线（看上去像抛物线）

验证：
取$\Bbb R^3$的平面$y+z=0$的正交基$(1,0,0),\frac1{\sqrt2}(0,1,-1)$. 把上面算出的A2表示为
$$A_2=a(1,0,0)+b\frac1{\sqrt2}(0,1,-1)$$
算出\begin{aligned}a&=\frac{\sin (\text{angleXW})+\cos (\text{angleXW})}{\sin (\text{angleXW})-\cos (\text{angleXW})+\text{distance}}\\b&=\frac{\sqrt{2}}{\sin (\text{angleXW})-\cos (\text{angleXW})+\text{distance}}\end{aligned}消去angleXW得\[2 a^2+b^2 \left(\text{distance}^2-2\right)-2 \sqrt{2} \text{distance}\ b+2=0\tag1\label1\]所以当$\text{distance}=\sqrt2$时，$A_2,A_{11},A_{10},A_3$共同的轨迹是一条抛物线
\[a^2-2 b+1=0\]
式\eqref{1}整理得\[\frac{a^2}{\frac{2}{\text{distance}^2-2}}+\frac{\left(b-\frac{\sqrt{2} \text{distance}}{\text{distance}^2-2}\right)^2}{\frac{4}{\left(\text{distance}^2-2\right)^2}}=1\]
当$0<\text{distance}<2$时，离心率为$\sqrt{2-\frac{\text{distance}^2}{2}}$；
当$\text{distance}>2$时，离心率为$\sqrt{\frac{\text{distance}^2-4}{\text{distance}^2-2}}$.
Plot[Piecewise[{{Sqrt[2 - x^2/2], x < 2}, {Sqrt[(-4 + x^2)/(-2 + x^2)], x > 2}}], {x, 0, 4}]

hbghlyj

hbghlyj 发表于 2024-3-12 10:35
且$A_2A_{11}$与$A_{10}A_3$的交点$\left(0,\frac{1}{\text{distance}},-\frac{1}{\text{distance}}\right)$是不变的

计算$A_2A_{11}$与$A_{10}A_3$的交点：

1. FullSimplify[
2. k A2 + (1 - k) A11 /.
3.   Solve[k A2[[1]] + (1 - k) A11[[1]] == 0, k][[1]]]

Copy the Code

$$\left(0,\frac{1}{\text{distance}},-\frac{1}{\text{distance}}\right)=0(1,0,0)+\frac{\sqrt{2}}{\text{distance}}\frac1{\sqrt2}(0,1,-1)$$
所以该点在上面取的平面$y+z=0$的正交基中的坐标为$(0,\frac{\sqrt{2}}{\text{distance}})$，当$\text{distance}=\sqrt{2}$时恰好为抛物线$a^2-2 b+1=0$的焦点$(0,1)$

1. ResourceFunction["ConicProperties"][(2+2 a^2-2 Sqrt[2] b distance+b^2 (-2+distance^2)/.distance->Sqrt[2])==0,{a,b}]["Focus"]

Copy the Code

import graph;

size(200);

real f(real x) {
    return (x^2 + 1)/2;
}

real xmax = 4.5;
real xmin = -xmax;
real ymin = -.5;
real ymax = f(xmax);

xaxis(Label("$a$"), xmin, xmax, above=true, Arrow);
yaxis(Label("$b$"), ymin, ymax, above=true, Arrow);

path g=graph(f, xmin, xmax);
pair P = (0, 1);
real angle1 = 60;
real angle2 = 150;
real[] x=intersections(g, P, P + dir(angle1));
pair A11=point(g,x[0]);
pair A2=point(g,x[1]);
x=intersections(g, P, P + dir(angle2));
pair A3=point(g,x[0]);
pair A10=point(g,x[1]);
draw(g);
dot("$A_2$", A2, NE);
dot("$A_{11}$", A11, SW);
dot("$A_{10}$", A10, SE);
dot("$A_3$", A3, NW);
draw(A2--A11, red);
draw(A10--A3, red);
dot((0,1));

hbghlyj

hbghlyj 发表于 2024-3-12 11:00
是不是与保持原点且保持顶点在原点的角度不变的射影变换有关？
把外接圆投影成了一条二次曲线（看上去 ...

透视投影限制在$P_2,P_3,P_{10},P_{11}$所在的平面上是一个射影变换，
取distance = Sqrt[2]我们把这个射影变换算出来：P是$P_2,P_3,P_{10},P_{11}$所在的平面上任意点，它被xw-旋转然后透视投影到$A_2,A_3,A_{10},A_{11}$所在的平面上的点A：

1. distance = Sqrt[2]
2. P = {{0}, {1}, {-1}, {0}} + {{u}, {0}, {0}, {v}}
3. Q = rotXW . P
4. W = 1/(distance - Q[[4, 1]])
5. M = {{W, 0, 0, 0}, {0, W, 0, 0}, {0, 0, W, 0}}
6. MP = M . Q
7. A = Transpose[MP][[1]]
8. FullSimplify[{a, b} /.
9.   Solve[A == a {1, 0, 0} + (b + 1)/Sqrt[2] {0, 1, -1}, {a, b}][[1]]]

Copy the Code

输出A的(a,b)坐标：\[\left\{\frac{u \cos (\text{angleXW})+v \sin (\text{angleXW})}{u \sin (\text{angleXW})-v \cos (\text{angleXW})+\sqrt{2}},\frac{\sqrt{2}}{u \sin (\text{angleXW})-v \cos (\text{angleXW})+\sqrt{2}}-1\right\}\]（这里的b是上面帖子中的b减1，平移坐标系是为了使这个变换把$u=0,v=0$的点映射到$a=0,b=0$的点，符合那帖的条件。）
写成矩阵
\[\left(
\begin{array}{ccc}
\cos (\text{angleXW}) & \sin (\text{angleXW}) & 0 \\
-\sin (\text{angleXW}) & \cos (\text{angleXW}) & 0 \\
\sin (\text{angleXW}) & -\cos (\text{angleXW}) & \sqrt{2} \\
\end{array}
\right)\]符合保持原点且保持顶点在原点的角度不变的射影变换的形式$\pmatrix{\cosα&-\sinα&0\\\sinα&\cosα&0\\A&B&C}$.

hbghlyj

hbghlyj 发表于 2024-3-11 21:43
英壬画板的作图方法非常简短：
在Oxy平面作$AC$和过A的AC垂线
$AC$、过A的AC垂线交椭圆于$E,F,E_1,F_1$

从上面的计算发现它这个作法还是不太对的，因为如果要保持直角应该是作一条抛物线，不是椭圆。正确作法是在适当的位置作一条抛物线，再过抛物线焦点作两条垂直线，交抛物线于四点，再关于x轴转90°,180°,270°，这样得到16个点是$[0,1]^4$关于xw面旋转后的透视投影。