The world’s leading publication for data science, AI, and ML professionals.

Create 3-D Galactic Art with Matplotlib

More Fun with Logarithmic Spirals

In a previous article, I demonstrated how you could use Python’s Tkinter GUI module to make 2-D galactic art using the equation for a logarithmic spiral [1]. In this article, we’ll take things a step further and use Python’s primary plotting library, matplotlib, to produce attractive 3-D simulations of spiral galaxies. Besides being fun, this project is a great way to present the concept of logarithmic spirals to students.

The Polar Equation for a Logarithmic Spiral

Logarithmic spirals are produced when natural objects grow geometrically. Some examples are galaxies, hurricanes, and pinecones. In the case of galaxies, these spirals are represented by the star-packed spiral arms.

Because spirals radiate from a central point or pole, they can be graphed with polar coordinates. In this system, the (x, y) coordinates used in the more familiar Cartesian coordinate system are replaced by (r, Ɵ), where r is the distance from the pole (at 0, 0), and Ɵ (theta) is the angle measured counterclockwise from the x-axis.

Using these terms, the polar equation for a logarithmic spiral is:

where e is the base of natural logarithms, a is a scaling factor that controls the size, and b is a growth factor that controls the spiral’s "openness" and direction of growth.

The Programming Strategy

As in the previous article, we’ll model a four-armed spiral galaxy by using the logarithmic equation to draw a single spiral and then rotate and redraw the spiral three more times. To model the central core of the galaxy we’ll randomly distribute points in a sphere. To capture the "background glow" of the galaxy, we’ll randomly distribute small "stars" around the galactic disc.

Finally, to make it all three-dimensional, we’ll randomly move stars in the galactic disc up and down a little, and then use matplotlib to display the results in a 3-D chart.

Installing Matplotlib

We’ll need matplotlib for this project. You can install it with either:

conda install matplotlib

or

pip install matplotlib

The Code

The following code was written in JupyterLab and is described by cell. You can download the full script from this Gist.

Importing Libraries and Setting Up the Display

For working with the logarithmic spiral equation and randomly shifting star locations we’ll use the math, random, and numpy libraries. For plotting we’ll use matplotlib.

We’ll display the results in a Qt window, which will quickly update and let us zoom, pan, save, and so on. And since we’re working in outer space, we’ll set the plot’s background style to dark.

Finally, we’ll represent the radius of our galactic disc with a constant called SCALE. Values between 200 and 700 will give the best results. You can relate these to real-world units by using the Milky Way galaxy as a reference. The radius of the Milky Way is 50,000 light-years. So, using a SCALE value of 500 means each pixel on your screen will equal 50,000/500 = 100 light-years.

import math
from random import randint, uniform, random
import numpy as np
import matplotlib.pyplot as plt

plt.style.use('dark_background')
%matplotlib qt

# Set the radius of the galactic disc (scaling factor):
SCALE = 350  # Use range of 200 - 700.

Defining a Function to Build a Spiral Arm

Next, we’ll define a function to build a single spiral arm using the logarithmic spiral equation. This equation will draw a single curving line, but we’ll "flesh it out" by varying the size and location of stars ("fuzzing"), and by duplicating the spiral for each arm and shifting it slightly backward while shrinking its stars. This will create a bright "leading" edge and a dim "trailing" edge.

The logarithmic spiral equation outputs points in two dimensions. To add 3-D "height" to our galaxy, we’ll create a z variable and assign it a random value scaled to the size of the galaxy (based on the SCALE constant).

The function will return a list of (x, y, z) locations for a single spiral arm.

def build_spiral_stars(b, r, rot_fac, fuz_fac):
    """Return list of (x,y,z) points for a logarithmic spiral.

    b = constant for spiral direction and "openness"
    r = scale factor (galactic disc radius)
    rot_fac = factor to rotate each spiral arm
    fuz_fac = randomly shift star position; applied to 'fuzz' variable
    """
    fuzz = int(0.030 * abs(r))  # Scalable initial amount to shift locations.
    num_stars = 1000
    spiral_stars = []
    for i in range(0, num_stars):
        theta = math.radians(i)
        x = r * math.exp(b*theta) * math.cos(theta - math.pi * rot_fac)
            - randint(-fuzz, fuzz) * fuz_fac
        y = r * math.exp(b*theta) * math.sin(theta - math.pi * rot_fac)
            - randint(-fuzz, fuzz) * fuz_fac
        z = uniform((-SCALE / (SCALE * 3)), (SCALE / (SCALE * 3)))
        spiral_stars.append((x, y, z))
    return spiral_stars

One thing to note here is that spiral galaxies are broad in the x and y directions but thin in the z-direction. Matplotlib’s default settings won’t handle this well, so, before plotting, we’ll manually limit the z-axis values to between -15 and 15 (-20 to 20 also yields good results).

I’m taking this exaggerated z-axis into account when scaling in the z-direction. Consequently, if you manually force all the axes to be the same scale, you may see some distortion.

Assigning Spiral Arm Parameters

To make a four-armed spiral galaxy, we’ll need to call the previous build_spiral_stars()function eight times (twice for each spiral to draw the leading and trailing edges). To prepare for this, let’s build a list of tuples that contain the last three arguments that we need to pass. Since the b parameter stays the same for each arm, we’ll just pass it directly when we call the function. The arguments I’m using were determined through trial and error.

# Assign scale factor, rotation factor, and fuzz factor for spiral arms.
# Each arm is a pair: leading arm + trailing arm:
arms_info = [(SCALE, 1, 1.5), (SCALE, 0.91, 1.5), 
             (-SCALE, 1, 1.5), (-SCALE, -1.09, 1.5),
             (-SCALE, 0.5, 1.5), (-SCALE, 0.4, 1.5), 
             (-SCALE, -0.5, 1.5), (-SCALE, -0.6, 1.5)]

Defining a Function to Build the Spiral Arms

The next function takes the arms_info list we just made, plus a value for b, and returns a list of all the stars in the leading edges of each of the four spirals, and a similar list for the trailing edges. We need separate lists as we’ll use a smaller marker size when plotting the trailing edges.

The function discriminates between leading and trailing edges by checking whether the arms_info index is odd or even.

def build_spiral_arms(b, arms_info):
    """Return lists of point coordinates for galactic spiral arms.

    b = constant for spiral direction and "openness"
    arms_info = list of scale, rotation, and fuzz factors
    """
    leading_arms = []
    trailing_arms = []
    for i, arm_info in enumerate(arms_info):
        arm = build_spiral_stars(b=b, 
                                 r=arm_info[0], 
                                 rot_fac=arm_info[1], 
                                 fuz_fac=arm_info[2])
        if i % 2 != 0:
            leading_arms.extend(arm)
        else:
            trailing_arms.extend(arm)            
    return leading_arms, trailing_arms

Defining a Function to Generate Spherical Coordinates

While most of a spiral galaxy consists of a flattened disc, the center tends to "bulge" out [3].

To create a bulge in our simulation, we’ll define a function that generates a spherical cloud of points. It will accept arguments for the number of points we want plus a radius, and will then generate random x, y, and z points and return them as a list of tuples.

We’ll generate the points using NumPy’s random.normal() function, which draws random samples from a normal distribution. It takes arguments for the mean, standard deviation, and output shape (dimensions). We’ll then multiply these points by the radius argument to scale them to the simulation.

To account for the vertical exaggeration in the matplotlib display, we’ll reduce the z values, found at index [2] in each coords tuple, by multiplying them by 0.02. (This value is somewhat arbitrary; feel free to play with it).

def spherical_coords(num_pts, radius):
    """Return list of uniformly distributed points in a sphere."""
    position_list = []
    for _ in range(num_pts):
        coords = np.random.normal(0, 1, 3)
        coords *= radius
        coords[2] *= 0.02  # Reduce z range for matplotlib default z-scale.
        position_list.append(list(coords))
    return position_list

Defining a Function to Build the Core Stars

Next, we’ll define a function that calls the previous function to build and return a list of coordinates for the outer and inner core. We’ll divide the core into two volumes to provide a visual gradation into the dense central area.

def build_core_stars(scale_factor):
    """Return lists of point coordinates for galactic core stars."""
    core_radius = scale_factor / 15
    num_rim_stars = 3000
    outer_stars = spherical_coords(num_rim_stars, core_radius)
    inner_stars = spherical_coords(int(num_rim_stars/4), core_radius/2.5)
    return (outer_stars + inner_stars)

Defining a Function for Scattering Star Haze

The space between the spiral arms isn’t devoid of stars, so we need a function to randomly cast points across the whole galactic model. Think of this as the "haze" or "glow" you see in photographs of galaxies.

We’ll include scalars for the radius (r_mult) and z-value (z_mult). These will let us "thicken" the galactic disc towards the center by calling the function twice to generate "inner" and "outer" sets of points.

def haze(scale_factor, r_mult, z_mult, density):
    """Generate uniform random (x,y,z) points within a disc for 2-D display.

    scale_factor = galactic disc radius
    r_mult = scalar for radius of disc
    z_mult = scalar for z values
    density = multiplier to vary the number of stars posted
    """
    haze_coords = []
    for _ in range(0, scale_factor * density):
        n = random()
        theta = uniform(0, 2 * math.pi)
        x = round(math.sqrt(n) * math.cos(theta) * scale_factor) / r_mult
        y = round(math.sqrt(n) * math.sin(theta) * scale_factor) / r_mult
        z = np.random.uniform(-1, 1) * z_mult
        haze_coords.append((x, y, z))
    return haze_coords

Here’s the impact of star haze versus no star haze:

Generating Star Coordinates and Plotting the Galaxy

The following code completes the program by calling the functions that build the star coordinate lists and then passing them to matplotlib’s scatter() method.

To prepare the data for plotting, we use the splat operator (*) such as in *zip(*leading_arm). This syntax unzips a list of tuples, effectively transposing rows into columns. Tuples are unpacked, and each element is passed as a separate argument to the zip() function. This creates separate tuples for the x, y, and z values, with the correct corresponding values in each tuple. Here’s an example:

a_list = [(x1, y1, z1), (x2, y2, z2)]

zip(a_list) => (x1, x2), (y1, y2), (z1, z2)

# Create lists of star positions for galaxy:
leading_arm, trailing_arm = build_spiral_arms(b=-0.3, arms_info=arms_info)
core_stars = build_core_stars(SCALE)
inner_haze_stars = haze(SCALE, r_mult=2, z_mult=0.5, density=5)
outer_haze_stars = haze(SCALE, r_mult=1, z_mult=0.3, density=5) 

# Plot stars in 3D using matplotlib:
fig, ax = plt.subplots(1, 1, 
                       subplot_kw={'projection': '3d'}, 
                       figsize=(12, 12))
ax.set_axis_off()
ax.set_zlim (-15, 15)

ax.scatter(*zip(*leading_arm), c='w', marker='.', s=5)
ax.scatter(*zip(*trailing_arm), c='w', marker='.', s=2)
ax.scatter(*zip(*core_stars), c='w', marker='.', s=1)
ax.scatter(*zip(*inner_haze_stars), c='w', marker='.', s=1)
ax.scatter(*zip(*outer_haze_stars), c='lightgrey', marker='.', s=1)

That’s a nice-looking galaxy. I want to live there!

As you may have noticed, our code includes several parameters that you can use to tweak the output and indulge your inner artist.

For instance, changing the b parameter in the build_spiral_arms() function changes the appearance of the leading and trailing arms. Here’s an example where b is set to -0.5:

Changing the SCALE constant affects the size of the galactic bulge. Here are some examples using different values:

In the previous figure, notice how using light grey for the outer haze stars produces "dust belts" similar to those you see in profiles of real galaxies.

You can also change the marker styles. For example, if you want to use actual star shapes, add this snippet:

import matplotlib.path as mpath
star = mpath.Path.unit_regular_star(6)

Then pass star as the marker argument in the call to ax.scatter(). The (6) indicates the number of vertices. A value of 5 produces a five-pointed star and 6 produces six points.

Bonus: Generating Globular Clusters

We’ve finished the spiral galaxy project, but I would be remiss if I didn’t show you how to build an astronomical oddity called a _globular cluster_ [5]. __ These are spherical collections of stars that orbit most spiral galaxies such as our Milky Way. They are among the oldest features in a galaxy and can contain millions of tightly packed stars [6].

The function we defined to build the spherical core stars is perfect for building globular clusters, so why not use it? Here’s the code:

def spherical_coords(num_pts, radius):
    """Return list of uniformly distributed points in a sphere."""
    position_list = []
    for _ in range(num_pts):
        coords = np.random.normal(0, 1, 3) 
        coords *= radius 
        position_list.append(list(coords))
    return position_list

rim_radius = 2
num_rim_stars = 3000
rim_stars = spherical_coords(num_rim_stars, rim_radius)
core_stars = spherical_coords(int(num_rim_stars/4), rim_radius/2.5)

fig, ax = plt.subplots(1, 1, subplot_kw={'projection':'3d'})
ax.axis('off')
ax.scatter(*zip(*core_stars), s=0.5, c='white')
ax.scatter(*zip(*rim_stars), s=0.1, c='white')
ax.set_xlim(-(rim_radius * 4), (rim_radius * 4))
ax.set_ylim(-(rim_radius * 4), (rim_radius * 4))
ax.set_zlim(-(rim_radius * 3), (rim_radius * 3))
ax.set_aspect('auto')

And here’s the result:

Summary

With Python, matplotlib, and a simple equation, we were able to easily generate 3-D simulations of a spiral galaxy. Besides producing some neat digital art, this project is a great way to introduce logarithmic spirals to students in diverse fields of study including mathematics, Astronomy, and programming.

Citations

  1. Vaughan, L. (2023), "Create Galactic Art with Tkinter: Model Mother Nature with Logarithmic Spirals," Towards Data Science, https://towardsdatascience.com/create-galactic-art-with-tkinter-e0418a59b215.
  2. Vaughan, Lee, 2018, Impractical Python Projects: Playful Programming Activities to Make You Smarter, No Starch Press, San Francisco.
  3. Wikipedia contributors. (2023, October 3). Galactic bulge. In Wikipedia, The Free Encyclopedia. Retrieved 16:23, October 3, 2023, from https://en.wikipedia.org/w/index.php?title=Galactic_bulge&oldid=1178373003.
  4. ESO/NASA/JPL-Caltech/M. Kornmesser/R. Hurt, CC BY 4.0 <https://creativecommons.org/licenses/by/4.0>, via Wikimedia Commons, File:Artist’s impression of the central bulge of the Milky Way.jpg – Wikimedia Commons.
  5. Wikipedia contributors. (2023, September 24). Globular Cluster. In Wikipedia, The Free Encyclopedia. Retrieved 00:56, October 4, 2023, from https://en.wikipedia.org/w/index.php?title=Globular_cluster&oldid=1176809907.
  6. Vaughan, Lee, 2023, Python Tools for Scientists: An Introduction to Using Anaconda, JupyterLab, and Python’s Scientific Libraries, No Starch Press, San Francisco.

Thanks!

Thanks for reading and please follow me for more Quick Success Data Science projects in the future.


Related Articles