Create random number within an annulus

Question

I am trying generate a random number that is within an annulus, i.e. we have a max and min radius. I tried doing:

while True:
    x=random.uniform(-maxR, maxR)
    y=random.uniform(-maxR, maxR)
    R=math.sqrt(x**2 + y**2)
    if R <= maxRadius and R >= minRadius:
        if x>= -maxRadius and x <= maxRadius and x<=-minRadius and x>= minRadius:
            print "passed x"
            if y>= -maxRadius and y <= maxRadius and y<=-minRadius and y>= minRadius: 
                break

But this is very slow. Is it possible to feed more contraints into random.uniform or is there another method?

What are the last two if-statements for? Why is simply checking if `R` (the radius) is within range not sufficient? — Rob Wouters, Jan 28 '12 at 20:04
You _could_ always generate the numbers beforehand and store them in a dict mapping. — Joel Cornett, Jan 28 '12 at 20:25
There are already endless questions about drawing from non-unifrom PDFs, and from circles in particular. The better answers to those question contain the method for drawing from arbitrary invertable PDFS. — dmckee --- ex-moderator kitten, Jan 28 '12 at 20:26

score 22 · Answer 1 · edited May 23 '17 at 11:59

In general you can either draw the correct distribution directly or use rejection.

To draw directly use

draw theta uniformly on [0,2pi): theta = random.uniform(0,2*pi)
draw r from the power-law distribution r^1.

The only complexity compared to doing this for a circle is that you PDF runs from [r_min,r_max] not [0,r_max]. This leads to

CDF = A \int_{r_min}^{r} r' dr' = A (r^2 - r_min^2)/2

for A the normalizing constant
```
A = 2/(r_max*r_max - r_min*r_min)
```
implying that
```
r = sqrt(2*random.uniform(0,1)/A + r_min*r_min)
```
and you can simplify slightly.
then compute (x,y) by the usual transformation from radial coordinates
x = r * cos(theta)
y = r * sin(theta)

This method of integrating the PDF, normalizing the CDF and inverting is general and is sometimes called the "Fundamental Theorem of Sampling".

Rejection

Draw (x,y) on a box big enough to contain the annulus, then reject all cases where `r = sqrt(xx + yy) exceeds r_max or is less than r_min.

This is reasonably efficient if the hole in the middle is small, and very inefficient if the hole is large.

I have drawn this solution in Mathematica: http://i.imgur.com/7L04w.png It's looking nice and uniform! :D — Ricardo Sanchez-Saez, Oct 25 '12 at 10:36
"and you can simplify slightly." Am I correct in thinking that in the general case (assuming your language of choice supports it) it can be simplified to this: `r = sqrt(random.uniform(r_min*r_min,r_max*r_max))` ? — dgnuff, Jan 21 '16 at 23:52

Create random number within an annulus

1 Answers1

To draw directly use

Rejection

Linked

Related