Feature request: Custom atom form factors

[DanPorter]

A user has requested a way to add custom anomalous form factors for atoms.

I should investigate a way of adding custom atoms with any properties, but this will mean quite a big change.

[DanPorter]

You could probably hack the current code by adding  your pseudo-atom name to ELEMENT_LIST, then adding the parameters to for this to "data/Dans Element Properties.txt" and "data/atomic_scattering_factors.npy".

The following untested code should be roughly right:

"""
Dans_Diffraction example
Alter the element properties files

14/4/2025
"""
import os
import numpy as np
import Dans_Diffraction as dif
from Dans_Diffraction.functions_crystallography import ATOMFILE, datadir


# 1. Add pseudo-atom label to ELEMENT_LIST in functions_crystallography.py:
#    'c', 'i', 'p', 'h', 's', 'n', 'o',
#    'D', 'd',  # add Deuterium
#    'ClO',  # this must match other labels for this atom below.
#]

# 2. Edit Element Properties
print(f"Location of element properties file: {ATOMFILE}")

# all of these parameters must be filled in, but most don't need to mean anything
new_properties = {
    'Z': 27,
    'Element': 'ClO',  # this one is used to select the atom
    'Name': 'myClO4',  # this one will be printed
    'Group': 9,
    'Period': 4,
    'Block': 'd',
    'ValenceE': 2,
    'Config': '1s2.2s2.2p6.3s2.3p6.3d7.4s2',
    'Radii': 135.,
    'Weight': 58.9332,
    'Coh_b': 2.49,  # used for neutron scattering
    'Inc_b': -6.2,
    'Nabs': 37.18,
    'Nscat': 5.6,
    'a1': 12.2841,  # a1-c are used for x-ray scattering, if unsure, put a1==sum(Z) and the rest to zero.
    'b1': 4.2791,
    'a2': 7.3409,
    'b2': 0.2784,
    'a3': 4.0034,
    'b3': 13.5359,
    'a4': 2.3488,
    'b4': 71.1692,
    'c': 1.0118,
    'j0_A': 0.4139,
    'j0_a': 16.1616,
    'j0_B': 0.6013,
    'j0_b': 4.7805,
    'j0_C': -0.1518,
    'j0_c': 0.021,
    'j0_D': 0.1345,
    'K': 7.709,
    'L1': 0.925,
    'L2': 0.793,
    'L3': 0.778,
    'M1': 0.101,
    'M2': 0.059,
    'M3': 0.06,
    'M4': 0.,
    'M5': 0.,
    'N1': 0.,
    'N2': 0.,
    'N3': 0.,
}

# edit atom file
with open(ATOMFILE, 'r') as f:
    atomfile_str = f.read()

# make a backup
backup = ATOMFILE.replace('.txt', '_backup.txt')
if not os.path.isfile(backup):
    with open(backup, 'w') as f:
        f.write(atomfile_str)

# add extra line
atomfile_str += ' '.join(new_properties.values()) + '\n'
with open(ATOMFILE, 'w') as f:
    f.write(atomfile_str)


# 3. add values to atomic scattering factors file
asf_file = os.path.join(datadir, 'atomic_scattering_factors.npy')
asf = np.load(asf_file, allow_pickle=True)
asf = asf.item()

# make a backup
asf_backup = asf_file.replace('.npy', '_backup.npy')
if not os.path.isfile(backup):
    np.save(asf_backup, asf, allow_pickle=True)

new_asf = {
    'Z': 27,
    'energy': [10.0, 10.1, 30000.0],  # energy values in eV
    'f1': [1, 1, 1],  # all lists should be the same length. Typically, about 500 values
    'f2': [0, 0, 0],
}

# add extra line
asf['ClO'] = new_asf  # should have the same name as 'Element' above
np.save(asf, asf, allow_pickle=True)


# 4. Make new crystal
xtl = dif.Crystal()
xtl.Atoms.addatom(0, 0, 0, 'ClO')
xtl.generate_structure()

print(xtl)  # if all goes well - this will work!

[DanPorter]

Added custom form factor functions in current feature branch: adbeb52

Example usage:

import numpy as np
import Dans_Diffraction as dif

xtl = dif.Crystal('file.cif')

# add custom element 
name = 'Co3+'
# custom analytical coefficients for the form factor (a1, b1, a2, b2, a3, b3, a4, b4, c, 0)
coefs = (12.2841,  4.5791,  7.4409,  0.6784,  4.2034, 12.5359,  2.2488, 72.1692,  1.1118,  0.)
# custom dispersion corrections as a table of [energy_kev, f1, f2]
# the dispersion corrections are only used when calculating structure factors with energy defined
table = np.array([
    [7.7, 7.71, 7.72, 7.73, 7.74, 7.75, 7.76, 7.77, 7.78, 7.79],  # keV
    [-9.6866, -11.0362, -10.1969, -9.3576, -8.5184, -7.6791, -6.8398, -6.0005, -5.1613, -4.322],  # f1
    [-0.4833, -3.7933, -3.7862, -3.7791, -3.772, -3.765, -3.7579, -3.7508, -3.7437, -3.7367]  # f2
])
# add the new element to the internal table (this is lost when python is closed)
dif.fc.add_custom_form_factor_coefs(name, *coefs, dispersion_table=table)
# Change one of the atoms in the crystal to the new atom
xtl.Structure.changeatom(3, type=name)

# Calculate custom structure factors
en = np.arange(7.7, 7.8, 0.01)  # this must be within the range of the table
ref = [0, 0, 3]
inten = xtl.Scatter.intensity(ref, 'custom', energy_kev=en)