Scripts to produce publication-ready figures.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

173 lines
6.6 KiB

  1. import os, csv
  2. from pathlib import Path
  3. import numpy as np
  4. import pandas as pd
  5. import seaborn as sns
  6. import matplotlib.pyplot as plt
  7. sns.set_theme(style="ticks")
  8. def parse_NIST_data(path, min_J, max_J, max_wavenumber):
  9. data_list = []
  10. with open(path) as f:
  11. reader = csv.reader(f, delimiter="\t")
  12. data_list = list(reader)
  13. #collect data
  14. Parity = np.zeros(len(data_list))
  15. J = np.zeros(len(data_list))
  16. Wavenumber = np.zeros(len(data_list))
  17. for i in range(1, len(data_list)):
  18. try:
  19. tmp = data_list[:][i]
  20. if not tmp[0] == '' and tmp[0].find('*') != -1:
  21. Parity[i] = 1
  22. elif not tmp[0] == '':
  23. Parity[i] = 0
  24. J[i] = int(tmp[1])
  25. Wavenumber[i] = float(tmp[3])
  26. except ValueError:
  27. J[i] = np.nan
  28. Wavenumber[i] = np.nan
  29. remove_idxs = []
  30. for i in range(1, len(data_list)):
  31. p = Parity[i]
  32. j = J[i]
  33. wn = Wavenumber[i]
  34. if np.isnan(p) or np.isnan(j) or np.isnan(wn):
  35. remove_idxs.append(i)
  36. Parity = np.delete(Parity, remove_idxs)
  37. J = np.delete(J, remove_idxs)
  38. Wavenumber = np.delete(Wavenumber, remove_idxs)
  39. #sort data
  40. sorting_indices = np.argsort(J)
  41. Parity = Parity[sorting_indices]
  42. J = J[sorting_indices]
  43. Wavenumber = Wavenumber[sorting_indices]
  44. # splice data to within user-defined range of Js
  45. splice_idx_start = np.where(J==min_J)[0][0]
  46. splice_idx_stop = len(J) - 1 - np.where(J[::-1]==max_J)[0][0]
  47. Parity = Parity[splice_idx_start:splice_idx_stop]
  48. J = J[splice_idx_start:splice_idx_stop]
  49. Wavenumber = Wavenumber[splice_idx_start:splice_idx_stop]
  50. # splice data to within user-defined range of Wavenumbers
  51. splice_idxs = [i for i in range(len(Wavenumber)) if Wavenumber[i] > max_wavenumber]
  52. Parity = [ele for idx, ele in enumerate(Parity) if idx not in splice_idxs]
  53. J = [ele for idx, ele in enumerate(J) if idx not in splice_idxs]
  54. Wavenumber = [ele for idx, ele in enumerate(Wavenumber) if idx not in splice_idxs]
  55. # Create a Pandas data frame with the data
  56. dataset = pd.DataFrame(np.array(list(zip(Parity, J, Wavenumber))), columns=['Parity', 'J', 'Wavenumber'])
  57. return dataset
  58. def plot_level_structure_with_red_and_blue_transitions(*args, **kwargs):
  59. #start plotting
  60. f, ax = plt.subplots(figsize=(4, 8))
  61. named_colors = ['r', 'm']
  62. Red_Blue_colors = ['#ab162a', '#cf5246', '#eb9172', '#fac8af', '#faeae1', '#e6eff4', '#bbdaea', '#7bb6d6', '#3c8abe', '#1e61a5']
  63. #draw levels
  64. plot_handle = sns.scatterplot(x='J', y='Wavenumber', data = dataframe, s=2000, hue = 'Parity', palette = sns.color_palette(named_colors), marker = '_', linewidth=1.5, legend=False)
  65. #write electronic configuration for GS
  66. ax.text(gs_J + 0.15, gs_wavenumber + 400, '$6s^2$')
  67. #draw guide line for GS
  68. #plt.axhline(y=gs_wavenumber, color='m', linestyle='--', linewidth=1, alpha=0.5)
  69. #write wavelength of red transition
  70. ax.text(red_J - 0.4, red_wavenumber * 0.5, '$626.082 ~ \mathrm{nm}$', color = '#db2929')
  71. ax.text(red_J - 0.4, red_wavenumber * 0.46, '$(\\Gamma = 2\\pi\\times 136 ~ \mathrm{kHz})$', fontsize = 8, color = '#db2929')
  72. #draw red transition arrow
  73. ax.annotate('',
  74. xy=(red_J, red_wavenumber),
  75. xytext=(gs_J, gs_wavenumber),
  76. arrowprops=dict(color='#db2929', alpha=0.8, width=1.5),
  77. horizontalalignment='right',
  78. verticalalignment='top')
  79. #write electronic configuration for triplet excited state
  80. ax.text(red_J + 0.35, red_wavenumber + 200, '$6s6p(^3P_1)$', fontsize = 10)
  81. #draw guide line for triplet excited state
  82. plt.axhline(y=red_wavenumber, color='m', linestyle='--', linewidth=1, alpha=0.5)
  83. #write wavelength of red transition
  84. ax.text(blue_J - 1.5, blue_wavenumber * 0.55, '$421.291~ \mathrm{nm}$', color = '#2630ea')
  85. ax.text(blue_J - 1.55, blue_wavenumber * 0.52, '$(\\Gamma = 2\\pi\\times 32.2 ~ \mathrm{MHz})$', fontsize = 8, color = '#2630ea')
  86. #draw blue transition arrow
  87. ax.annotate('',
  88. xy=(blue_J, blue_wavenumber),
  89. xytext=(gs_J, gs_wavenumber),
  90. arrowprops=dict(color='#2630ea', alpha=0.8, width=3.5),
  91. horizontalalignment='right',
  92. verticalalignment='top')
  93. #write electronic configuration for singlet excited state
  94. ax.text(blue_J + 0.35, blue_wavenumber + 200, '$6s6p(^1P_1)$', fontsize = 10)
  95. #draw guide line for singlet excited state
  96. plt.axhline(y=blue_wavenumber, color='m', linestyle='--', linewidth=1, alpha=0.5)
  97. #figure options
  98. f.canvas.draw()
  99. plt.xlabel('$J$', fontsize=16)
  100. plt.ylabel('$\\tilde{v}~(cm^{-1})$', fontsize=16)
  101. plt.ylabel('$\\lambda~(nm)$', fontsize=16)
  102. plot_handle.set_xticks(range(min_J-1, max_J+2))
  103. ax.get_xticklabels()[0].set_visible(False)
  104. ax.get_xticklabels()[-1].set_visible(False)
  105. ax.get_xticklines()[0].set_visible(False)
  106. ax.get_xticklines()[-2].set_visible(False)
  107. yticklabels = [item.get_text() for item in ax.get_yticklabels()]
  108. yticklabels = ['' if item.startswith('') or item.startswith('0') else item for item in yticklabels]
  109. yticks = [float(item) if item != '' else 0.0 for item in yticklabels]
  110. new_yticks = np.arange(min(yticks), max(yticks), 4000)
  111. plot_handle.set_yticks(new_yticks)
  112. new_yticklabels = [round(1e7/item) if item != 0 else item for item in new_yticks]
  113. ax.set_yticklabels(new_yticklabels)
  114. ax.get_yticklabels()[0].set_visible(False)
  115. ax.get_yticklabels()[-1].set_visible(False)
  116. ax.get_yticklines()[0].set_visible(False)
  117. ax.get_yticklines()[-2].set_visible(False)
  118. plt.tick_params(axis='both', which='major', labelsize=14)
  119. #plt.show()
  120. f.savefig(Path(home_path + os.sep + 'result.pdf'), format='pdf', bbox_inches = "tight")
  121. if __name__ == '__main__':
  122. min_J = 8
  123. max_J = 9
  124. max_wavenumber = 25500.0
  125. #Ground State: [Xe]4f^{10} 6s^2(^5I_8)
  126. gs_J = 8
  127. gs_wavenumber = 0.0
  128. #626 nm transition GS ---> [Xe]4f^{10}(^5I_8) 6s6p(^3P_1)(8,1)_9
  129. red_J = 9
  130. red_wavenumber = 15972.35
  131. #421 nm transition GS ---> [Xe]4f^{10}(^5I_8) 6s6p(^1P_1)(8,1)_9
  132. blue_J = 9
  133. blue_wavenumber = 23736.610
  134. NIST_Dy_level_data_filename= 'Dylevels.txt'
  135. home_path = str(Path(__file__).parent.resolve())
  136. dataframe = parse_NIST_data(home_path + os.sep + NIST_Dy_level_data_filename, min_J, max_J, max_wavenumber)
  137. plot_level_structure_with_red_and_blue_transitions(dataframe, gs_J, gs_wavenumber, red_J, red_wavenumber, blue_J, blue_wavenumber, home_path)