beta_shape_pnpi/notebooks/decay_paceENSDF.ipynb

paceENSDF methods for accessing ENSDF radioactive-decay data sets

Import the Python package paceENSDF and load the JSON-formatted ENSDF-decay data sets to a list object.

In [ ]:

import paceENSDF as pe
e = pe.ENSDF()
edata = e.load_ensdf()

In [ ]:

import numpy as np
import matplotlib.pyplot as plt
%matplotlib notebook

Find all radioactive-decay data pairs in ENSDF

Explore the available data sets and radioactive-decay data pairs, and associated meta data, using the ensdf_pairs method.

In [ ]:

# Find all radioactive-decay pairs in ENSDF
alpha_pairs = e.ensdf_pairs(edata,"A")
bm_pairs = e.ensdf_pairs(edata,"BM")
ecbp_pairs = e.ensdf_pairs(edata,"ECBP")
for k,v, in ecbp_pairs.items(): print(k,v)

Plot the nuclear chart relevant to the decay data

In [ ]:

# Gather the Z, A information from each decay-mode dictionary

# Alpha decay
Z_A = [k[1] for (k,v) in alpha_pairs.items()]
A_A = [k[2] for (k,v) in alpha_pairs.items()]

# Beta-minus decay
Z_BM = [k[1] for (k,v) in bm_pairs.items()]
A_BM = [k[2] for (k,v) in bm_pairs.items()]

# Electron-capture/beta-plus decay
Z_ECBP = [k[1] for (k,v) in ecbp_pairs.items()]
A_ECBP = [k[2] for (k,v) in ecbp_pairs.items()]

In [ ]:

%matplotlib notebook
fig, ax = plt.subplots(figsize=(9,6))

ax.scatter(np.array(A_A), np.array(Z_A), color='y', marker='s', alpha=0.5, label=r'$\alpha$')
ax.scatter(np.array(A_BM), np.array(Z_BM), color='r', marker='s', alpha=0.10, label=r'$\beta^{-}$')
ax.scatter(np.array(A_ECBP), np.array(Z_ECBP), color='b', marker='s', alpha=0.10, label=r'$\epsilon/\beta^{+}$')

ax.set_xlabel(r'$A$', size=20)
ax.set_ylabel(r'$Z$', size=20)
ax.tick_params(axis='both', which='major', labelsize=15)
ax.legend(loc='best', fontsize=20)

plt.tight_layout()
plt.savefig("decay_nuclides.png", dpi=fig.dpi)
plt.show()

pyENSDF: Manipulation of all primary- and continuation-record data from ENSDF

pyENSDF modules contain methods to access all data from:

• Parent record
• Normalization record
• Level record (+ continuation)
• Gamma record (+ continuation)
• Beta-minus record (+ continuation)
• Electron-capture/beta-plus record (+ continuation)
• Alpha record

Parent isotope properties

Methods associated with the properties of the parent isotope

(i) Find all isomeric parent isotopes in the ENSDF database associated with a particular decay mode

The returned quantities are explained in the docstrings:

>>> help(e.find_all_parent_isomers)

In [ ]:

# Find all alpha-decay isomers
e.find_all_parent_isomers(edata, mode="A")

In [ ]:

# Find all beta-minus-decay isomers
e.find_all_parent_isomers(edata, mode="BM")

In [ ]:

# Find all electron-capture/beta-plus-decay isomers
e.find_all_parent_isomers(edata, mode="ECBP")

(ii) Find all parents with multiple spin-parity permutations for a given decay mode

The returned quantities are explained in the docstrings:

>>> help(e.find_parents_multiple_jpi)

In [ ]:

# Find all alpha-decay parents with multiple JPi
e.find_parents_multiple_jpi(edata, mode="A")

In [ ]:

# Find all beta-minus-decay parents with multiple JPi
e.find_parents_multiple_jpi(edata, mode="BM")

In [ ]:

# Find all electron-capture/beta-plus-decay parents with multiple JPi
e.find_parents_multiple_jpi(edata, mode="ECBP")

