1

I would like to sample points around the surface of a unit sphare without repetition in python. phi - [0-2pi] and theta [0, pi]

I have this code:

def get_random_angles(n):
    phis = np.random.uniform(0, 2*np.pi, n)
    thetas = np.arccos(1 - np.random.uniform(0, 2, n))
    
    return thetas, phis

How can I make sure the generated numbers do not reapeat?

  • 2
    The probability that two points are equal is 0 (or extremely close to 0 since the computer random numbers have a finite number of digits). Hence no need to check for that. – RandomGuy Apr 26 '23 at 14:58
  • If you absolutely wanted to verify that no phi-theta pair is repeated, you could put the zipped results into a set, and verify the length of the set matches n. If not, then generate n - len(set) additional angles, check that for repeats, and then merge the results. However, as RandomGuy said, the probability of repeats is so low its probably not worth checking for. – nigh_anxiety Apr 26 '23 at 15:29

2 Answers2

1

Something like this did the trick:

def get_random_angles(n):

    result=set()

    while len(result) < n:
        phi = np.random.uniform(0, 2 * np.pi, size=n-len(result))
        theta = np.arccos(1 - np.random.uniform(0, 2, size=n-len(result)))
        new_tuples = set(zip(theta, phi))
        result.update(new_tuples)

    return result
  • This will only attempt to reach `n` values one extra time. In my edited answer, there's a recursive while loop that will guarantee the function outputs `n` coordinates – bc1155 Apr 26 '23 at 15:42
  • sorry, but I did not get it. What do you mean by one extra time? –  Apr 26 '23 at 15:49
  • On second look, I was wrong. Your function is recursive too, so it will continue to run until you get the right number—my bad! The if statement, `if (length is not n)` should be tweaked to `if length != n`, though. See [here](https://stackoverflow.com/questions/132988/is-there-a-difference-between-and-is) for an explanation as to why. – bc1155 Apr 26 '23 at 15:58
  • 1
    It's recursive, but the recursive calls don't update the caller's set object - they just return their result to the caller to update its set. So the code does indeed have the problem you identified @bc1155 – slothrop Apr 26 '23 at 16:12
  • 1
    I redid the implementation. Now it should be OK –  Apr 26 '23 at 16:31
0

There are 3 options I can think of:

  1. you can pre-generate non-redundant coordinates to a suitable precision and use the sample function from the random built-in library
  2. you can save the coordinates to a list and do if coord in coord_list conditional before adding to the list
  3. you can use a set(), which cannot contain any duplicate objects, i.e. set(1,1,2,3,4,4,5) becomes {1,2,3,4,5}.

Given your function though, I think the number of repetitions will be negligible unless you're generating a lot of coordinates.

Applying the most performant option (3) to your code would look something like:

import numpy as np

def get_random_angles(n):
    phis = np.random.uniform(0, 2*np.pi, n)
    thetas = np.arccos(1 - np.random.uniform(0, 2, n))
    while l := len(set(zip(thetas,phis))) != n:
        t,p=get_random_angles(n-l)
        np.concatenate(thetas,t)
        np.concatenate(phis,p)
    return thetas, phis
bc1155
  • 263
  • 7
  • 1
    (3) has a performance advantage over (2), since checking for the presence of an item in an N-element list is O(N), whereas in the corresponding set it's O(1). – slothrop Apr 26 '23 at 15:04
  • The `set` approach seems interesting. Is there a one-liner to make work with the current function? Such that the size `n` is mantained? –  Apr 26 '23 at 15:10
  • 1
    Assuming that you're looking for unique pairs of `thetas` and `phis`, it would be something like `set(zip(thetas, phis))`, otherwise `set(thetas)` and `set(phis)`. My strategy would be to add that to the end of the function and check its length. If the length < n, then generate n minus length new coordinates, update the set, and check again, repeat if necessary, and return when length == n. For set documentation see [here](https://docs.python.org/3/tutorial/datastructures.html#sets) – bc1155 Apr 26 '23 at 15:19
  • This worked. I will post the final version –  Apr 26 '23 at 15:37