Exploring the Visible Spectrum in Python

The electromagnetic spectrum covers a vast range of wavelengths and frequencies, only a tiny fraction of which is visible to the human eye. Wavelengths and frequencies are inversely proportional and the relationship between the two can easily be plotted. In this post I will write code in Python to do just that, using the actual colours each wavelength/frequency combination represents.

Wavelengths of electromagnetic radiation range from 1 picometre (one trillionth of a metre) to 100,000 km (about a quarter of the way to the moon). Each wavelength has a corresponding frequency, and it is best to think of wavelengths and frequencies as two ways to measure the same thing.

The visible part of the electromagnetic spectrum has wavelengths from approximately* 380nm (nanometres, or billionths of a metre) to 780nm. The corresponding frequencies go from 789THz (terahertz, 1 trillion Hertz) to 384Thz. As I mentioned above wavelengths and frequencies are inversely proportional so as wavelengths increase frequencies decrease.

* There is no fixed definition of visible light; 380nm seems widely accepted as the starting point of the violet end of the spectrum but estimates of the finishing point of the red end range from just 700nm to the 780nm I will be using.

The following graphic should make the wavelength/frequency relationship clear.

The violet wave has a short wavelength, ie the distance between two points at the same height. This means we can cram 4 waves into this small graph: low wavelength = high frequency.

The red wave has a wavelength twice as long as the violet wave so we can only get two of them into the same space: high wavelength = low frequency.

Calculating frequency from wavelength is simple, we just divide the speed of light by the wavelength. We need to be careful with units though: with loads of nanos, teras etc. floating around we can easily get errors of many orders of magnitude. I'll deal with this when we get to the code.

A particular wavelength/frequency also has an associated energy. I won't plot those here but will calculate them for display in a table. It is a tiny constant called Planck's constant multiplied by the frequency, the result being in Joules.

The project consists of the following two files which you can download as a zip from the Downloads page, or you can clone/download the Github repository if you prefer. You will also need the Pillow imaging library, documentation for which is here https://pillow.readthedocs.io/en/stable/ .

• visiblespectrum.py
• visiblespectrum_demo.py

The first file, visiblespectrum.py, contains a set of functions to generate a data structure containing data on the integer wavelengths of visible light, and to generate plots of this data.

visiblespectrum.py