(iii) Decay-level energy and Q-value for given parent and defined decay mode and index

NB: Ground-state decays have index 0; isomeric decays have index \geq1.

The returned quantities are explained in the docstrings:

>>> help(e.get_parent_decay)

In [ ]:

# 173Ir ground-state alpha decay
e.get_parent_decay(edata, "Ir173", 0, mode="A")

In [ ]:

# 102Tc isomer beta-minus decay
e.get_parent_decay(edata, "Tc102", 1, mode="BM")

In [ ]:

# 50V ground-state electron-capture/beta-plus decay
#e.get_parent_decay(edata, "V50", 0, mode="ECBP")
e.get_parent_decay(edata, "Dy155", 0, mode="ECBP")

(iv) Get the parent halflife for a defined isotope, decay mode and index

The returned quantities are explained in the docstrings:

>>> help(e.get_parent_halflife)

In [ ]:

# Beta-minus ground-state decay of 60Co parent (T1/2 in best units)
e.get_parent_halflife(edata, "Co60", 0, mode="BM", units="seconds")

In [ ]:

# Beta-minus isomer decay of 60Co parent (T1/w in best units)
e.get_parent_halflife(edata, "Co60", 1, mode="BM", units="best")

In [ ]:

# Beta-minus isomer decay of 60Co parent (T1/w in seconds units)
e.get_parent_halflife(edata, "Co60", 1, mode="BM", units="seconds")

(v) Get spin and parity of a given parent isotope for a defined decay mode and index

The returned quantities are explained in the docstrings:

>>> help(e.get_parent_jpi)

In [ ]:

# 258Db ground state electron-capture/beta-plus decay (parent with two JPi values)
#e.get_parent_jpi(edata, "Db258", 0, mode="ECBP")
e.get_parent_jpi(edata, "Dy155", 0, mode="ECBP")

In [ ]:

# 60Co isomer beta-minus decay (parent with single JPi value)
e.get_parent_jpi(edata, "Co60", 1, mode="BM")

In [ ]:

# 60Co ground state beta-minus decay (parent with single JPi value)
e.get_parent_jpi(edata, "Co60", 0, mode="BM")

Normalization and Production Normalization records

Refer to the docstrings for more information regarding the returned qauntities:

>>> help(e.norm_record)
>>> help(e.prod_norm_record)

In [ ]:

# Normalization record for 60Co (g.s.)
#e.norm_record(edata, "Co60", 0, mode="BM")
e.norm_record(edata, "Dy155", 0, mode="ECBP")

In [ ]:

# Production Normalization record for 60Co (g.s.)
#e.prod_norm_record(edata, "Co60", 0, mode="BM")
e.prod_norm_record(edata, "Dy155", 0, mode="ECBP")

Decay scheme levels and gammas

Methods for accessing the levels and gammas, together with their associated properties, populated in the daughter nucleus following radioactive decay.

The following keyword arguments are accepted and passed to callables requiring the radioactive-decay mode:

mode = "A" # alpha decay
mode = "BM" # beta-minus decay
mode = "ECBP" # electron-capture/beta-plus decay

(i) Decay scheme of daughter nucleus populated in decay

The returned quantities are explained in the docstrings:

>>> help(e.get_levels)

In [ ]:

# Levels populated in daughter nucleus following radioactive decay
#e.get_levels(edata,"Ra226",0,mode="A" )
e.get_levels(edata,"Dy155",0,mode="ECBP" )

(ii) Find all levels of a particular nucleus with unique spin-parity assignments

The returned quantities are explained in the docstrings:

>>> help(e.find_unique_jpi)

In [ ]:

# For the levels with unique JPi assignments
#e.find_unique_jpi(edata, "Ar45", 0, mode="BM")
e.find_unique_jpi(edata, "Dy155", 0, mode="ECBP")

(iii) Find all levels of a particular nucleus with unique spin-parity assignments

