Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Aurora processing of spectra #45

Open
fsciortino opened this issue Dec 21, 2021 · 0 comments
Open

Aurora processing of spectra #45

fsciortino opened this issue Dec 21, 2021 · 0 comments
Assignees
Labels
tutorials Issues intended to demonstrate Aurora features

Comments

@fsciortino
Copy link
Owner

This is a short tutorial intended to answer questions from @Didou09, but posted here publicly hoping that it may be useful to other users in the future.

Here we demonstrate how to plot the He-like Ar x-ray spectrum using Aurora. We will make use of some photon emissivity coefficients (PECs) that are not distributed by ADAS (it is originally from atomDB, another database), but that can be requested by emailing me. Here we focus on the He-like Ar spectrum, discussed for example in J E Rice et al 2015 J. Phys. B: At. Mol. Opt. Phys. 48 144013 (https://iopscience.iop.org/article/10.1088/0953-4075/48/14/144013), which is analogous to the He-like Ca spectrum discussed in F. Sciortino et al 2021 Nucl. Fusion 61 126060 (https://iopscience.iop.org/article/10.1088/1741-4326/ac32f2). The same Aurora functionality can be used for any spectra for which ADAS data is available. Important: it is up to the user to assess the atomic data quality. Ask an expert for help if you're in doubt.

Let's begin with a few basic setup tasks, run in IPython:

import numpy as np, sys, os
import matplotlib.pyplot as plt
plt.ion()
from scipy.interpolate import interp1d
import aurora

# just some preferences to better see the spectra later, but not necessary
import matplotlib as mpl
mpl.rcParams['axes.titlesize'] = 24
mpl.rcParams['axes.labelsize'] = 20
mpl.rcParams['lines.linewidth'] = 1.5
mpl.rcParams['lines.markersize'] = 10
mpl.rcParams['xtick.labelsize'] = 16
mpl.rcParams['ytick.labelsize'] = 16
mpl.rcParams['legend.fontsize'] = 16

Next, we specify the PEC files of interest -- modify these appropriately to point to your files of interest. These must be in an ADF15 ADAS format.

Z_ion=18
filepath_h = 'pec#ar17.dat'
filepath_he = 'pec#ar16.dat'
filepath_li = 'pec#ar15.dat'

Here we referred to each file as "_h", "_he" and "_li" because these files correspond to the H-like, He-like and Li-like states of Ar ions. We then specify values of electron density and temperatures at which we want to compute the spectrum. We are going to pick the following values here, with units clearly specified in the variables' names:

Te_eV = 3.5e3
ne_cm3 = 1e14

To balance the abundance of Li-like, He-like and H-like Ar charge states, we must choose a model. In general, the density of these charge states will depend on both atomic physics (ionization and recombination balance) and transport physics. Here we are going to ignore the latter component for simplicity (although Aurora can simulate impurity transport in 1.5D) and just consider fractional abundances at ionization equilibrium, computed from ADAS data as follows:

# first load atomic data for Ar:
atom_data = aurora.get_atom_data('Ar',['scd','acd'])  # accounting for effective ionization and recombination rates (ignore CX for simplicity)

# now compute fractional abundances -- this can be done for arrays of ne and Te (also including neutrals, if useful),
# but we are going to do it only for single values of ne and Te here
logTe, fz = aurora.get_frac_abundances(atom_data, np.array([ne_cm3,]), np.array([Te_eV,]), plot=False)

Then, we set up a figure and add the contributions to our spectrum of interest from each charge state:

fig, ax = plt.figure()
fig.set_size_inches(9,6, forward=True)

# value of Doppler shift -- here set to 0.15 mA for example
dlam_A = 0.00015

# now add spectra for each charge state -- see `get_local_spectrum` docs for details
out = aurora.get_local_spectrum(filepath_li, ion, ne_cm3, Te_eV, n0_cm3=0.0,
                                ion_exc_rec_dens=[fz[0,-5], fz[0,-4], fz[0,-3]], # Be-like, Li-like, He-like
                                dlam_A = dlam_A, plot_spec_tot=False, no_leg=True, plot_all_lines=True, ax=ax)
wave_final_li, spec_ion_li, spec_exc_li, spec_rec_li, spec_dr_li, spec_cx_li, ax = out

out= aurora.get_local_spectrum(filepath_he, ion, ne_cm3, Te_eV, n0_cm3=0.0,
                               ion_exc_rec_dens=[fz[0,-4], fz[0,-3], fz[0,-2]], # Li-like, He-like, H-like
                               dlam_A = dlam_A, plot_spec_tot=False, no_leg=True, plot_all_lines=True, ax=ax)
wave_final_he, spec_ion_he, spec_exc_he, spec_rec_he, spec_dr_he, spec_cx_he, ax = out

out = aurora.get_local_spectrum(filepath_h, ion, ne_cm3, Te_eV, n0_cm3=0.0,
                                ion_exc_rec_dens=[fz[0,-3], fz[0,-2], fz[0,-1]], # He-like, H-like, fully stripped
                                dlam_A = dlam_A, plot_spec_tot=False, no_leg=True, plot_all_lines=True, ax=ax)
wave_final_h, spec_ion_h, spec_exc_h, spec_rec_h, spec_dr_h, spec_cx_h, ax = out

# total spectra from all ions
spec_tot_li = spec_ion_li + spec_exc_li + spec_rec_li + spec_dr_li + spec_cx_li
spec_tot_he = spec_ion_he + spec_exc_he + spec_rec_he + spec_dr_he + spec_cx_he
spec_tot_h = spec_ion_h + spec_exc_h + spec_rec_h + spec_dr_h + spec_cx_h

# add plot of total spectrum
wave_all = np.linspace(plt.gca().get_xlim()[0], plt.gca().get_xlim()[1], 10000) #A
spec_all = interp1d(wave_final_li, spec_tot_li, bounds_error=False, fill_value=0.0)(wave_all)
spec_all += interp1d(wave_final_he, spec_tot_he, bounds_error=False, fill_value=0.0)(wave_all)
spec_all += interp1d(wave_final_h, spec_tot_h, bounds_error=False, fill_value=0.0)(wave_all)
plt.gca().plot(wave_all, spec_all, 'k', label='total')

# just to have nice legends:
ax2.plot([], [], c='r', ls='--', label='ionization')
ax2.plot([], [], c='b', ls='--', label='excitation')
ax2.plot([], [], c='g', ls='--', label='radiative recomb')
ax2.plot([], [], c='m', ls='--', label='dielectronic recomb')

This gives us a busy image like this:
drawing
Zooming into the region near 3.95 A, we see the familiar K-alpha spectrum of Ar (other spectral regions are also interesting, e.g. one can find the Ly-alpha lines using these files):

drawing

Here I also added some vertical dashed line that show the location and nomenclature of some standard lines (using Gabriel notation). See for example Sciortino et al RSI 92, 053508 (2021) (https://doi.org/10.1063/5.0043765) for a discussion of these components. The code above doesn't plot the vertical red-dashed lines, I only added them in my figure to help visualization.
@fsciortino fsciortino added the tutorials Issues intended to demonstrate Aurora features label Dec 22, 2021
@fsciortino fsciortino self-assigned this Dec 22, 2021
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
tutorials Issues intended to demonstrate Aurora features
Projects
None yet
Development

No branches or pull requests

1 participant