from PIL import Image, ImageDraw, ImageFontdef generate_data():data = []c = 3 * 10**8   # speed of light in m/s h = 6.62607015 * 10**-34   # Planck's constant in Js for nm in range(380, 781):item = {}item["nm"] = nmitem["Hz"] = (c / (nm * 10**-9))item["THz"] = int(item["Hz"] * 10**-12)item["rgb"] = _wavelength_to_rgb(nm)item["J"] = h * item["Hz"]item["E"] = item["J"] * 10**19data.append(item)return datadef print_data(data):width = 37pow_minus19 = chr(8315) + chr(185) + chr(8313)print("-" * width)print(f"|{chr(955)}(nm)|f(THz)|E(J)      |R  |G  |B  |")print("-" * width)for item in data:print(f"|{item['nm']:>5.0f}|", end="")print(f"{item['THz']:>6.0f}|", end="")print(f"{item['E']:>4.2f}x10{pow_minus19}|", end="")print(f"{item['rgb'][0]:>3.0f}|", end="")print(f"{item['rgb'][1]:>3.0f}|", end="")print(f"{item['rgb'][2]:>3.0f}|")print("-" * width)def plot_wavelength_frequency(data, filename):border_width = 70width_scaling = 2height_scaling = 0.5image_width = int(401 * width_scaling) + (border_width * 2)image_height = int(800 * height_scaling) + (border_width * 2)column_width = width_scalingcolumn_bottom = image_height - border_widthx = border_widthimage = Image.new("RGB", (image_width, image_height), (32, 32, 32))image = _draw_labels(image, "Visible Spectrum", "Wavelength (nm)", "Frequency (THz)")_draw_x_axes(image, border_width, 380, 780, 50)_draw_y_axes(image, border_width, 0, 800, 50)draw = ImageDraw.Draw(image)for item in data:column_top = column_bottom - (item["THz"] * height_scaling)draw.rectangle(xy=[x, column_bottom, x - column_width, column_top],fill=item["rgb"])x += column_widthtry:image.save(filename, "PNG")except IOError as e:print(e)def plot_frequency_wavelength(data, filename):border_width = 70width_scaling = 2height_scaling = 0.5image_width = int(401 * width_scaling) + (border_width * 2)image_height = int(800 * height_scaling) + (border_width * 2)column_width = width_scalingcolumn_bottom = image_height - border_widthx = image_width - border_widthimage = Image.new("RGB", (image_width, image_height), (32, 32, 32))image = _draw_labels(image, "Visible Spectrum", "Frequency (THz)", "Wavelength (nm)")_draw_x_axes(image, border_width, 384, 789, 50)_draw_y_axes(image, border_width, 0, 800, 50)draw = ImageDraw.Draw(image)for item in data:column_top = column_bottom - (item["nm"] * height_scaling)draw.rectangle(xy=[x, column_bottom, x - column_width, column_top],fill=item["rgb"])x -= column_widthtry:image.save(filename, "PNG")except IOError as e:print(e)def _draw_labels(image, heading_text, x_axis_text, y_axis_text):heading_font = ImageFont.truetype('Pillow/Tests/fonts/FreeSans.ttf', 32)axis_font = ImageFont.truetype('Pillow/Tests/fonts/FreeSans.ttf', 16)draw = ImageDraw.Draw(image)heading_text_size = draw.textsize(text=heading_text, font=heading_font)draw.text(xy=((image.width / 2)-(heading_text_size[0] / 2), 8),text=heading_text,align="center",font=heading_font,fill=(255, 255, 255))x_axis_text_size = draw.textsize(text=x_axis_text, font=axis_font)draw.text(xy=((image.width / 2) - (x_axis_text_size[0] / 2), image.height - 24),text=x_axis_text,font=axis_font,fill=(255, 255, 255))image = image.rotate(270, expand=1)draw = ImageDraw.Draw(image)y_axis_text_size = draw.textsize(text=y_axis_text, font=axis_font)draw.text(xy=((image.width / 2)-(y_axis_text_size[0] / 2), 8),text=y_axis_text,font=axis_font,fill=(255, 255, 255))image = image.rotate(90, expand=1)return imagedef _draw_y_axes(image, border_width, y_axis_start, y_axis_end, y_axis_interval):y_axis_indices_x_left = border_width - 8y_axis_indices_x_right = border_widthy = image.height - border_widthy_distance = ((image.height - (border_width * 2)) / (y_axis_end - y_axis_start)) * y_axis_intervalindex_font = ImageFont.truetype('Pillow/Tests/fonts/FreeSans.ttf', 12)draw = ImageDraw.Draw(image)for v in range(y_axis_start, y_axis_end + 1, y_axis_interval):draw.line(xy=[y_axis_indices_x_left, y, y_axis_indices_x_right, y],fill=(255, 255, 255),width=1)v_str = str(v)v_str_size = draw.textsize(text=v_str, font=index_font)draw.text(xy=[y_axis_indices_x_left - 2 - (v_str_size[0]), y - (v_str_size[1] / 2)],text=v_str,font=index_font,fill=(255, 255, 255))y -= y_distancedef _draw_x_axes(image, border_width, x_axis_start, x_axis_end, x_axis_interval):x_axis_indices_y_top = image.height - border_widthx_axis_indices_y_bottom = x_axis_indices_y_top + 8x = border_widthx_distance = ((image.width - (border_width * 2)) / (x_axis_end - x_axis_start)) * x_axis_intervalindex_font = ImageFont.truetype('Pillow/Tests/fonts/FreeSans.ttf', 12)draw = ImageDraw.Draw(image)for v in range(x_axis_start, x_axis_end + 1, x_axis_interval):draw.line(xy=[x, x_axis_indices_y_bottom, x, x_axis_indices_y_top],fill=(255, 255, 255),width=1)v_str = str(v)v_str_size = draw.textsize(text=v_str, font=index_font)draw.text(xy=[x - (v_str_size[0] / 2), x_axis_indices_y_bottom + 2],text=v_str,font=index_font,fill=(255, 255, 255))x += x_distancedef _wavelength_to_rgb(nm):gamma = 0.8max_intensity = 255factor = 0rgb = {"R": 0, "G": 0, "B": 0}if 380 <= nm <= 439:rgb["R"] = -(nm - 440) / (440 - 380)rgb["G"] = 0.0rgb["B"] = 1.0elif 440 <= nm <= 489:rgb["R"] = 0.0rgb["G"] = (nm - 440) / (490 - 440)rgb["B"] = 1.0elif 490 <= nm <= 509:rgb["R"] = 0.0rgb["G"] = 1.0rgb["B"] = -(nm - 510) / (510 - 490)elif 510 <= nm <= 579:rgb["R"] = (nm - 510) / (580 - 510)rgb["G"] = 1.0rgb["B"] = 0.0elif 580 <= nm <= 644:rgb["R"] = 1.0rgb["G"] = -(nm - 645) / (645 - 580)rgb["B"] = 0.0elif 645 <= nm <= 780:rgb["R"] = 1.0rgb["G"] = 0.0rgb["B"] = 0.0if 380 <= nm <= 419:factor = 0.3 + 0.7 * (nm - 380) / (420 - 380)elif 420 <= nm <= 700:factor = 1.0elif 701 <= nm <= 780:factor = 0.3 + 0.7 * (780 - nm) / (780 - 700)if rgb["R"] > 0:rgb["R"] = int(max_intensity * ((rgb["R"] * factor) ** gamma))else:rgb["R"] = 0if rgb["G"] > 0:rgb["G"] = int(max_intensity * ((rgb["G"] * factor) ** gamma))else:rgb["G"] = 0if rgb["B"] > 0:rgb["B"] = int(max_intensity * ((rgb["B"] * factor) ** gamma))else:rgb["B"] = 0return (rgb["R"], rgb["G"], rgb["B"])