The returned quantities are explained in the docstrings:

>>> help(e.find_multiple_jpi)

In [ ]:

# Find the levels that have multiple JPi assignments
levels = e.find_multiple_jpi(edata, "Dy155", 0, mode="ECBP")
for k,v in levels.items():
    print("Parent = {0} ({1} keV); Daughter = {2}".format(k[0],k[4],k[5]))
    for vv in v:
        print("index = {0}; energy = {1} \xb1 {2} keV; J = {3}; Pi = {4}".
              format(vv[0],vv[1],vv[2],vv[5],vv[6]))

(iv) Find all isomeric levels populated in the daughter nucleus

The returned quantities are explained in the docstrings:

>>> help(e.find_isomers)

The halflife information is returned in its 'best' units or in units of 'seconds' depending on the keyword argument:

units = "best"
units = "seconds"

In [ ]:

# Find all isomers in daughter nucleus following radioactive decay ('best' or 's' units for halflife)
isomers = e.find_isomers(edata, "Dy155", 0, mode="ECBP", units="best")
for k,v in isomers.items():
    print("Parent = {0} ({1} keV); Daughter = {2}".format(k[0],k[4],k[5]))
    for vv in v: 
        print("index = {0}; energy = {1} \xb1 {2} keV; T1/2 = {3} \xb1 {4} {5}".
              format(vv[0], vv[1], vv[2], vv[3], vv[4], vv[5]))

(v) Get all levels, gammas, and associated properties of the daughter nucleus

The returned quantities are explained in the docstrings:

>>> help(e.get_levels_and_gammas)

In [ ]:

# Levels and gammas associated with daughter following radioactive decay
lg = e.get_levels_and_gammas(edata,"Dy155",0,mode="ECBP")
#for i in lg: print(i)
print("Transition: {0} -> {1}; {2} keV -> {3} keV".format(lg[0][0],lg[0][1],lg[0][2],lg[0][3]))
print("Gamma-ray energy: {0} \xb1 {1} keV".format(lg[0][6],lg[0][7]))
print("Gamma-ray multipolarity: {0}".format(lg[0][10]))
print("Gamma-ray mixing ratio: {0} \xb1 {1}".format(lg[0][11],lg[0][12]))
print("Raw gamma-ray intensity: {0} \xb1 {1}".format(lg[0][8],lg[0][9]))
print("Total ICC: {0}".format(lg[0][14]))

Atomic subshell information can also be retrieved by passing one of the following keyword arguments to the get_levels_and_gammas_subshells method:

subshell = 'calc'
subshell = 'sumcalc'
subshell = 'ratio'
subshell = 'expt'
subshell = 'electron'

Refer to the docstrings for more information:

>>> help(e.get_levels_and_gammas_subshells)

In [ ]:

# Level and gamma-decay information as before, with additional information for the calculated atomic 
# subshell internal conversion coefficients (K to Q shell)
lgs = e.get_levels_and_gammas_subshells(edata,"Dy155",0,mode="ECBP",subshell='calc')
#for i in lgs: print(i)
print("Gamma-ray energy: {0} \xb1 {1} keV".format(lgs[0][6],lgs[0][7]))
print("Total ICC: {0}".format(lgs[0][14]))
print("K-shell ICC: {0} \xb1 {1}".format(lgs[0][16],lgs[0][17]))
print("L-shell ICC: {0} \xb1 {1}".format(lgs[0][18],lgs[0][19]))
print("M-shell ICC: {0} \xb1 {1}".format(lgs[0][20],lgs[0][21]))
print("N-shell ICC: {0} \xb1 {1}".format(lgs[0][22],lgs[0][23]))
print("O-shell ICC: {0} \xb1 {1}".format(lgs[0][24],lgs[0][25]))
print("P-shell ICC: {0} \xb1 {1}".format(lgs[0][26],lgs[0][27]))
print("Q-shell ICC: {0} \xb1 {1}".format(lgs[0][28],lgs[0][29]))

