本来想做绕轴旋转,绕轴旋转是在垂直于轴向量的空间里旋转,但是n维空间里与某个向量垂直的空间为n-1维,而旋转只在二维空间里有定义。所以这里就改成了,向任意方向旋转。

计算单位向量 X = ( x 1 , x 2 , ⋯ , x n ) X=(x_1,x_2,\cdots,x_n) X=(x1​,x2​,⋯,xn​)旋转到单位向量 V = ( v 1 , v 2 , ⋯ , v n ) V=(v_1,v_2,\cdots,v_n) V=(v1​,v2​,⋯,vn​)的旋转矩阵 R R R。

旋转坐标系

任意的旋转都可以看作绕着一个轴,在某个平面上的旋转。
不失一般性,假定向量 V = ( v 1 , v 2 , ⋯ , v n ) V=(v_1,v_2,\cdots,v_n) V=(v1​,v2​,⋯,vn​)在 v 1 × v 2 v_1\times v_2 v1​×v2​张成的平面上旋转,用矩阵乘法可表示为: V ′ = R 1 V V'=R_1V V′=R1​V其中旋转矩阵 R 1 R_1 R1​定义为
R 1 = ( cos ⁡ α − sin ⁡ α sin ⁡ α cos ⁡ α 1 ⋱ 1 ) R_1=\left( \begin{array}{ccc} \cos\alpha& -\sin\alpha& & & \\ \sin\alpha& \cos\alpha& & & \\ & &1 & & \\ & & & \ddots & \\ & & & & 1 \\ \end{array} \right) R1​=⎝⎜⎜⎜⎜⎛​cosαsinα​−sinαcosα​1​⋱​1​⎠⎟⎟⎟⎟⎞​
其中 cos ⁡ α = v 2 v 1 2 + v 2 2 \cos\alpha=\frac{v_2}{\sqrt{v_1^2+v_2^2}} cosα=v12​+v22​ ​v2​​, sin ⁡ α = v 1 v 1 2 + v 2 2 \sin\alpha=\frac{v_1}{\sqrt{v_1^2+v_2^2}} sinα=v12​+v22​ ​v1​​。
由此我们将向量 V V V旋转到与 ( v 1 , 0 , ⋯ , 0 ) (v_1,0,\cdots,0) (v1​,0,⋯,0)正交。得 V ′ = ( 0 , v 1 2 + v 2 2 , v 3 , v 4 , ⋯ , v n ) V'=(0,\sqrt{v_1^2+v_2^2}, v_3,v_4,\cdots,v_n) V′=(0,v12​+v22​ ​,v3​,v4​,⋯,vn​)。
以此类推,将 V V V旋转到 ( 0 , ⋯ , 0 , 1 ) (0,\cdots,0,1) (0,⋯,0,1)需要 n − 1 n-1 n−1次旋转。将这些旋转变换组合为一个旋转矩阵 R 0 = R n − 1 R n − 2 ⋯ R 1 R_0=R_{n-1}R_{n-2}\cdots R_1 R0​=Rn−1​Rn−2​⋯R1​(注意,乘法不满足交换律,顺序不能颠倒)。顺便我们可以计算出这个旋转的逆矩阵 R 0 − 1 = R 1 ⊤ ⋯ R n − 2 ⊤ R n − 1 ⊤ R_0^{-1}=R_1^\top\cdots R_{n-2}^\top R_{n-1}^\top R0−1​=R1⊤​⋯Rn−2⊤​Rn−1⊤​。

旋转向量

现在我们要做的是将 R 0 X R_0X R0​X旋转到坐标轴 ( 0 , ⋯ , 0 , 1 ) (0,\cdots,0,1) (0,⋯,0,1),并求旋转矩阵。做法与上面相同。最后可以求得,在旋转后的坐标系下,旋转矩阵为 R X R_X RX​: R 0 V = R X R 0 X R_0V=R_XR_0X R0​V=RX​R0​X

虽然 R R R和 R X R_X RX​表示的旋转角度相同,但是这两个旋转所在的坐标系是不同的,还需要将 R X R_X RX​的坐标轴旋转回去:
R = R 0 − 1 R X R 0 R=R_0^{-1}R_XR_0 R=R0−1​RX​R0​

代码实现

import numpy as npdef Rot_map(V):assert len(V.shape) == 1assert np.linalg.norm(V) - 1 < 1e-8n_dim = V.shape[0]Rot = np.eye(n_dim)Rot_inv = np.eye(n_dim)for rotate in range(n_dim-1):rot_mat = np.eye(n_dim)rot_norm = np.sqrt(V[rotate]**2 + V[rotate+1]**2)cos_theta = V[rotate+1]/rot_normsin_theta = V[rotate]/rot_normrot_mat[rotate,rotate] = cos_thetarot_mat[rotate,rotate+1] = - sin_thetarot_mat[rotate+1,rotate] = sin_thetarot_mat[rotate+1,rotate+1] = cos_thetaV = np.dot(rot_mat, V)Rot = np.dot(rot_mat, Rot)Rot_inv = np.dot(Rot_inv,rot_mat.transpose())return Rot, Rot_invn_dim = 512X = np.random.rand(n_dim)
X = X / np.linalg.norm(X)
V = np.random.rand(n_dim)
V = V / np.linalg.norm(V)R_0, R_0_inv = Rot_map(V)
R_X, _ = Rot_map(np.dot(R_0, X))R = np.dot(np.dot(R_0_inv, R_X), R_0)assert np.linalg.norm(np.dot(R,X) - V) < 1e-8

扩展

这个方法可以用来在正交于某个向量的空间中随机采样。例如,在 n n n维空间中,采样与 V = ( v 1 , v 2 , ⋯ , v n ) V=(v_1,v_2,\cdots,v_n) V=(v1​,v2​,⋯,vn​)正交的向量,则代码如下:

import numpy as npdef Rot_map(V):assert len(V.shape) == 1assert np.linalg.norm(V) - 1 < 1e-8z_dim = V.shape[0]Rot = np.eye(z_dim)Rot_inv = np.eye(z_dim)for rotate in range(z_dim-1):rot_mat = np.eye(z_dim)rot_norm = np.sqrt(V[rotate]**2 + V[rotate+1]**2)cos_theta = V[rotate+1]/rot_normsin_theta = V[rotate]/rot_normrot_mat[rotate,rotate] = cos_thetarot_mat[rotate,rotate+1] = - sin_thetarot_mat[rotate+1,rotate] = sin_thetarot_mat[rotate+1,rotate+1] = cos_thetaV = np.dot(rot_mat, V)Rot = np.dot(rot_mat, Rot)Rot_inv = np.dot(Rot_inv,rot_mat.transpose())return Rot, Rot_invz_dim = 512
n_samples = 10000V = np.random.rand(z_dim)
V = V / np.linalg.norm(V)R_0, R_0_inv = Rot_map(V)samples = np.random.rand(z_dim-1, n_samples)
samples = np.concatenate((samples, np.zeros((1, n_samples))))
samples = np.dot(R_0_inv, samples)assert np.mean(np.abs(np.dot(V, samples))) < 1e-8

进一步,我们采样与一组向量 V s Vs Vs正交的向量:

import numpy as npdef gs(X, row_vecs=True, norm = True):'''https://gist.github.com/iizukak/1287876/edad3c337844fac34f7e56ec09f9cb27d4907cc7'''if not row_vecs:X = X.TY = X[0:1,:].copy()for i in range(1, X.shape[0]):proj = np.diag((X[i,:].dot(Y.T)/np.linalg.norm(Y,axis=1)**2).flat).dot(Y)Y = np.vstack((Y, X[i,:] - proj.sum(0)))if norm:Y = np.diag(1/np.linalg.norm(Y,axis=1)).dot(Y)if row_vecs:return Yelse:return Y.Tdef Rot_map(V):assert len(V.shape) == 1assert np.linalg.norm(V) - 1 < 1e-8z_dim = V.shape[0]Rot = np.eye(z_dim)Rot_inv = np.eye(z_dim)for rotate in range(z_dim-1):rot_mat = np.eye(z_dim)rot_norm = np.sqrt(V[rotate]**2 + V[rotate+1]**2)cos_theta = V[rotate+1]/rot_normsin_theta = V[rotate]/rot_normrot_mat[rotate,rotate] = cos_thetarot_mat[rotate,rotate+1] = - sin_thetarot_mat[rotate+1,rotate] = sin_thetarot_mat[rotate+1,rotate+1] = cos_thetaV = np.dot(rot_mat, V)Rot = np.dot(rot_mat, Rot)Rot_inv = np.dot(Rot_inv,rot_mat.transpose())return Rot, Rot_invz_dim = 10 # number of dimensions
n_vectors= 2 # number of Vs
assert n_vectors < z_dim
n_samples = 10 # number of samples to orthogonal VsVs = np.random.rand(n_vectors, z_dim)
# Gram-Schmidt Orthogonization
orth_Vs = gs(Vs)
orth_Vs = orth_Vs.TRots = np.eye(z_dim)
Rot_invs = np.eye(z_dim)
for idx in range(n_vectors):V = orth_Vs[:, idx]R_0 = np.eye(z_dim)R_0_inv = np.eye(z_dim)R_0[0:z_dim-idx,0:z_dim-idx], R_0_inv[0:z_dim-idx,0:z_dim-idx] = Rot_map(V[range(z_dim-idx)])orth_Vs = R_0.dot(orth_Vs)Rots = np.dot(R_0, Rots)Rot_invs = np.dot(Rot_invs, R_0_inv)samples = np.random.rand(z_dim-n_vectors, n_samples)
samples = np.concatenate((samples, np.zeros((n_vectors, n_samples))))
samples = np.dot(Rot_invs, samples)assert np.mean(np.abs(np.dot(Vs, samples))) < 1e-8

用矩阵旋转太慢了,上述采样操作可以用矩阵分解解决(在任意向量的正交空间中采样):

import numpy as npz_dim = 10 # number of dimensions
n_vectors= 2 # number of Vs
assert n_vectors < z_dim
n_samples = 20 # number of samples to orthogonal VsVs = np.random.rand(n_vectors, z_dim)
Vs = np.concatenate((Vs,np.zeros((z_dim-n_vectors, z_dim))))U, s, V = np.linalg.svd(Vs)directions = np.eye(z_dim)
directions[range(n_vectors),range(n_vectors)] = 0
directions = np.dot(np.dot(U, directions), V)samples = np.random.rand(z_dim, n_samples)
samples = np.dot(directions.T, samples)assert np.mean(np.abs(np.dot(Vs, samples))) < 1e-8

python n维向量向任意方向旋转,求旋转矩阵相关推荐

  1. python二维向量运算模拟_python二维向量运算_[VB.NET][C#]二维向量的基本运算

    前言 在数学中,几何向量是指具有大小和方向的几何对象. 在编程中,向量有着广泛的应用,其作用在图形编程和游戏物理引擎方面尤为突出. 第一节 构造函数 通过创建一个二维向量的类(或结构体),实现向量的表 ...

  2. python二维向量公式坐标叉乘_Python 或 sage 有没有哪个函数可以求指定向量的模?进行点乘运算,叉乘运算。...

    Hlianbobo:python 或 sage 有没有哪个函数可以求指定向量的模?进行点乘运算,叉乘运算?求相关的库名称已经代码示例,谢谢! zhzy:numpy Tony042:numpy,或者自己 ...

  3. python多维向量聚类_机器学习:Python实现聚类算法(三)之总结

    考虑到学习知识的顺序及效率问题,所以后续的几种聚类方法不再详细讲解原理,也不再写python实现的源代码,只介绍下算法的基本思路,使大家对每种算法有个直观的印象,从而可以更好的理解函数中参数的意义及作 ...

  4. python二维向量运算模拟_Python数学基础之向量定义与向量运算(附代码)

    患难与困苦是磨练人格的最高学府.--苏格拉底(公元前470年-公元前399年) Adversity and pinch are the highest institution of higher le ...

  5. python二维向量运算模拟_【转载 Python】Numpy基础:数组和矢量运算

    Numpy基础:数组和矢量运算 目录: 尊重原创,转载出处:https://blog.csdn.net/sunyaowu315/article/details/82733249?utm_medium= ...

  6. python二维向量运算_python中二维数组的Elementwise与or或运算

    numpy有logical_or.logical_xor和{},它们有一个reduce方法>> np.logical_or.reduce(a, axis=0) array([ True, ...

  7. python二维向量运算,二维数组python中的矢量化计数

    像这样的东西呢,它是矢量化的,没有for循环:def moving_count(a, value, axis=0): """Return sequential count ...

  8. python运维是做什么的_Python运维工程师主要干什么?

    Python运维工程师主要干什么? 更新时间:2020年05月09日 作者:spoto 广义运维工程师的职责介绍 Python运维工程师主要干什么?其实无论是做什么运维,运维工程师最基本的职责都是负责 ...

  9. 从某一点出发沿任意一方向旋转矩阵计算思考与实现

    欢迎关注更多精彩 关注我,学习常用算法与数据结构,一题多解,降维打击. 上期讲到 绕任一向量旋转矩阵计算思考与实现 点击前往 点击前往 问题提出 之前讲到绕任一向量旋转矩阵实现,原来的向量都是从原点出 ...

最新文章

  1. 每个公司都会用的短 URL 服务,怎么设计与实现?
  2. Transaction And Lock--事务中使用return会回滚事务吗?
  3. 归并排序执行次数_肯定能懂的常见算法讲解(1)——排序算法
  4. 思科Webex与下一代视频会议
  5. cognos report在做同比时遇到的问题解决方法
  6. Spring中的计时器StopWatch
  7. linux mariadb 升级,linux mariadb
  8. 【飞秋】使用C++语言创建Silverlight中的VisualState
  9. 2020计算机大纲,计算机专业2020考试大纲.doc
  10. java有用还是c有用_这周的有用资源
  11. 苹果真伪查询_拆解报告:山寨版苹果AirPods Pro
  12. unix bsd linux shell bash GNU之间的联系,歪讲Linux(一)
  13. scratch英语计算机,scratch 2怎么做计算器
  14. android rom打包解包工具,Android刷机包解包打包
  15. LeetCode第9题 回文数(Palindrome Number)
  16. 提升网站关键词的排名
  17. 掌握盲打?写个脚本陪自己练,每天练一遍,基本上一个星期就成了!
  18. 如何成为优秀的数据人经验分享
  19. 软件开发的43款可视化工具
  20. 计算机过程控制系统教材,计算机过程控制系统

热门文章

  1. 当鼠标移动到表单上,改变表单的颜色
  2. FTP测试手机软件图标素材可爱,手把手教你测试FTP
  3. 反思中国十大失败企业
  4. 如何修改wincc服务器画面,关于OS站的wincc画面修改的问题-工业支持中心-西门子中国...
  5. 全志T7/T507 Qt5.12.5移植记录
  6. 逻辑思维强的人适合学计算机不,逻辑思维强的人适合做什么工作?
  7. 编写dylib_iOS逆向-无需越狱注入动态库
  8. 前端切页快速上手系列
  9. AnyRTC将携互动直播连线2016杭州·云栖大会
  10. Spring cloud 多种限流方案