在单位半球的表面上快速均匀分布的随机点 [英] Fast uniformly distributed random points on the surface of a unit hemisphere

问题描述

我试图在蒙特卡罗射线追踪程序的单位球体表面上产生均匀的随机点。当我说均匀时，我的意思是点相对于表面积均匀分布。我当前的方法是在指向正z轴的半球上计算均匀随机点，并在xy平面中计算基底。

半球上的随机点表示方向

$ b
$ b

注意：dsfmt *将返回一个介于0和1之间的随机数。

 code> azimuthal = 2 * PI * dsfmt_genrand_close_open（& dsfmtt）; 
 zenith = asin（sqrt（dsfmt_genrand_close_open（& dsfmtt）））; 
 
 //计算笛卡尔点
 osRay.c._x = sin（zenith）* cos（azimuthal）; 
 osRay.c._y = sin（zenith）* sin（azimuthal）; 
 osRay.c._z = cos（zenith）; 
  

但是这是很慢的，并且分析表明它占用了大部分的运行时间。因此，我寻找了一些替代方法：

Marsaglia 1972拒绝方法

  do {
 x1 = 2.0 * dsfmt_genrand_open_open（& dsfmtt）-1.0; 
 x2 = 2.0 * dsfmt_genrand_open_open（& dsfmtt）-1.0 
 S = x1 * x1 + x2 * x2; 
} while（S> 1.0f）; 
 
 
 osRay.c._x = 2.0 * x1 * sqrt（1.0-S）; 
 osRay.c._y = 2.0 * x2 * sqrt（1.0-S）; 
 osRay.c._z = abs（1.0-2.0 * S）; 
  

分析笛卡尔坐标计算

  azimuthal = 2 * PI * dsfmt_genrand_close_open（& dsfmtt）; 
 u = 2 * dsfmt_genrand_close_open（& dsfmtt）-1; 
 w = sqrt（1-u * u）; 
 
 osRay.c._x = w * cos（azimuthal）; 
 osRay.c._y = w * sin（azimuthal）; 
 osRay.c._z = abs（u）; 
  

虽然最后两种方法的运行速度比第一种快，但当我使用它们时，表示它们不在球体表面上产生均匀的随机点，而是产生偏好赤道的分布。

此外，最后两种方法给出相同的最终结果但我确定它们是不正确的，因为我正在与一个分析解决方案比较。

我发现的每个引用表明这些方法产生均匀分布，但我不实现正确的结果。

我的实现中是否有错误，或者我错过了第二和第三种方法中的一个基本想法？

$实际上，例如在维度3中，e ^（ - x ^ 2/2）e ^（ - y ^ 2/2）e ^

这是很快的，如果你使用一个快速正态分布发生器（Ziggurat或Ratio-Of-Uniforms）和快速归一化程序（快速反平方根的谷歌）。不需要超越函数调用。

此外，Marsaglia在半球上不统一。由于在半球上的2D盘上的对应点不是等距的，因此在赤道附近将具有更多的点。最后一个似乎是正确的（但我没有做计算，以确保这一点）。

I am trying to generate uniform random points on the surface of a unit sphere for a Monte Carlo ray tracing program. When I say uniform I mean the points are uniformly distributed with respect to surface area. My current methodology is to calculate uniform random points on a hemisphere pointing in the positive z axis and base in the x-y plane.

The random point on the hemisphere represents the direction of emission of thermal radiation for a diffuse grey emitter.

I achieve the correct result when I use the following calculation :

Note : dsfmt* is will return a random number between 0 and 1.

azimuthal = 2*PI*dsfmt_genrand_close_open(&dsfmtt);
zenith = asin(sqrt(dsfmt_genrand_close_open(&dsfmtt)));

// Calculate the cartesian point
osRay.c._x = sin(zenith)*cos(azimuthal); 
osRay.c._y = sin(zenith)*sin(azimuthal);
osRay.c._z = cos(zenith);

However this is quite slow and profiling suggests that it takes up a large proportion of run time. Therefore I sought out some alternative methods:

The Marsaglia 1972 rejection method

do {
   x1 = 2.0*dsfmt_genrand_open_open(&dsfmtt)-1.0;
   x2 = 2.0*dsfmt_genrand_open_open(&dsfmtt)-1.0;
   S = x1*x1 + x2*x2;
} while(S > 1.0f);


osRay.c._x = 2.0*x1*sqrt(1.0-S);
osRay.c._y = 2.0*x2*sqrt(1.0-S);
osRay.c._z = abs(1.0-2.0*S);

Analytical cartesian coordinate calculation

azimuthal = 2*PI*dsfmt_genrand_close_open(&dsfmtt);
u = 2*dsfmt_genrand_close_open(&dsfmtt) -1;
w = sqrt(1-u*u);

osRay.c._x = w*cos(azimuthal);
osRay.c._y = w*sin(azimuthal);
osRay.c._z = abs(u);

While these last two methods run serval times faster than the first, when I use them I get results which indicate that they are not generating uniform random points on the surface of a sphere but rather are giving a distribution which favours the equator.

Additionally the last two methods give identical final results however I am certain that they are incorrect as I am comparing against an analytical solution.

Every reference I have found indicates that these methods do produce uniform distributions however I do not achieve the correct result.

Is there an error in my implementation or have I missed a fundamental idea in the second and third methods?

解决方案

The simplest way to generate a uniform distribution on the unit sphere (whatever its dimension is) is to draw independent normal distributions and normalize the resulting vector.

Indeed, for example in dimension 3, e^(-x^2/2) e^(-y^2/2) e^(-z^2/2) = e^(-(x^2 + y^2 + z^2)/2) so the joint distribution is invariant by rotations.

This is fast if you use a fast normal distribution generator (either Ziggurat or Ratio-Of-Uniforms) and a fast normalization routine (google for "fast inverse square root). No transcendental function call is required.

Also, the Marsaglia is not uniform on the half sphere. You'll have more points near the equator since the correspondence point on the 2D disc <-> point on the half sphere is not isometric. The last one seems correct though (however I didn't make the calculation to ensure this).

这篇关于在单位半球的表面上快速均匀分布的随机点的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！