$\beta^{-}$-decay properties

(i) Get $\beta^{-}$-decay information to all states populated in the decay

The returned quantities are explained in the docstrings:

>>> help(e.get_beta_minus)

In [ ]:

# Properties of states observed following beta-minus decay of 60Co ground state (best halflife units):
x=e.get_beta_minus(edata, "Co60", 0, units='best')
print(x)

In [ ]:

# Properties of states observed following beta-minus decay of 60Co ground state (seconds halflife units):
x=e.get_beta_minus(edata, "Co60", 0, units='seconds')
print(x)

(ii) Get log(ft) values from \beta^{-} decay or \epsilon/\beta^{+} decay

The returned quantities are explained in the docstrings:

>>> help(e.get_logft)

In [ ]:

# Retrieve log(ft) data following beta-minus decay of 60Co g.s.
e.get_logft(edata,"Co60",0,mode='BM')

In [ ]:

# Retrieve log(ft) data following electron-capture/beta-plus decay of 86Y g.s.
e.get_logft(edata,"Y86",0,mode='ECBP')

(iii) Find nuclides of a given forbiddenness in \beta^{-} decay and \epsilon/\beta^{+} decay

The forbiddenness classification needs to be passed as a string argument in these functions:

'0A': Allowed (S=0)
'1A' : Allowed (S=1)
'1F': First forbidden
'1UF': First-forbidden unique
'2F': Second forbidden
'2UF': Second-forbidden unique
'3F': Third forbidden
'3UF': Third-forbidden unique
'4F': Fourth forbidden
'4UF': Fourth-forbidden unique
'5F': Fifth forbidden
'5UF': Fifth-forbidden unique

The callable find_forbidden takes either three arguments to scan the entire database, or five arguments to search for a given nucleus. For more information refer to the docstrings:

>>> help(e.find_forbidden)

In [ ]:

# To find all '2UF' transitions throughout the ENSDF beta-minus decay data sets
x=e.find_forbidden(edata,'2UF',mode='ECBP')
print(x)

In [ ]:

# To find all '3F' transitions throughout the ENSDF electron-capture/beta-plus decay data sets
x=e.find_forbidden(edata,'3F',mode='ECBP')
print(x)

In [ ]:

# Find all '2UF' transitions in 60Co beta-minus decay
x=e.find_forbidden(edata,"Co60",0,'2UF',mode='BM')
print(x)

In [ ]:

# Find all '1UF' transitions in 158Tm electron-capture/beta-plus decay
x=e.find_forbidden(edata,"Tm158",0,'1UF',mode='ECBP')
print(x)

$\epsilon/\beta^{+}$-decay properties

(i) $\epsilon/\beta^{+}$-decay properties to all states populated following the decay of the parent nucleus

The returned quantities are explained in the docstrings:

>>> help(e.get_ecbp)

In [ ]:

# Get electon-capture/beta-plus decay properties from the ground-state decay of 86Y
#x=e.get_ecbp(edata, "Y86", 0, units='best')
x=e.get_ecbp(edata, "Dy155", 0, units='best')
for k,v in x.items():
    print(v)

(ii) $\epsilon$-decay fractions to all states populated following the decat of the parent nucleus

The returned quantities are explained in the docstrings:

>>> help(e.get_ec_fractions)

NB: One of the following keyword arguments is needed for the desired atomic-shell \epsilon data:

subshell = 'calc'
subshell = 'sumcalc'
subshell = 'expt'

In [ ]:

# Get calculated electron-capture decay fractions to states populated in 86Sr following the g.s. decay of 86Y
#x=e.get_ec_fractions(edata, "Y86", 0, subshell='calc')
x=e.get_ec_fractions(edata, "Dy155", 0, subshell='calc')
print(x)

$\alpha$-decay properties

The returned quantities are explained in the docstrings:

>>> help(e.get_alpha)

In [ ]:

# Get alpha-decay properties to all states populated in Rn222 following the ground-state decay of Ra226
alpha=e.get_alpha(edata, "Ra226", 0, units='best')
print(alpha)

