RE: Random spherical troubles
- To: mathgroup at smc.vnet.net
- Subject: [mg25198] RE: [mg25170] Random spherical troubles
- From: PNichols at cornell-iowa.edu
- Date: Tue, 12 Sep 2000 21:24:43 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
Barbara-- I think this will do what you need; some explanation follows later. ============================================ phi := Random[Real, {0, 2 Pi}; (* uniform distribution, for longitude *) theta := ArcCos[1 - 2 Random[]]; (* non-uniform distribution, for latitude *) Note the use of :=, not =, for then we can make a list of random points thus: sphericalData = Table[{phi, theta}, {5000}]; ListPlot[sphericalData] The ListPlot will clearly show the sparseness of the points near theta=0 or theta=Pi. To map the points into the sphere, I'll use sphereMap[{phi_, theta_}] := {Sin[theta] Cos[phi], Sin[theta] Sin[phi], Cos[theta]} mappedSphericalData=sphereMap /@ sphericalData; To display this as a "cloud" of 3D points, use Show[Graphics3D[Point /@ mappedSphericalData]]; Here's the contrasting "wrong" case -- uniform distribution on the rectangle -- in similar style: rectangularData = Table[{Random[Real, {0, Pi}], Random[Real, {0, Pi}]}, {5000}]; ListPlot[rectangularData]; (Note uniform distribution on the rectangle.) mappedRectangularData = sphereMap /@ rectangularData; Show[Graphics3D[Point /@ mappedRectangularData]]; If you have Mathematica version 4, evaluate <<RealTime3D` before the Show statements above, and then you can move the pictures around with your mouse. That makes the 3D pictures easier to understand. (Use <<Default3D` to revert to ordinary 3D display.) ============================================ DISCUSSION (more mathematics than Mathematica) The formula theta := ArcCos[1 - 2 Random[]]; needs explanation. I'm not an expert in probability theory, but here goes. Because the Jacobian of the sphere map measures the distortion of the area, we need random values of theta with a probability density f[theta_]:=Sin[theta]/2 on [0, Pi]. The Jacobian is independent of phi, reflecting the fact that the "pull-back" area-measure from the sphere to the rectangle is uniform along the phi-direction. (There is also a factor of rho^2, where rho is the radius of the sphere; we use rho=1. The division by 2 makes the total probability 1.) There is a theorem (not terribly difficult to prove) which states that if you have a random variable X uniformly distributed on [0,1] (like Random[]), and you want a random variable Y with density f, then you can do this (warning: pseudo-Mathematica pseudo-code follows!): (1) Find the "cumulative distribution function", cdf[y_]=Integrate[f[t],{t,0,y}]; (2) solve x==cdf[y] for y (i.e. find the inverse function of the cdf); (3) Substitute the random variable X for the dummy variable x in the inverse function just found. (All this works well if f is non-negative and zero only at isolated points, so that the cdf is continuous and strictly increasing. More general cases require a generalized type of "inverse function", but that's not relevant for your problem.) Since the functions in this problem are so nice, Mathematica can do nearly everything for us: Solve[x == Integrate[Sin[y]/2, {y, 0, y0}], y0] {{y0 -> -ArcCos[1 - 2*x]}, {y0 -> ArcCos[1 - 2*x]}} There, in the second solution, is the function your problem needs. I have made a notebook containing more detailed notes on this problem; let me know if you want a copy. Best regards, Preston Nichols Department of Mathematics Cornell College ============================================ -----Original Message----- From: Barbara DaVinci [mailto:barbara_79_f at yahoo.it] To: mathgroup at smc.vnet.net Subject: [mg25198] [mg25170] Random spherical troubles Hi MathGrouppisti This time, my problem is to generate a set of directions randomly distributed over the whole solid angle. This simple approach is incorrect (spherical coordinates are assumed) : Table[{Pi Random[], 2 Pi Random[]} , {100}] because this way we obtain a set of point uniformly distributed over the [0 Pi] x [0 2Pi] rectangle NOT over a spherical surface :-( If you try doing so and plot the points {1, random_theta , random_phi} you will see them gathering around the poles because that simple transformation from rectangle to sphere isn't "area-preserving" . Such a set is involved in a simulation in statistical mechanics ... and I can't get out this trouble. May be mapping [0 Pi] x [0 2Pi] in itself , using an suitable "non-identity" transformation, can spread points in a way balancing the poles clustering effect. ==================================================================== While I was brooding over that, an intuition flashed trought my mind : since spherical to cartesian transformation is x = rho Sin[ theta ] Cos[ phi ] y = rho Sin[ theta ] Sin[ phi ] z = rho Cos[ theta ] perhaps the right quantities to randomly spread around are Cos[ theta ] and Cos[ phi ] rather than theta and phi for itself. Give a glance at this : Table[{ ArcCos[ Random[] ], ArcCos[ Random[] Sign[ 0.5 - Random[] ] } , {100}] Do you think it is close to the right ? Do you see a better way ? Have you just done the job in the past ? Should I reinvent the wheel ? ==================================================================== I thanks you all for prior replies and in advance this time. Distinti Saluti (read : "Faithfully yours") Barbara Da Vinci barbara_79_f at yahoo.it ______________________________________________________________________ Do You Yahoo!? Il tuo indirizzo gratis e per sempre @yahoo.it su http://mail.yahoo.it