generate_data

This function creates a list of dictionaries, one for each of the integer wavelengths between 380nm and 780nm. Each dictionary holds the wavelength, frequency in Hertz and TeraHertz, the RGB value of the colour, and the energy in Joules and 10 -19 Joules. The list can then be passed to the following functions which print the data to the terminal, and create graphs from it.

Firstly we create an empty list and a couple of constants, c for the speed of light and h for Planck's Constant.

Then we iterate through the required wavelengths, creating a new dictionary and adding it to the list. The frequencies and energies are calculated as described above, and the colours as RGB values are calculated by the _wavelength_to_rgb function which I'll get to in a moment.

print_data

The print_data function takes a list as created by generate_data and prints it out in a neat table format. The peculiar pow_minus19 variable is initialised to -19 for use after the E values. I have done this as otherwise Python would use the ugly e-19 notation.

plot_wavelength_frequency

The first plotting function takes our data structure and plots it with wavelengths on the x-axis and frequencies on the y-axis. The wavelengths fall between 380 and 780, and the frequencies are plotted from between 0 and 800. To give a more appealing shape I am stretching it to twice the width and half the height using scaling factors which can be applied to the wavelengths and frequencies to get pixel values.

After setting up the necessary variables we create a Pillow image, and then pass it to three functions to draw the labels and axes. The images is then used to get an ImageDraw object; this provides a set of methods to draw on the image.

Next we iterate the data, calculating the y coordinate of the top of each column before drawing it in its correct colour. Finally the image is saved.

plot_frequency_wavelength

This functions works in the same way except that frequencies are along the x-axis and wavelengths are on the y-axis.

_draw_labels

This function draws a heading at the top, as well as labelling the x and y axes.

A major drawback with Pillow's drawing functionality is that it does not use installed system fonts and we have to provide the location of a font file. I copied the file FreeSans.ttf to the project directory and then hard-coded the filename. You might wish to do the same, or hard code the path to an installed font. I didn't discover this until I had written most of the graphing code; if I'd found out sooner I would have abandoned Pillow and used something else, probably SVG. Anyhow, I won't be using Pillow for anything that requires drawing text for any further projects.

Having created a couple of fonts we draw the specified pieces of text. Another drawback with Pillow is that there is no built-in way of drawing text at an angle so I had to resort to rotating the image, drawing the text, and then rotating it round to its correct orientation.

_draw_y_axes

_draw_x_axes

These functions draw the little lines along each axis, along with their corresponding numeric values.

_wavelength_to_rgb

This function is based on code originally written in Fortran by Dan Bruton, and now available here: http://www.midnightkite.com/color.html . There are a number of implementations floating around in various languages (Pascal, R etc.), often with the logic tweaked slightly.

Firstly we create some constants and a dictionary for the RGB values. Next comes a pile of if/elifs, setting the RGB values depending on which range the wavelength falls in.

Next the factor variable is set - this is used to make the values trail off at each end of the spectrum while leaving those in the middle unchanged with a factor of 1.0.

Finally the RGB values, which are currently between 0 and 1, are multiplied up to give values between 0 and 255, the factor is applied, and the result "toned down" by gamma .

This implementation is a straight translation from a Pascal version and is frankly a bit rough. Maybe one day I will go back to the code and tidy it up a bit...

Now we just need a short program to demo our module.

visiblespectrum_demo.py

import visiblespectrumdef main():print("--------------------")print("| codedrome.com    |")print("| Visible Spectrum |")print("--------------------\n")data = visiblespectrum.generate_data()visiblespectrum.print_data(data)visiblespectrum.plot_wavelength_frequency(data, "wavelength_frequency.png")visiblespectrum.plot_frequency_wavelength(data, "frequency_wavelength.png")main()

As you can see we just call generate_data and then pass the result to print_data , plot_wavelength_frequency and plot_frequency_wavelength . Now run the program.

Run

python3.7 visiblespectrum_demo.py

The console output is

Program Output (partial)

--------------------| codedrome.com    || Visible Spectrum |---------------------------------------------------------|λ(nm)|f(THz)|E(J)      |R  |G  |B  |-------------------------------------|  380|   789|5.23x10⁻¹⁹| 97|  0| 97||  381|   787|5.22x10⁻¹⁹|100|  0|101||  382|   785|5.20x10⁻¹⁹|103|  0|106||  383|   783|5.19x10⁻¹⁹|106|  0|110||  384|   781|5.18x10⁻¹⁹|108|  0|115||  385|   779|5.16x10⁻¹⁹|111|  0|119||  386|   777|5.15x10⁻¹⁹|113|  0|123||  387|   775|5.14x10⁻¹⁹|115|  0|127||  388|   773|5.12x10⁻¹⁹|117|  0|132||  389|   771|5.11x10⁻¹⁹|119|  0|136||  390|   769|5.10x10⁻¹⁹|121|  0|140||  391|   767|5.08x10⁻¹⁹|123|  0|144||  392|   765|5.07x10⁻¹⁹|124|  0|148||  393|   763|5.06x10⁻¹⁹|125|  0|152||  394|   761|5.05x10⁻¹⁹|126|  0|156||  395|   759|5.03x10⁻¹⁹|127|  0|160||  396|   757|5.02x10⁻¹⁹|128|  0|164||  397|   755|5.01x10⁻¹⁹|129|  0|168||  398|   753|4.99x10⁻¹⁹|129|  0|172||  399|   751|4.98x10⁻¹⁹|130|  0|176||  400|   749|4.97x10⁻¹⁹|130|  0|180||  401|   748|4.96x10⁻¹⁹|130|  0|184||  402|   746|4.94x10⁻¹⁹|130|  0|188||  403|   744|4.93x10⁻¹⁹|130|  0|192|...

It might strike you that the energy values in Joules are pretty small, but remember that they are the energy for a single photon, the smallest possible amount of electromagnetic energy. I think maybe I should have converted these value to attojoules, ie 10 -18 J.

Wikipedia has a few examples to give an intuitive idea of a Joule, my favourite being "The energy required to lift a medium-sized tomato up 1 metre". Not surprisingly, even a relatively powerful violet photon can only lift a tomato a very short distance.

If you open the folder where you saved the source code you'll find these two PNG files.

This one shows wavelengths increasing left to right and their corresponding colours: violet, indigo, blue, green, yellow, orange and red. Note that the frequencies decrease as the wavelengths increase.

This is the same data but with frequencies on the x-axis, increasing left to right. This of course means the colours are reversed.

Follow @CodeDrome