In [ ]:

# Put alpha intensity and HF data into lists for plotting
level_energy = []
IA, dIA = [], []
HF, dHF = [], []
for k,v in alpha.items():
    print("Daughter = {0}".format(k[3]))
    for vv in v:
        print("level energy = {0}; IA = {1} \xb1 {2}; HF = {3} \xb1 {4}".format(vv[1],vv[9],vv[10],vv[11],vv[12]))
        level_energy.append(vv[1])
        IA.append(vv[9])
        dIA.append(vv[10])
        HF.append(vv[11])
        dHF.append(vv[12])

In [ ]:

# Plot IA vs. E and IA vs. HF
%matplotlib notebook
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(9,6))

ax1.errorbar(level_energy, IA, xerr=None, yerr=dIA, color='b', fmt='o', capsize=5, label=r'$^{226}$Ra $\rightarrow$ $^{222}$Rn + $\alpha$')
ax1.legend(loc='best', fontsize=15)
ax1.set_xlabel(r'$E$ [keV]', size=20)
ax1.set_ylabel(r'$I_{\alpha}$', size=20)

#ax1.set_xscale('log')
ax1.set_yscale('log')

ax2.errorbar(HF, IA, xerr=dHF, yerr=dIA, color='b', fmt='o', capsize=5, label=r'$^{226}$Ra $\rightarrow$ $^{222}$Rn + $\alpha$')
ax2.legend(loc='best', fontsize=15)
ax2.set_xlabel(r'$HF$', size=20)
ax2.set_ylabel(r'$I_{\alpha}$', size=20)

#ax2.set_xscale('log')
ax2.set_yscale('log')

ax1.grid()
ax2.grid()

plt.tight_layout()
plt.savefig("E_hf_plot.png", dpi=fig.dpi)
plt.show()

In [ ]:

np.corrcoef(np.array(IA), np.array(HF))

In [ ]:

# Plot IA vs. HF
%matplotlib notebook
fig, (ax1,ax2) = plt.subplots(2,figsize=(9,6),sharex=True)

ax2.errorbar(level_energy, IA, xerr=None, yerr=dIA, color='b', fmt='o', capsize=5, label=r'$^{226}$Ra $\rightarrow$ $^{222}$Rn + $\alpha$')
#ax2.legend(loc='best', fontsize=15)
ax2.set_xlabel(r'$E$ [keV]', size=20)
ax2.set_ylabel(r'$I_{\alpha}$', size=20)
ax2.set_yscale('log')

ax1.errorbar(level_energy, HF, xerr=None, yerr=dHF, color='b', fmt='o', capsize=5, label=r'$^{226}$Ra $\rightarrow$ $^{222}$Rn + $\alpha$')
ax1.legend(loc='best', fontsize=15)
ax1.set_ylabel(r'$HF$', size=20)
ax1.set_yscale('log')

ax1.grid()
ax2.grid()

plt.tight_layout()
plt.savefig("HF_IA_vs_E_plot.png", dpi=fig.dpi)
plt.show()

Search for superallowed decays in \epsilon/\beta^{+} decay

In [ ]:

# Superallowed 0+->0+ decays in EC/BP nuclides

superallowed = []
for k,v, in ecbp_pairs.items():
    decay = e.get_ecbp(edata, str(k[0]), int(k[3]), units='best')
    for k_decay, v_decay in decay.items():
        # Firm J and pi
        if int(k_decay[9])==1 and int(k_decay[11])==1:
            # Find parents with Jpi=0+
            if int(2*k_decay[8])==0 and int(k_decay[10])==1:
                # Find daughter with Jpi=0+
                for daughter in v_decay:
                    # Search for ground states only
                    #if int(daughter[0])==0:
                    # Firm J and pi:
                    if int(daughter[2]==1) and int(daughter[4])==1 and int(daughter[6])==1:
                        # Find J=0 and pi=+:
                        try:
                            if int(2*daughter[3])==0 and int(daughter[5])==1 and daughter[17]<4.0:
                                print(k[0], k[1], k[2], k[3], k[4], k_decay[3], k_decay[4], \
                                      k_decay[5], daughter[0], daughter[1], daughter[3], \
                                      daughter[5], daughter[17], daughter[18], daughter[21], \
                                      daughter[22], daughter[23])
                                superallowed.append([k[0], k[1], k[2], k[3], k[4], k_decay[3], k_decay[4], \
                                      k_decay[5], daughter[0], daughter[1], daughter[3], \
                                      daughter[5], daughter[17], daughter[18], daughter[21], \
                                      daughter[22], daughter[23]])
                        except TypeError:
                            pass

