implemented uniformly random rotation quaternion

This commit is contained in:
Jiri Pittner 2022-11-02 14:23:26 +01:00
parent 4a8e008811
commit f368a6dd50
2 changed files with 28 additions and 0 deletions

View File

@ -408,6 +408,33 @@ r *= ::pow(xnorm,y);
return r; return r;
} }
#ifndef QUAT_NO_RANDOM
//method to efficiently generate a uniformly random point on a unit 4-sphere by G. Marsaglia, Ann. Math. Stat. 43, 645 (1972)
template<typename T>
void Quaternion<T>::random_rotation()
{
T s1,s2,s;
do
{
q[0]=2.*random()/(1. + RAND_MAX) - 1.;
q[1]=2.*random()/(1. + RAND_MAX) - 1.;
s1 = q[0]*q[0] + q[1]*q[1];
}
while(s1>1);
do
{
q[2]=2.*random()/(1. + RAND_MAX) - 1.;
q[3]=2.*random()/(1. + RAND_MAX) - 1.;
s2 = q[2]*q[2] + q[3]*q[3];
}
while(s2>1);
s = sqrt((1-s1)/s2);
q[2] *= s;
q[3] *= s;
}
#endif
//force instantization //force instantization
#define INSTANTIZE(T) \ #define INSTANTIZE(T) \
template class Quaternion<T>; \ template class Quaternion<T>; \

View File

@ -124,6 +124,7 @@ public:
void axis2normquat(const T *axis, const T &angle); void axis2normquat(const T *axis, const T &angle);
void normquat2axis(T *axis, T &angle) const; void normquat2axis(T *axis, T &angle) const;
void random_rotation(); //generate uniformly random unit quaternion
//C-style IO //C-style IO
int fprintf(FILE *f, const char *format) const {return ::fprintf(f,format,q[0],q[1],q[2],q[3]);}; int fprintf(FILE *f, const char *format) const {return ::fprintf(f,format,q[0],q[1],q[2],q[3]);};