1

This is a follow up question to this one: pyephem problem with alt/az to ra/dec and back

I want to distribute stars randomly above an observer. I was able to to do this in the last question.

The issue I now get is that the stars seems clustered towards the center.

How can I get rid of this?

import matplotlib.pyplot as plt
import numpy as np
import  math
import ephem
from datetime import datetime
import random

home = ephem.Observer()

home.lon = '-70.4'  # +E
home.lat = '-24.62'  # +N
home.elevation = 0  # meters
home.date = datetime.utcnow()

theta_plot = []
r_plot=[]

stars_pre_selected = []

for n in range(4000):
    ra, dec = home.radec_of(az=random.uniform(0, 2 * math.pi),
                            alt=random.uniform(0, math.pi))

    body = ephem.FixedBody()
    body._epoch = ephem.J2000
    body._ra = ephem.degrees(ra)
    body._dec = ephem.degrees(dec)

    body.compute(home)


    if math.degrees(body.alt) > 0:
        stars_pre_selected.append(body)


    else:
        print(math.degrees(body.alt))


for body in stars_pre_selected:
    body.compute(home)
    theta_plot.append(body.az)
    r_plot.append(math.cos(body.alt))
    print(body.alt,body.az)

ax = plt.subplot(111, polar=True)
#
ax.set_theta_direction(-1)
ax.set_theta_offset(np.pi / 2)
ax.grid(True)
ax.scatter(theta_plot, r_plot, c='blue', s=2)
ax.set_rmax(1.0)

plt.pause(0.1)

plt.show()

enter image description here

UPDATE:

    ra, dec = home.radec_of(az=random.uniform(0, 2 * math.pi),
    alt=math.asin(2*random.random()-1))
    

this will work too, but then the code hast to be changed to math.sin instead of math.cos :

    ra, dec = home.radec_of(az=random.uniform(0, 2 * math.pi),
    alt=math.acos(2*random.random()-1))

enter image description here

Taiko
  • 51
  • 7
  • 1
    Does this answer your question? [Fast uniformly distributed random points on the surface of a unit hemisphere](https://stackoverflow.com/questions/7280184/fast-uniformly-distributed-random-points-on-the-surface-of-a-unit-hemisphere) – Peter O. Jun 14 '21 at 13:17
  • sorry , i cannot read c++ – Taiko Jun 14 '21 at 13:21
  • i also found this [link](https://meyavuz.wordpress.com/2018/11/15/generate-uniform-random-points-within-a-circle/) which seems to be what i want but i cannot get it to work : / – Taiko Jun 14 '21 at 13:24
  • i found another [link](https://mathworld.wolfram.com/SpherePointPicking.html) that explains how it should work – Taiko Jun 14 '21 at 21:41

4 Answers4

3

You responded to Peter O's suggestion by saying that you don't read C, but while the question is in C, the answer https://stackoverflow.com/a/7280536/8544123 is language agnostic. In Python, it would be

    import numpy as np
    #create random three dimension vector
    vector = np.array([np.random.normal() for _ in range(3)])
    #find length
    length = np.dot(vector,vector)
    #normalize vector
    unit_vector = np.divide(vector, length)

This will create a vector on the unit sphere. You can then convert it to azimuth and altitude. Since you want a hemisphere, you can then take the absolute value of the altitude.

Acccumulation
  • 3,491
  • 1
  • 8
  • 12
1

I guess you can generate random position in a square, but then check, if this position is inside circle - place new point. This will let you avoid clustering in the center.

ukra174
  • 13
  • 1
  • 3
1

If I understand your application correctly, this is because a 360 degree 'ring' near the top of your hemisphere covers a smaller area than a 'ring' at the bottom of your sphere. Because of this, if you have the same number of stars in each ring, they will appear denser in the smaller rings (the center in your image).

The alt argument to home.radec_of() distributes uniformly around 0 to pi, although note that since you sweep 'az' by a full circle, you could actually just have an alt distribution of 0 to pi/2, I'll do that in the rest of my answer.

You want to have fewer stars in smaller rings, so you can't use a uniform distribution from 0 to pi/2. Instead you want more close to 0 (the larger ring) and fewer by pi/2.

My best pointer is to look up how to do random distributions in numpy where you specify the distribution characteristics.

As for the distribution itself, find the size of a ring at each value of alt and have the number of values in your distribution proportional to that size.

MichaelCG8
  • 579
  • 2
  • 14
0
ra, dec = home.radec_of(az=random.uniform(0, 2 * math.pi),
alt=math.asin(2*random.random()-1))
Taiko
  • 51
  • 7