eproducer package
Tools for EPR g-value calculations and spectral analysis.
Subpackages
Submodules
eproducer.g_valsetup
Molecular structure utilities for EPR g-value calculations: connectivity analysis, radical site location, and structure extraction from npz datasets.
- eproducer.g_valsetup.HFC_relevantA(atoms, radicals, depth=2)[source]
Extract the H, O and N atoms in the neighborhood of the radicals up to certain depth (neighbors of neighbors, moving only through heavy atoms). It also includes the H atoms bonded to those O and N atoms.
- Parameters:
atoms (ase.Atoms) – Atoms object containing the molecule.
radicals (list) – posible positions of the radicals, implying that that neighborhood is important for the HFC. Indexes starting from 1.
depth (int. Default=2) – number of times that it searches the neighbors of the neighbors.
- Returns:
(dict) {element (list of indexes} with the indexes of the H, O and N)
atoms grouped by element, e.g. {‘H’ ([3, 4], ‘O’: [7]}. The indices start)
with 1.
- eproducer.g_valsetup.ase_to_nx_graph(atoms, engine='vmd')[source]
Converts a molecule into a networkx graph whose nodes are the atoms, with their atomic number in the attribute ‘atomic_number’, and whose edges are the bonds.
- Parameters:
atoms (ase.Atoms) – molecule to be converted.
engine (str. Default='vmd') – engine used to get the connectivity, ‘vmd’ or ‘ase’. See get_connectivity.
- Return type:
(networkx.Graph) graph of the molecule.
- eproducer.g_valsetup.ext_xyz_from_npz(npz_file, selection=None, create_xyz_files=True)[source]
Extracts structures in xyz files from a given npz file.
- Parameters:
npz_file (str) – file containing the structures.
selection (list. Default=None) – list of indexes of the structures to extract. If None, all structures are extracted.
create_xyz_files (bool. Default=True) – whether to create xyz files for the extracted structures or not.
- Return type:
(list) list of ASE Atoms objects corresponding to the extracted structures.
- eproducer.g_valsetup.get_connectivity(atoms, engine)[source]
Computes the connectivity of a molecule.
- Parameters:
atoms (ase.Atoms) – molecule to be analyzed.
engine (str) – ‘vmd’ to use the bonds detected by VMD (VMD must be in the PATH) or ‘ase’ to use the natural cutoffs of ase.neighborlist.
- Returns:
(list) [n_atoms] list where the element i contains the indices of the
atoms bonded to the atom i. Indices start from 0.
- eproducer.g_valsetup.get_isomorphic_permutation(nx_graph1, nx_graph2)[source]
Returns a permutation list to reorder the atoms of nx_graph2 to match nx_graph1 based on their isomorphic structure (nodes are matched by their ‘atomic_number’), such that e.g.:
permutation = get_isomorphic_permutation(g1, g2) numbers2[permutation] == numbers1
- Parameters:
nx_graph1 (networkx.Graph) – first molecular graph (see ase_to_nx_graph).
nx_graph2 (networkx.Graph) – second molecular graph (see ase_to_nx_graph).
- Returns:
(list) permutation such that the node i of nx_graph1 corresponds to the
node permutation[i] of nx_graph2. Raises ValueError if the graphs are not
isomorphic.
- eproducer.g_valsetup.get_permutation_of_heavy_atoms(heavy_ref_mol, heavy_mch_mol, engine='vmd')[source]
Finds the permutation that maps the heavy atoms of the reference molecule onto the heavy atoms of the molecule to match, using the isomorphism of their molecular graphs.
- Parameters:
heavy_ref_mol (ase.Atoms) – reference molecule (only heavy atoms).
heavy_mch_mol (ase.Atoms) – molecule to match (only heavy atoms).
engine (str. Default='vmd') – engine used to get the connectivity, ‘vmd’ or ‘ase’. See get_connectivity.
- Returns:
(list) permutation such that the atom i of heavy_ref_mol corresponds to
the atom permutation[i] of heavy_mch_mol. Raises ValueError if the
molecular graphs are not isomorphic.
- eproducer.g_valsetup.heavya_idx_from_ref(ref_mol_xyz: str, rad_mol_xyz: str, index: int, engine: str = 'vmd') ndarray[source]
Find the index, in the radical, of the heavy atom that corresponds to a given heavy atom of the reference molecule. The heavy atoms of both molecules are matched by graph isomorphism.
- Parameters:
ref_mol_xyz (str) – path to the file describing the molecule with the hydrogen included.
rad_mol_xyz (str) – path to the file describing the molecule without the hydrogen.
index (int) – index of the heavy atom in the reference molecule to match in the new molecule. Starting from 1.
engine (str. Default='vmd') – engine to use to get the connectivity. It can be ‘vmd’ or ‘ase’.
- Returns:
(int) index of the corresponding heavy atom in rad_mol_xyz. Index
counting starting from 1.
- eproducer.g_valsetup.iHFC_fromxyz(file, radicals, depth)[source]
Extract the H, O and N atoms in the neighborhood of the radicals up to certain depth (see HFC_relevantA). Those atoms are the relevant ones for HyperFine Calculations (HFC). This function receives strings as inputs.
- Parameters:
file (str) – path to the molecule that you want to analyse. ASE must be able to read this file.
radicals (str) – index of the atoms of posible location of the radicals, implying that that neighborhoods of this atoms are important for the HFC. 1-based indexes. It should be a list of integers in string format, e.g.: “[1, 5, 7]”.
depth (str) – number of times that it searches the neighbors of the neighbors.
- Returns:
(dict) {element (list of indexes} with the indexes of the H, O and N)
atoms grouped by element. The indices start from 1.
- eproducer.g_valsetup.nrand_rad(npz_file, n=None, Nmin=None, entry='energy_MACE', exclude=None)[source]
Extracts a random subset of radical structures from a given npz file. It takes the Nmin structures with the lowest values of entry and randomly selects n of them.
- Parameters:
npz_file (str) – file containing the radicals.
n (int. Default=None) – number if radical structures to extract from the subset. If None, all the structures in the subset are selected.
Nmin (int. Default=None) – number of structures in the subset that minimizes the entry parameter. If None, all the structures are included in the subset.
entry (str. Default='energy_MACE') – entry of the npz file from where the Nmin are selected.
exclude (list. Default=None) – indexes of the structures that cannot be selected, e.g. because they were already computed.
- Return type:
(list) indexes of the selected entries in the npz file.
- eproducer.g_valsetup.rad_loc(ref_mol_xyz: str, rad_mol_xyz: str, engine='vmd') int[source]
Find the index of the atom where a hydrogen atom was attached. The heavy atoms of both molecules are matched by graph isomorphism and the first heavy atom with a different number of bonds is returned.
- Parameters:
- Returns:
(int) index, in rad_mol_xyz, of the atom where the hydrogen was attached.
Index counting starting from 1. None if no atom lost a bond.
- eproducer.g_valsetup.rad_loc_deprecated(ref_mol: str, radical: str) ndarray[source]
Find the index of the atom where a hydrogen atom was attached.
Parametes
- ref_mol: str
path to the file describing the molecule with the hydrogen included.
- radical: str
path to the file describing the molecule without the hydrogen.
- returns:
(int) index of the atom where the hydrogen was attached. Index counting
starting from 1.
eproducer.postpro
Post-processing of ORCA EPR outputs: g-value extraction, tabulation, spectral comparison, and bond dissociation enthalpy.
- class eproducer.postpro.GvalComparison(molecules)[source]
Bases:
objectTable (pandas.DataFrame stored in GvalComparison.df) to compare the g-values and spectra of several candidate molecules with each other and with experiments. Each method adds one or more columns to the table.
- Parameters:
molecules (list) – paths to the directories of the candidate molecules. They are stored in the column ‘molecule’.
- add_gvals(output_files)[source]
Extracts the g-values of each molecule from its orca output and stores them in the column ‘gval’.
- add_peak_loc()[source]
Stores the field where each spectrum has its maximum in the column ‘peak_loc’. Requires the column ‘spectrum’ (see GvalComparison.add_spectra).
- Return type:
(list) field of the maximum of each spectrum.
- add_spectra(spect_files, experiment=None)[source]
Loads the computed spectrum of each molecule, interpolated onto the field of the experiment (see TheoMatchExpe), and stores it in the column ‘spectrum’. It also adds the column ‘peak_loc’.
- Parameters:
spect_files (str or list) – .dat files (field vs intensity) of the computed spectra, one per molecule. If a str is given, it is taken as the name of the file inside the directory of each molecule, i.e. <molecule>/<spect_files>.
experiment (str. Default=None) – .dat file of the experimental spectrum whose field is used for the interpolation. If None, the field of the first spectrum is used.
- Return type:
(numpy.array) [n_molecules x n_points] interpolated spectra.
- d_target_gval(gval_ref, column_name)[source]
Computes the euclidean distance between the g-values of each molecule and a reference, and stores it in a new column. Requires the column ‘gval’ (see GvalComparison.add_gvals).
- fit_to(experiment_file, **kwargs)[source]
Fits the spectra of the molecules to an experimental spectrum, keeping only the molecules that contribute the most (see TheoMatchExpe.gradual_cleaning). Requires GvalComparison.add_spectra to be called first.
- individual_corr(experiment, column_name)[source]
Computes the Pearson correlation between each spectrum and an experimental spectrum, and stores it in a new column. Requires the column ‘spectrum’ (see GvalComparison.add_spectra).
- lc_max_corr(experiment, column_name)[source]
Finds the coefficients of the linear combination of the spectra that maximizes the Pearson correlation with an experimental spectrum, and stores them in a new column. Requires the column ‘spectrum’ (see GvalComparison.add_spectra).
- Parameters:
experiment (numpy.array) – experimental intensities, defined on the same field values as the spectra in the column ‘spectrum’.
column_name (str) – name of the new column with the coefficients.
- Return type:
(tuple) coefficients and the maximum correlation reached.
Note
The optimization starts from random coefficients, so the result can change between runs.
- peak_loc(x, y)[source]
Finds the maximum of a spectrum.
- Parameters:
x (numpy.array) – field values.
y (numpy.array) – intensities.
- Return type:
(tuple) field where the maximum is located and the maximum intensity.
- eproducer.postpro.add_gval2npz(path)[source]
Add g-values to the npz files in the entry ‘g-values’, based on the files with the shape <name>_<index>_epr.out next to <name>.npz. The g-values of each file are stored in the row <index> of ‘g-values’, which is created with zeros (one row per radical in ‘heavy_atom_missing_idxs’) if it does not exist.
- Parameters:
path (str) – path to the directory containing the npz files and the epr.out files, the npz file itself or a .dat file with a list of npz files.
- eproducer.postpro.compute_bde(reactants=None, products=None)[source]
Compute the bond dissociation enthalpy (BDE) for a reaction as the difference between the total enthalpy of products (P) and reactants (R), such that BDE = H_P - H_R.
- Parameters:
reactants (iterable) – An iterable (e.g., list or tuple) of reactant species. Each element should be an orca output file path with computed enthalpy.
products (iterable) – An iterable (e.g., list or tuple) of product species. Each element should be an orca output file path with computed enthalpy.
- Return type:
(float) The BDE computed. The units are the same as in the orca outputs.
- eproducer.postpro.create_table_per_method(data: ndarray, experiment: list, methods: list) str[source]
For a given system, creates a table with the different methods used to compute the g-values.
- Parameters:
data (numpy.array) – g-values and errors in the shape [meth][gx, gy, gz, max absolute error, mean absolute error], as returned by several_gvals.
experiment (numpy.array) – expected values obtained from the exteriment. [gx, gy, gz] in decreassing order.
methods (list) – list of names of the methods used. they has to be in the same order than data.
- Returns:
(str) Table in LaTeX format of g-values with the different methods and
the associated max error and average error. The table is also printed.
- eproducer.postpro.create_table_per_system(experiment: list, systems: list, orca_file: str = 'epr_info.out') str[source]
Creates a markdown table with the g-values of several systems computed with the same method, and their errors with respect to the experiment.
- Parameters:
experiment (list) – expected values obtained from the exteriment. [gx, gy, gz] in decreassing order.
systems (list) – paths to the directories of the systems, ending with ‘/’. Each directory contains orca_file, and the name of the system in the table links to <system>vmd_image.md.
orca_file (str. Default='epr_info.out') – name of the orca file containing the g-values.
- Returns:
(str) Table of g-values with the different systems and the associated
max error and average error.
- eproducer.postpro.extract_gvals(orca_out: str) ndarray[source]
This function extracts the g-values calculated by orca.
- Parameters:
orca_out (str) – path to the output file from orca.
- Return type:
(np.array) [gx, gy, gz] in decreasing order.
- eproducer.postpro.extract_system_info(sys_path: str, exp_values: str, orca_output: str = 'epr_info.out')[source]
Creates gvalues_table.md and gvalues.png in sys_path.
- Parameters:
sys_path (str) – path to the folder containing the candidates. The candidates are the subdirectories whose name contains ‘-‘.
exp_values (str) – String of the list of expected experimental values in shape of python list, f.e. ‘[2.0062, 2.0055, 2.0022]’. Use ‘’ if there are no experimental values.
orca_output (str. Default='epr_info.out') – name of the orca output file with the g-values in each candidate directory.
- Return type:
(StandardPlotter) StandardPlotter used to plot the gvalues.
- eproducer.postpro.extract_values(systems: list, orca_files: list, experiment: list)[source]
This function extracts the information of the gvalues.
- Parameters:
systems (list) – paths to the directories containing the orca files. Each path should correspond to a different system.
orca_files (str) – names of the orca files containing the g-values separated by commas, e.g. ‘cc-pVDZ.out,epr_info.out’. Each one should correspond to the calculation of the gvalues with different methods.
experiment (list) – expected values obtained from the experiment. [gx, gy, gz] in decreasing order.
- Returns:
(dict) {<system_path> (list with one array per orca file [gx, gy, gz,)
max absolute error, mean absolute error]}, see several_gvals.
- eproducer.postpro.several_gvals(names: str, experiment: ndarray, path: str = './')[source]
Extract several gvals and compare them with experimental values.
- Parameters:
names (str) – names of the orca output files separated by commas. Eg: ‘cc-pVDZ.out,epr_info.out’
experiment (numpy.array) – expected values according to the experiment [gx, gy, gz] in decreasing order.
path (str. Default='./') – path to the directory containing the output files from orca.
- Returns:
(list) one array per orca file with [gx, gy, gz, max absolute error,
mean absolute error] as strings.
eproducer.best_fit
Spectral fitting of computed EPR spectra against experimental data using non-negative least squares.
- class eproducer.best_fit.TheoMatchExpe(files, experiment_file=None, **kwargs)[source]
Bases:
objectCollects the computed data and tries to fit it to the experiment spectrum
- Parameters:
files (list) – paths to the .dat files (field vs intensity) of the computed spectra of the candidates.
experiment_file (str. Default=None) – .dat field vs intensity of the experiment. If None, the first file of files is used as experiment.
**kwargs – passed to TheoMatchExpe.fit_shifting, e.g. fieldrange or shifting_width.
Note
The computed spectra are interpolated onto the field values of the experiment, and the experimental intensities are normalized to a maximum of 1.
- clean_candidates(threshold=1, **kwargs)[source]
Removes the candidates that are expected to be less than certain percentage. It means, TheoMatchExpe.intensities and TheoMatchExpe.files might be reduced. The fit is repeated until all the remaining candidates reach the threshold.
- Parameters:
threshold (float. Default=1) – minimum percentage that a candidate should have to be considered.
**kwargs – passed to TheoMatchExpe.fit_shifting, e.g. shifting_width.
- Return type:
(np.array) new set of intensities.
- fit_shifting(files=None, shifting_width=None, shift_step=0.2, **kwargs)[source]
Shifts the experimental field within [-shifting_width, shifting_width] and keeps the shift whose fit (see TheoMatchExpe.fit_theoretical) has the smallest residual. The results are stored in the attributes coeffs, intensities, intensfit, error, corr and shift.
- Parameters:
files (list. Default=None) – paths to the .dat files of the computed spectra. If None, TheoMatchExpe.files is used.
shifting_width (float. Default=None) – maximum shift of the field in both directions. If None, the last value given is used (0 at the beginning, i.e. no shifting).
shift_step (float. Default=0.2) – step between the tested shifts.
**kwargs – fieldrange: list. [min, max] field range to fit. If not given, the last value given is used. None (the initial value) means the whole range of the experiment.
- Returns:
(tuple) coeffs, intensities, intensfit, error and corr of the best shift.
See TheoMatchExpe.fit_theoretical.
- fit_theoretical(files=None, fieldexp=None, fieldrange=None)[source]
Fits the experimental spectrum with a non-negative linear combination of the computed spectra (non-negative least squares).
- Parameters:
files (list. Default=None) – paths to the .dat files of the computed spectra. If None, TheoMatchExpe.files is used.
fieldexp (numpy.array. Default=None) – field values where the computed spectra are interpolated. If None, TheoMatchExpe.fieldexp is used.
fieldrange (list. Default=None) – [min, max] field range to fit. The experimental intensities outside this range (shifted by TheoMatchExpe.shift) are set to zero before fitting. If None, the whole range is used.
- Returns:
(tuple) coeffs [n_files x 1] of the linear combination, intensities
[n_files x n_points] interpolated computed spectra, intensfit [n_points]
fitted spectrum, error (residual of the fit) and corr (Pearson
correlation between the fitted spectrum and the experiment).
- gradual_cleaning(threshold=1, steps=0.1, **kwargs)[source]
Removes gradually the candidates that do not weight much in the fitting up to having all the candidates with a minimum percentage in the fitting. It starts again from all the candidates given initially.
- Parameters:
threshold (float. Default=1) – minimum percentage per candidate to be considered into the fitting.
steps (float. Default=0.1) – size of the steps in the filtering.
**kwargs – passed to TheoMatchExpe.fit_shifting, e.g. shifting_width.
- Return type:
(np.array) new set of intensities.
- eproducer.best_fit.best_fit(experiment, output, fieldrange, *argv)[source]
Fits an experimental spectrum with the computed spectra of several candidates, keeps the candidates that contribute at least 5 % to the fit, prints their percentages and saves a plot of the fit.
- Parameters:
experiment (str) – path to the .dat file containing the experimental field and the intensity.
output (str) – name of the png file where the plot is stored.
fieldrange (str) – [min, max] field range to fit in python list format, e.g. ‘[640, 650]’. Use ‘’ to fit the whole range. The range is marked with dotted lines in the plot.
*argv (str) – paths to the .dat files containing the computed spectra of the candidates.
- Returns:
(None) it prints the percentage of each candidate and the shift of the
field, and saves the plot in output.
Note
The experimental field is shifted up to 10 units in both directions (in steps of 0.2) to find the best fit. See TheoMatchExpe.fit_shifting.
- eproducer.best_fit.create_fit_file(output, experiment_file='../experiments/field_vs_spectrum.dat', basis_pattern='spectrum_wo_hyFiCorr.dat', basis_location='../')[source]
Creates a .dat file with field vs intensity of the best fit to an experiment with some given calculated spectrums.
- Parameters:
output (str) – name of the .dat file where you want to store the output.
experiment_file (str. Default='../experiments/field_vs_spectrum.dat') – .dat field vs intensity of the experiment.
basis_pattern (str. Default='spectrum_wo_hyFiCorr.dat') – pattern of the .dat files containing the spectrum of the computed spectrum.
basis_location (str. Default='../') – location where to search for the files with the basis_pattern.
- Return type:
(None) it does not return anything but creates the .dat file with the fit.
Note
The basis files are searched recursively in basis_location. They are interpolated onto the field of the experiment (see TheoMatchExpe).
- eproducer.best_fit.fit_experiment(experiment_file='../experiments/field_vs_spectrum.dat', basis_pattern='spectrum_wo_hyFiCorr.dat')[source]
Collects the computed data and tries to fit it to the experiment spectrum
- Parameters:
experiment_file (str. Default='../experiments/field_vs_spectrum.dat') – .dat field vs intensity of the experiment.
basis_pattern (str. Default='spectrum_wo_hyFiCorr.dat') – pattern of the .dat files containing the spectrum of the computed spectrum.
- Returns:
(tuple) fieldexp and intensexp of the experiment, field of the computed
spectra and best_fit_intensity [n_files x n_points], i.e. each computed
spectrum multiplied by its coefficient in the fit.
Note
This script assumes that the basis files has a common pattern and they are all in subdirectories of ‘../’. It also assumes that all the values of the computed intensities are defined for the same values of the field.
- eproducer.best_fit.mat2dat(file_field: str, file_spectrum: str, output: str)[source]
Converts an experimental spectrum stored in MATLAB .mat files into a .dat file with two columns: field and spectrum.
- Parameters:
- Return type:
(None) it creates the .dat file.
eproducer.bonds2npz
Utilities for adding bond connectivity and molecular ID information to npz dataset files.
- eproducer.bonds2npz.addfile2npz(npz_file_ref, npz_file_add, file_name)[source]
Copies one entry of an npz file into another npz file. If the entry already exists in the second file, it is overwritten.
- eproducer.bonds2npz.bonds2npz(npz_file)[source]
Add bond information to npz files in the entry ‘bonds’. The connectivity of the original molecule (‘original_atomic_numbers’ and ‘original_xyz’) is obtained with VMD. Then, for each hydrogen in ‘original_hydrogen_idxs’, the bonds of the molecule without that hydrogen are stored (with the indices shifted accordingly). The npz files that fail are listed in failed.dat.
- Parameters:
npz_file (str) – Path to the npz file or directory containing npz files. It could also be a path to a .dat file that contains a list of npz files.
- Returns:
(list) [n_hydrogens][n_bonds x 2] bonds of each radical of the last npz
file processed.
- eproducer.bonds2npz.molid2npz(npz_file)[source]
Add molecular ID information to npz files in the entry ‘molid’. The ID of the structure i is ‘<npz path without extension>_<i>’, with i written with three digits, e.g. ‘dir/mol_007’.
- Parameters:
npz_file (str) – Path to the npz file or directory containing npz files.
- Return type:
(None) it modifies the npz files.
eproducer.miscellaneous
General-purpose utilities: shell command execution and introspection helpers.
- eproducer.miscellaneous.args_and_defaults(module, *args)[source]
Takes a function and returns its parameters with their default values.
- Parameters:
- Return type:
Notes
moved to pkgdeveloper
- eproducer.miscellaneous.function_doc(module, *args)[source]
Takes a function and returns its documentation.
- Parameters:
- Return type:
(str) function documentation
Notes
moved to pkgdeveloper
- eproducer.miscellaneous.optimized_e(file)[source]
This code finds the last energy in a log file of gaussian computed using RBMK functional. The output is given in eV.
- eproducer.miscellaneous.output_terminal(cmd, print_output=True, skip_error=False, print_cmd=False, **kwargs)[source]
Runs a command in a terminal and returns its output.
- Parameters:
cmd (str) – bash command to be executed in the terminal.
print_output (bool. Default=True) – True for printing the output (without empty lines) besides of returning it.
skip_error (bool. Default=False) – True for continuing running although the command fails. Otherwise, an AssertionError with the error message is raised.
print_cmd (bool. Default=False) – print the command before executing it.
**kwargs – additional options for subprocess.Popen
- Return type:
(str) standard output of the executed command.
- eproducer.miscellaneous.time_g09(logfile)[source]
Function that extracts the time spend for one gaussian simulation from a .log file.
- Parameters:
logfile (string) – .log file obtained during a gaussian simulation.
- Returns:
(float) Time in minutes. The time in seconds, minutes and hours is also
printed. If the simulation starts and ends in different months, it only
prints a message and returns None.
Shell scripts
The shell scripts run the calculations and the post-processing from the
terminal. Get the path of a script with $(eproducer <script> -path) to
run it, copy it or submit it with sbatch. The scripts with flags print
the help messages below when they are called with -h. Click
[script] to see the source code of a script.
g-value calculations
gval_workflow
[script] eproducer/gval_workflow.sh
gval_dataset
[script] eproducer/gval_dataset.sh
gval_data_restart
[script] eproducer/gval_data_restart.sh
Spectra
extract_EPRspec
[script] eproducer/extract_EPRspec.sh
spectrum_from_gval
[script] eproducer/spectrum_from_gval.sh
Post-processing
gval_postpro
[script] eproducer/gval_postpro.sh
gvals_template
[script] eproducer/gvals_template.sh
create_mol_png
[script] eproducer/create_mol_png.tcl
Utilities
basics
[script] eproducer/basics.sh
extract_bonds_vmd
[script] eproducer/extract_bonds_vmd.sh