In [ ]:

for l in superallowed: print(l)
print(len(superallowed))

In [ ]:

sa = [i for i in superallowed if i[12] is not None]

In [ ]:

print(len(sa))

In [ ]:

sa_array = np.array(sa)
A = sa_array[:,2].astype(int)
logft = sa_array[:,12].astype(float)
d_logft = sa_array[:,13].astype(float)

In [ ]:

%matplotlib notebook
fig, ax = plt.subplots(figsize=(9,6))

ax.errorbar(A, logft, yerr=d_logft, color='k', fmt='o', capsize=5, label=r'Superallowed $0^{+} \rightarrow 0^{+}$')
#ax.scatter(A, logft, color='k', label=r'Superallowed $0^{+} \rightarrow 0^{+}$')
ax.legend(loc='best', fontsize=15)
ax.set_xlabel(r'$A$', size=20)
ax.set_ylabel(r'log$(ft)$', size=20)

plt.grid()
plt.tight_layout()
plt.savefig("superallowed.png", dpi=fig.dpi)
plt.show()

In [ ]:

x=e.get_ecbp(edata, "S28", 0, units="best")
for k,v in x.items():
    print(k,v)

Search for first-forbidden transitions in \beta^{-} decay

In [ ]:

# Find 1st forbidden transitions based on a.m. selection rules and compare logft values
forbidden_first_list = []
forbidden_first=e.find_forbidden(edata,'1F',mode='ECBP')
for k,v in forbidden_first.items():
    for vv in v:
        if vv[7] is not None:
            forbidden_first_list.append([k[11], vv[7], vv[8]])

forbidden_firstu_list = []
forbidden_firstu=e.find_forbidden(edata,'1UF',mode='ECBP')
for k,v in forbidden_firstu.items():
    for vv in v:
        if vv[7] is not None:
            forbidden_firstu_list.append([k[11], vv[7], vv[8]])

In [ ]:

first_array = np.array(forbidden_first_list)
f_A = first_array[:,0].astype(int)
f_logft = first_array[:,1].astype(float)
f_dlogft = first_array[:,2].astype(float)

firstu_array = np.array(forbidden_firstu_list)
fu_A = firstu_array[:,0].astype(int)
fu_logft = firstu_array[:,1].astype(float)
fu_dlogft = firstu_array[:,2].astype(float)

In [ ]:

%matplotlib notebook
fig, ax = plt.subplots(figsize=(9,6))

ax.errorbar(f_A, f_logft, yerr=f_dlogft, color='k', fmt='o', capsize=5, alpha=0.5, label=r'First forbidden')
ax.errorbar(fu_A, fu_logft, yerr=fu_dlogft, color='r', fmt='o', capsize=5, alpha=0.5, label=r'First forbidden unique')
ax.legend(loc='best', fontsize=15)
ax.set_xlabel(r'$A$', size=20)
ax.set_ylabel(r'log$(ft)$', size=20)

ax.hlines(y=6, xmin=0, xmax=260, linewidth=2, color='b')
ax.hlines(y=9, xmin=0, xmax=260, linewidth=2, color='b')

plt.grid()
plt.tight_layout()
plt.savefig("forbidden.png", dpi=fig.dpi)
plt.show()

In [ ]: