给定一个任意单位向量，计算任意正交单位向量的最佳方法是啥？
Posted
技术标签:
【中文标题】给定一个任意单位向量，计算任意正交单位向量的最佳方法是啥？【英文标题】：Given a single arbitrary unit vector, what is the best method to compute an arbitrary orthogonal unit vector?给定一个任意单位向量，计算任意正交单位向量的最佳方法是什么？ 【发布时间】：2013-11-08 02:10:13 【问题描述】：here 基本上问了同样的问题，但在非编程环境中。建议的解决方案是采用 y, -x, 0 。这适用于所有具有 x 或 y 分量的向量，但如果向量等于 + 或 - 0, 0, 1 ，则会失败。在这种情况下，我们会得到 0, 0, 0 。
我目前的解决方案（c++）：
// floating point comparison utilizing epsilon
bool is_equal(float, float);
// ...
vec3 v = /* some unit length vector */
// ...
// Set as a non-parallel vector which we will use to find the
// orthogonal vector. Here we choose either the x or y axis.
vec3 orthog;
if( is_equal(v.x, 1.0f) )
orthog.set(1.0f, 0.0f, 0.0f);
else
orthog.set(0.0f, 1.0f, 0.0f);
// Find orthogonal vector
orthog = cross(v, orthog);
orthog.normalize();
此方法有效，但我觉得可能有更好的方法，我的搜索结果仅此而已。
[编辑]
只是为了好玩，我在 c++ 中对每个建议的答案的简单实现做了一个快速代码，并验证它们都有效（尽管有些并不总是自然地返回单位向量，我在需要的地方添加了一个 noramlize() 调用)。
我最初的想法：
vec3 orthog_peter(vec3 const& v)
vec3 arbitrary_non_parallel_vec = v.x != 1.0f ? vec3(1.0, 0.0f, 0.0f) : vec3(0.0f, 1.0f, 0.0f);
vec3 orthog = cross(v, arbitrary_non_parallel_vec);
return normalize( orthog );
https://***.com/a/19650362/2507444
vec3 orthog_robert(vec3 const& v)
vec3 orthog;
if(v.x == 0.0f && v.y == 0.0f)
orthog = vec3(1.0f, 1.0f, 0.0f);
else if(v.x == 0.0f)
orthog = vec3(1.0f, v.z / v.y, 1.0f);
else if(v.y == 0.0f)
orthog = vec3(-v.z / v.x, 1.0f, 1.0f);
else
orthog = vec3(-(v.z + v.y) / v.x, 1.0f, 1.0f);
return normalize(orthog);
https://***.com/a/19651668/2507444
// NOTE: u and v variable names are swapped from author's example
vec3 orthog_abhishek(vec3 const& v)
vec3 u(1.0f, 0.0f, 0.0f);
float u_dot_v = dot(u, v);
if(abs(u_dot_v) != 1.0f)
return normalize(u + (v * -u_dot_v));
else
return vec3(0.0f, 1.0f, 0.0f);
https://***.com/a/19658055/2507444
vec3 orthog_dmuir(vec3 const& v)
float length = hypotf( v.x, hypotf(v.y, v.z));
float dir_scalar = (v.x > 0.0) ? length : -length;
float xt = v.x + dir_scalar;
float dot = -v.y / (dir_scalar * xt);
return vec3(
dot * xt,
1.0f + dot * v.y,
dot * v.z);
;
【问题讨论】：
【参考方案1】：嗯，这是一种解决方法。给定一个向量 (a, b, c)。求解方程 (a, b, c) dot (aa, bb, cc) = 0 对于 aa、bb 和 cc（并确保 aa、bb 和 cc 不全为零），所以 (aa, bb, cc ) 与 (a, b, c) 正交。我已经使用 Maxima (http://maxima.sf.net) 来解决它。
(%i42) solve ([a, b, c] . [aa, bb, cc] = 0, [aa, bb, cc]), a=0, b=0;
(%o42) [[aa = %r19, bb = %r20, cc = 0]]
(%i43) solve ([a, b, c] . [aa, bb, cc] = 0, [aa, bb, cc]), a=0;
%r21 c
(%o43) [[aa = %r22, bb = - ------, cc = %r21]]
b
(%i44) solve ([a, b, c] . [aa, bb, cc] = 0, [aa, bb, cc]), b=0;
%r23 c
(%o44) [[aa = - ------, bb = %r24, cc = %r23]]
a
(%i45) solve ([a, b, c] . [aa, bb, cc] = 0, [aa, bb, cc]);
%r25 c + %r26 b
(%o45) [[aa = - ---------------, bb = %r26, cc = %r25]]
a
请注意，我首先解决了特殊情况（a = 0 和 b = 0，或 a = 0，或 b = 0），因为找到的解决方案并非对某些等于零的组件都有效。出现的 %r 变量是任意常量。我将它们设置为 1 以获得一些特定的解决方案。
(%i52) subst ([%r19 = 1, %r20 = 1], %o42);
(%o52) [[aa = 1, bb = 1, cc = 0]]
(%i53) subst ([%r21 = 1, %r22 = 1], %o43);
c
(%o53) [[aa = 1, bb = - -, cc = 1]]
b
(%i54) subst ([%r23 = 1, %r24 = 1], %o44);
c
(%o54) [[aa = - -, bb = 1, cc = 1]]
a
(%i55) subst ([%r25 = 1, %r26 = 1], %o45);
c + b
(%o55) [[aa = - -----, bb = 1, cc = 1]]
a
希望这会有所帮助。祝你好运，继续努力。
【讨论】：
我对 Maxima 不熟悉，但是从我从输出和您的笔记中收集到的信息来看，您是否还在求解一组非常具体的输入向量（其中 a = 0 和 b = 0)?我不是 100% 确定您所说的“请注意，我首先解决了特殊情况......因为找到的解决方案并非对某些等于零的组件都有效”。我需要求解任意输入单位向量（标题已更新以说明这一点），而这似乎仍然只是求解一个子集。 在这种情况下得到的解决方案中，如果 a = 0，则除以零。例如，如果输入向量是 0, 0, 1 ，你会得到 -(1/0), 1, 1 这不是一个有效的解决方案。可能是我误会了。 @PeterClark 您是正确的，给出的解决方案并非对所有输入都有效。每一个都涵盖了一组不同的条件，因此它们一起涵盖了所有输入。您会注意到其他受访者给出了类似的个案解决方案——例如“如果点 u 已经位于轴上，那么只需为点 v 选择任何其他轴”。我已经组织了解决方案，以便它们从最窄的特殊情况到最一般的情况。到了一般情况时，您已经排除了一般解决方案无效的特殊情况。 @PeterClark 当您有特殊情况要处理时，在这个问题或任何其他问题中，将它们从最具体到最一般组织起来很有用，例如if (special case 1) ... else if (special case 2) ... else /* 最一般的情况 */ .... 当你满足条件时，你排除了一个又一个的特殊情况，所以由当你走到最后时，你已经排除了所有的特殊情况，因此你处于最一般的情况。 哦，好吧，不是点击代码的第二个sn-p是4个唯一的解决方案。现在说得通了，谢谢！【参考方案2】：你需要选择一个不等于 0 且不在与给定单位向量 u 连接原点的线上的点 v。
如前所述，您可以在任何轴上选择一个单位向量，只要该点满足上述条件即可。如果点 u 已经位于轴上，则只需为点 v 选择任何其他轴。
然后你需要解方程(v + tu).u = 0。 （只需求解 t）
当然，您需要对其进行归一化以获得正交单位向量。
【讨论】：
为什么我们求解**(v + tu)**.u = 0 而不是**v.u** = 0？当你说我们把 v 作为一个点，不等于零，而不是在连接原点和 u 的线上（这意味着如果我们将 v 视为从原点开始的向量，它不是 u 的某个倍数） - 不v+tu 只需给出一条线的参数方程，其中 v 是线上的一个点，u 是线的方向（因此它平行于 u，这个方程永远不会正确）？
@PeterClark 查看我的更新答案。我已经上传了图片。【参考方案3】：
另一种方法是使用Householder reflectors。
我们可以找到一个反射器 Q 将我们的向量映射到 (1,0,0) 的倍数。将 Q 应用于 (0,1,0) 将给出一个垂直于我们的向量的向量。这种方法的一个优点是它适用于任意数量的维度。另一个是我们可以 获得垂直于原始向量和新向量的其他向量：将 Q 应用于 (0,0,1)。听起来可能很复杂，但这里是 3d 的 C 代码（xp,yp,zp 是必需的向量，长度为 1；正如所写，所有内容都是双精度数，但您可以使用 float 并使用 hypotf 而不是 hypot）
l = hypot( x, hypot(y,z));
s = (x > 0.0) ? l : -l;
xt = x + s;
dot = -y/(s*xt);
xp = dot*xt;
yp = 1.0 + dot*y;
zp = dot*z;
【讨论】：
也许 l = sqrt(x * x + y * y + z * z)，节省一个平方根，一个乘数和两个函数调用。【参考方案4】：这是一个 C 版本，它使用主轴来提供更具确定性的结果。
调用者需要对ortho_v3_v3的结果进行归一化处理。
inline int axis_dominant_v3_single(const float vec[3])
const float x = fabsf(vec[0]);
const float y = fabsf(vec[1]);
const float z = fabsf(vec[2]);
return ((x > y) ?
((x > z) ? 0 : 2) :
((y > z) ? 1 : 2));
/**
* Calculates \a p - a perpendicular vector to \a v
*
* \note return vector won't maintain same length.
*/
void ortho_v3_v3(float p[3], const float v[3])
const int axis = axis_dominant_v3_single(v);
assert(p != v);
switch (axis)
case 0:
p[0] = -v[1] - v[2];
p[1] = v[0];
p[2] = v[0];
break;
case 1:
p[0] = v[1];
p[1] = -v[0] - v[2];
p[2] = v[1];
break;
case 2:
p[0] = v[2];
p[1] = v[2];
p[2] = -v[0] - v[1];
break;
【讨论】：
以上是关于给定一个任意单位向量，计算任意正交单位向量的最佳方法是啥？的主要内容，如果未能解决你的问题，请参考以下文章