旋转球面系统的点
Rotating a Spherical System of Points
假设我在一个球体上有10个点(随机分布),我想旋转整个系统以确保有一个点位于北极。我该如何使用c++实现这一点呢?
我通过观察3D旋转矩阵来解决这个问题:
http://en.wikipedia.org/wiki/Rotation_matrix我绕x轴旋转我的点直到y分量为零,然后绕y轴旋转直到x分量为零。这应该使问题的点在北极或南极,对吧?
我的代码是这样的:
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <iostream>
#include <iomanip>
#include <fstream>
#include <time.h>
#include <stdlib.h>
#include <sstream>
using namespace std;
#define PI 3.14159265358979323846
int main()
{
int a,b,c,f,i,j,k,m,n,s;
double p,Time,Averagetime,Energy,energy,Distance,Length,DotProdForce,Forcemagnitude,
ForceMagnitude[101],Force[101][4],E[1000001],En[501],x[101][4],y[101][4],
Dist[101][101],Epsilon,z[101],theta,phi;
/* set the number of points */
n=10;
/* check that there are no more than 100 points */
if(n>100){
cout << n << " is too many points for me :-( n";
exit(0);
}
/* reset the random number generator */
srand((unsigned)time(0));
for (i=1;i<=n;i++){
x[i][1]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
x[i][2]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
x[i][3]=((rand()*1.0)/(1.0*RAND_MAX)-0.5)*2.0;
Length=sqrt(pow(x[i][1],2)+pow(x[i][2],2)+pow(x[i][3],2));
for (k=1;k<=3;k++){
x[i][k]=x[i][k]/Length;
}
}
/* calculate the energy */
Energy=0.0;
for(i=1;i<=n;i++){
for(j=i+1;j<=n;j++){
Distance=sqrt(pow(x[i][1]-x[j][1],2)+pow(x[i][2]-x[j][2],2)
+pow(x[i][3]-x[j][3],2));
Energy=Energy+1.0/Distance;
}
}
cout << fixed << setprecision(10) << "energy=" << Energy << "n";
/* Save Values so far */
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
y[i][j]=x[i][j];
}
}
/* Choose each point in turn and make it the north pole note this is what the while loop os for, but have set it to a<2 to just look at first point */
a=1;
b=0;
c=0;
while(a<2){
/* Find theta and phi to rotate points by */
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "n";
theta=x[a][2]/x[a][3];
theta=b*PI+atan(theta);
/* Rotate Points by theta around x axis and then by phi around y axis */
for(i=1;i<=n;i++){
x[i][1]=x[i][1];
x[i][2]=x[i][2]*cos(theta)-x[i][3]*sin(theta);
x[i][3]=x[i][2]*sin(theta)+x[i][3]*cos(theta);
}
phi=x[a][1]/x[a][3];
phi=c*PI+atan(phi);
for(i=1;i<=n;i++){
x[i][1]=x[i][1]*cos(phi)-x[i][3]*sin(phi);
x[i][2]=x[i][2];
x[i][3]=x[i][1]*sin(phi)+x[i][3]*cos(phi);
}
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "n";
if(x[a][3]==1.0 && x[a][2]==x[a][3]==0)
a=a+1;
else if(b==0 && c==0)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
c=1;
}
}
else if(b==0 && c==1)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
b=1;
c=0;
}
}
else if(b==1 && c==0)
for(i=1;i<=n;i++){
for(j=1;j<=3;j++){
x[i][j]=y[i][j];
c=1;
}
}
else if(b==1 && c==1)
break;
}
energy=0.0;
for(i=1;i<=n;i++){
for(j=i+1;j<=n;j++){
Distance=sqrt(pow(x[i][1]-x[j][1],2)+pow(x[i][2]-x[j][2],2)
+pow(x[i][3]-x[j][3],2));
energy=energy+1.0/Distance;
}
}
cout << fixed << setprecision(10) << "ENERGY=" << energy << "n";
cout << fixed << setprecision(5) << "x[" << a << "][1]=" << x[a][1] <<
" x[" << a << "][2]=" << x[a][2] << " x[" << a << "][3]=" << x[a][3] << "n";
/* Output to help with gnuin.txt */
ofstream File4 ("mypoints");
for(i=1;i<=n;i++){
File4 << x[i][1] << " " << x[i][2] << " " << x[i][3] << "n";
}
File4.close();
return 0;
}
好的,这里有很多问题,比如if语句(第103行)不应该有一个等于双精度的条件,因为它永远不会工作,但我可以稍后使用间接比较(一些epsilon的东西)来解决这个问题。我真正的问题是,为什么旋转即使作用于所有的点也会使球面上的点离开?(正如你所看到的,这些值已经被标准化,使它们都在第38行单位球体上)。
注意:b, c的东西是检查点是在北极还是南极
您的旋转代码有问题。例如:
x[i][1]=x[i][1];
x[i][2]=x[i][2]*cos(theta)-x[i][3]*sin(theta);
x[i][3]=x[i][2]*sin(theta)+x[i][3]*cos(theta);
在第二行修改x[i][2]
,然后在第三行使用它。对于中间结果,应该使用临时存储,以避免在引用完成之前修改值。
第一行是相当多余的,其余的应该看起来更像:
double new_y, new_z;
new_y=x[i][2]*cos(theta)-x[i][3]*sin(theta);
new_z=x[i][2]*sin(theta)+x[i][3]*cos(theta);
x[i][2] = new_y;
x[i][3] = new_z;
(在执行这样的计算时显然要这样做)
一个更好的方法来定位你的球体,可能是计算一个变换矩阵,以相同的方式作为一个"看"矩阵。在查看矩阵中,帧被旋转以使某个向量与z轴对齐。在您的情况下,您可能希望沿y轴对齐,但原理是一样的。
我还评论说,你似乎忽略了数组中的第0个元素…恕我直言,这是个坏习惯——您应该习惯数组从0开始的事实。迟早你会弄错索引,或者你将不得不接口到别人的代码。
相关文章:
- C++,系统无法执行指定的程序
- 与互斥锁相比,旋转锁可以保证上下文切换
- 在UNIX系统中使用DIR查找文件的字节大小
- 错误处理.将系统错误代码映射到泛型
- 当系统的卷被修改时,如何修改WASAPI环回捕获卷
- 绘制旋转的三角形
- 有什么好的方法可以让系统调用代理允许在单元测试中进行模拟
- 在C++游戏中与库存系统作斗争
- 旋转模型矩阵时的形状失真
- 四边形的 2D 旋转
- 文件系统:复制功能的速度秘诀是什么
- c++17文件系统::recursive_directory迭代器()在mac上没有给出这样的目录,但在windows上
- 在gtest.中使用fff.h模拟系统API
- 垂直方向的 Gtk3+ 旋转按钮 (c/c++)
- 如何制作无限制照明系统
- 系统.将数组移交给c#中动态加载的c++DLL时发生AccessViolationException
- 发布旋转矩阵(openGL/glm)
- 如何旋转粒子系统,使其不绕世界轴旋转
- Cocos2d粒子系统不跟随移动和旋转的摄像机
- 旋转球面系统的点