On my spectral analysis of a quasar project, I had cheated by already knowing the redshift of the quasar from professional research, and was able to reverse-discover what elements should be in my spectrum. It led me to wonder how I might discover the elements in my spectra naturally. The ratios of the observed lines is really all there is to go off of, but there are just too many combinations of possible elements observable in a galaxy to guess this by hand. So, I wrote a Python script to do this for me.
I used a table of emission lines observed in galaxies, converted to a spreadsheet. The script takes every combination of elements in this chart, computes the ratio of the wavelengths of the longer line to the shorter line, and compares it to the observed wavelengths entered in the command line. For example, using the data I had from my quasar project, I would use the script thusly,
% python3.8 redshift.py 445.333 654.016 730.710
You simply just run the script with every emission line mean wavelength from the command line. It will first show you the ratios of all the combinations of your wavelengths,
Observed wavelengths and ratios 1 4453.33000 2 6540.16000 6540.16000 / 4453.33000 = 1.468600 3 7307.10000 7307.10000 / 4453.33000 = 1.640817 7307.10000 / 6540.16000 = 1.117266
After it chugs through a large amount of possibilities from the reference data, it will show you the top ten matches based on the RMS of the errors between the observed and rest ratios, e.g., here are the top three from my data,
RMS: 0.000230 Rest Ratios: [[1.4684041497472178], [1.6410161306706743, 1.1175507308074357]] Elements: [Fe V]: 3839.270000 -> 4453.330000 [Fe VI]: 5637.600000 -> 6540.160000 [O I]: 6300.304000 -> 7307.100000 Redshift: 0.1599465122528052 RMS: 0.000426 Rest Ratios: [[1.4683580844685535], [1.6412941771876017, 1.1177751493646317]] Elements: C III]: 1908.734000 -> 4453.330000 Mg II]: 2802.705000 -> 6540.160000 O III: 3132.794000 -> 7307.100000 Redshift: 1.3330348690831522 RMS: 0.000835 Rest Ratios: [[1.467395237750597], [1.6406031906864216, 1.1180376959661793]] Elements: [Ne III]: 3868.760000 -> 4453.330000 [Fe VI]: 5677.000000 -> 6540.160000 Si II: 6347.100000 -> 7307.100000 Redshift: 0.15146512185773703
The Gaussian fits of my data were not great looking back at what I did, a noobie mistake, but the real emission lines come up as the second best fit! So if you emission lines are measured accurately enough, this script should give you the elements in your spectrum as well as the calculated mean redshift.
The script and data file can be downloaded here. You will need the
pandas library installed for the LibreOffice file operations. The full script is shown below,
# -*-coding: utf-8 -*- import sys import math import statistics from itertools import combinations import pandas as pd def rms(array): n = len(array) square = 0.0 root = 0.0 mean = 0.0 #calculating Squre for i in range(0, n): square += (array[i] ** 2) #Calculating Mean mean = (square/ (float)(n)) #Calculating Root root = math.sqrt(mean) return root def flatten_list(_2d_list): flat_list =  # Iterate through the outer list for element in _2d_list: if type(element) is list: # If the element is of type list, iterate through the sublist for item in element: flat_list.append(item) else: flat_list.append(element) return flat_list os = [float(a) for a in sys.argv[1:]] os.sort() ratios =  if os < 1000: os = [10 * o for o in os] print("Observed wavelengths and ratios") for i, w in enumerate(os): print("%i %0.5f" % (i + 1, w)) rs =  for j in range(0, i): rs.append(os[i] / os[j]) print(" %0.5f / %0.5f = %0.6f" % (os[i], os[j], os[i] / os[j])) if len(rs) > 0: ratios.append(rs) data = pd.read_excel('data.ods', engine='odf') print("") print("") errors =  max_row = 0 for i in range(0, len(data.index)): if os[-1] >= data.at[i, "λ (Å)"]: max_row = i combs = list(combinations(range(0, max_row), len(os))) for comb in combs: ratios_e =  for i in range(0, len(comb)): r_es =  for j in range(0, i): r_es.append(data.at[comb[i], "λ (Å)"] / data.at[comb[j], "λ (Å)"]) #print(" %0.5f / %0.5f = %0.6f" % (data.at[i, "λ (Å)"], data.at[j, "λ (Å)"], data.at[i, "λ (Å)"] / data.at[j, "λ (Å)"])) if len(r_es) > 0: ratios_e.append(r_es) # Use RMS of ratio errors. score = rms([o - e for o, e in zip(flatten_list(ratios), flatten_list(ratios_e))]) errors.append( (score, ratios_e, comb) ) errors = sorted(errors, key=lambda x: x) for i in range(0, 10): e = errors[i] print("RMS: %f" % (e)) print("Rest Ratios:") print(e) print("Elements:") comb = e lambda_es =  for j in range(0, len(comb)): lambda_es.append(data.at[comb[j], "λ (Å)"]) print("%s: %f -> %f" % (data.at[comb[j], "Ion"], data.at[comb[j], "λ (Å)"], os[j])) print("Redshift:") print(statistics.mean([(o - e) / e for o, e in zip(os, lambda_es)])) print("")