Experimental measurements are often obtained from a mixture of conformers. To accurately predict the experimental values, a series of calculations over a set of conformers must be performed, and the Boltzmann-weighted average of the quantity of interest obtained. The spectroscopy.py script automates such calculations, which include conformational searching of the input structure (with MacroModel, using the OPLS4 or OPLS_2005 force field), quantum mechanical refinement of the conformers (using Jaguar), elimination of redundant conformers (again using MacroModel), calculation of the property of interest for the remaining conformers (using Jaguar), and computation of Boltzmann weighting factors for the conformers (using Jaguar). This information can then be used to generate Boltzmann-averaged quantities. This script currently supports computation of VCD or ECD spectra, optical rotation calculations, NMR chemical shifts of carbon and hydrogen atoms, and NMR spin-spin couplings between atom pairs. IR spectra can be computed separately or as part of a VCD computation. UV-Vis can be computed separately or as part of a ECD computation. NOE distances can also be computed, either separately or as part of an NMR computation.

The command syntax is as follows:

jaguar run spectroscopy.py ([-ir] [-vcd] [uv_vis] [-ecd] |[-nmr] [-chem_shifts] [-ssc] [-noe] | [-optrot) [options] filename

where filename is a Maestro structure file. If the file contains multiple structures, all of the structures must be conformers of each other and conformer searching is skipped. If the input file contains a single structure, then a conformational search is run on it by default, although this may be suppressed with the appropriate option, -skip_csrch. For a description of the options, see spectroscopy.py Command Help or run the command with the -h option.

The script offers several choices of solvent. Only water and chloroform are available in MacroModel, so the conformational search for the other solvents is performed with either water or chloroform, depending on whether the chosen solvent is polar or nonpolar. The geometry optimization of the conformers is then done with the chosen solvent by Jaguar. When a solvent is requested, the geometry optimizations and spectroscopy property calculations are done using the PCM solvation model. The Boltzmann weights used to average properties are improved by recomputing the electronic energy using the PBF solvation model.

For a more complete conformational search, you can set -conf_search_acc=high.

You can set the temperature at which the Boltzmann weighting is calculated for the generated conformers with the --temperature TEMPERATURE option.

Below are details of important output files and quantities for CD, NMR, and NOE calculations.

VCD/ECD spectra calculations:

Boltzmann-Averaged VCD/ECD .spm and .csv files are generated and found in the jobname folder in the directory the script is executed in. They are named jobname + a descriptive suffix. The contents of the files are described in Vibrational and Electronic Circular Dichroism Spectra.

The computation of IR spectra is included in a VCD calculation (setting -vcd also sets -ir), but you can compute it on its own by setting the -ir property flag. The computation of UV/Visible spectra is included in an ECD calculation (setting -ecd also sets -uv_vis), but you can compute it on its own by setting the -uv_vis property flag.

Note that while the -vcd/-ir options can be used together with -ecd/-uv_vis, none can be used in conjunction with any nmr or noe property flags (-nmr/-chem_shifts/-ssc/-noe).

Optical Rotation calculations:

A jobname_optrot.csv file is generated after the calculation has finished, which lists the name, weight, specific optical rotation, and data type. The data type is "raw" for conformers and "averaged" for the Boltzmann-averaged structure. -optrot cannot be set in conjunction with the -nmr, -vcd, or -ecd options.

Calculate all NMR properties:

Set the -nmr option to calculate NMR chemical shifts, spin-spin couplings, and H-H/H-F vector distances for NOE analysis. This is equivalent to setting all three of the nmr properties: -chem_shifts, -ssc, and -noe. -nmr cannot be set in conjunction with the -vcd or -ecd options. For each NMR spectrum jobname_nmr.csv file generated in the calculation, a jobname_nmr_peaks.csv file reporting peak positions and heights is also generated.

NMR chemical shift calculations:

These calculations are currently supported for hydrogen, boron, carbon, nitrogen, fluorine, and phosphorus atoms in the gas phase or in chloroform.

The Boltzmann-averaged chemical shifts are reported for each atom, as well as further averaged over 3D equivalent atoms. The shifts are written as atom-level Maestro properties in a copy of the minimum-energy conformer, as Boltzmann-Avg NMR Chemical Shift and Boltzmann-Avg NMR Chemical Shift 3D Avg. The Boltzmann-averaged chemical shifts are reported in csv files that start with jobname_avg_.

Note that -chem_shifts cannot be used in conjunction with the -vcd or -ecd options.

NMR spin-spin coupling calculations:

The Boltzmann-averaged isotropic spin-spin couplings are reported for each atom pair, as well as further averaged over 3D equivalent pairs of atoms. By default, pairs whose members are 3D-equivalent to one another are not computed since their couplings have no physical meaning. To compute couplings for all pairs, add the -compute_all_ssc option. The couplings are written as structure-level Maestro properties in a copy of the minimum-energy conformer as, for example, Boltzmann-Avg F2 H2 Total SpinSpin Coupling iso and Boltzmann-Avg F2 H2 Total SpinSpin Coupling iso 3D Avg.

The Boltzmann- and 3D-averaged spin-spin couplings are printed in both the .out file and a csv file named jobname_SpinSpin_Couplings.csv.

Note that -ssc cannot be used in conjunction with the -vcd or -ecd options.

NOE distance calculations:

For Nuclear Overhauser Effect (NOE) analysis, we also calculate the Boltzmann-averaged distances of all H-H and H-F vectors.

H-H and H-F distances are written to a CSV file and as a table in the .out file, for the nuclear Overhauser effect (NOE). The distance values include Boltzmann-averaged, minimum, and maximum values, and values for each conformer.

The Boltzmann-averaged distances are printed in both the .out file and a csv file named jobname_HH_HF_distances.csv.

Note that -noe cannot be used in conjunction with the -vcd or -ecd options. NOE distances will still be calculated if any nmr property flag (-nmr/-chem_shifts/-ssc) is specified. Using the -noe option alone will skip calculation of other nmr properties.

Atropisomerism

Compounds exhibiting two atropisomers around one dihedral can be accounted for with the -check_for_atropisomers option. Setting this option alone will task the workflow with searching all dihedrals in the input compound for atropisomerism. If the structure of the atropisomers is known, searching around a specific dihedral can be set with the -target_dihedral option. The search involves rotating each dihedral of interest and constructing potential energy surfaces (PES) to look for barriers. If high enough barriers around a dihedral are found, atropisomers are recognized. The PES is constructed using QRNN or xTB by default, depending on the elements present in the input compound. One can specify a DFT method/basis to construct the PES with the -pes_method and -pes_basis options. Currently atropisomerism is not compatible with NMR spectra.

The below example shows how a check for atropisomers can be included as part of the command:

jaguar run spectroscopy.py ... -check_for_atropisomers -target_dihedral [1,4,5,8] -pes_method b3lyp-d3 -pes_basis 6-31G* filename

{ 
 s_m_m2io_version
 :::
 2.0.0 
} 

f_m_ct { 
 s_m_title
 i_m_Molecular_charge
 i_m_Spin_multiplicity
 i_m_ct_format
 :::
 Structure1 
  0
  1
  2
 m_atom[11] { 
  # First column is atom index #
  i_m_mmod_type
  r_m_x_coord
  r_m_y_coord
  r_m_z_coord
  i_m_residue_number
  i_m_color
  i_m_atomic_number
  s_m_color_rgb
  s_m_atom_name
  i_m_minimize_atom_index
  :::
  1 3 -1.593359 -1.026400 -2.195457 1 24 6 808080  "C1"  1
  2 3 -0.565628 -0.480762 -1.202775 1 2 6 A0A0A0  "C2"  2
  3 41 -1.560211 -0.493430 -3.146939 1 21 1 FFFFFF  ""  3
  4 41 -2.606114 -0.944213 -1.799475 1 21 1 FFFFFF  ""  4
  5 41 -1.413367 -2.081553 -2.404048 1 21 1 FFFFFF  ""  5
  6 41 -0.588336 -1.025736 -0.256676 1 21 1 FFFFFF  ""  6
  7 26 0.749607 -0.520029 -1.779849 1 43 7 5757FF  "N7"  7
  8 16 -0.772232 0.882337 -0.969039 1 70 8 FF2F2F  "O8"  8
  9 42 -1.525474 1.167314 -1.465735 1 21 1 FFFFFF  ""  9
  10 43 1.354087 0.054715 -1.209072 1 21 1 FFFFFF  ""  10
  11 43 1.121596 -1.457737 -1.735619 1 21 1 FFFFFF  ""  11
  :::
 } 
 m_bond[10] { 
  # First column is bond index #
  i_m_from
  i_m_to
  i_m_order
  :::
  1 1 2 1
  2 1 3 1
  3 1 4 1
  4 1 5 1
  5 2 6 1
  6 2 7 1
  7 2 8 1
  8 7 10 1
  9 7 11 1
  10 8 9 1
  :::
 } 
} 

$SCHRODINGER/jaguar run spectroscopy.py -jobname=wflow_spectroscopy_vcd_0001 -vcd -forcefield OPLS4 -conf_search_method mixed -conf_search_steps 200 -mm_energy_window 5.0 -outconfs 100 -max_atom_dev 0.5 -qm_energy_window 5.0 -solvent None -k=dftname=b3lyp-d3 '-k=basis=lacvp**' -k=iacc=1 wflow-spectroscopy_vcd-0001.mae

{ 
 s_m_m2io_version
 :::
 2.0.0 
} 

f_m_ct { 
 s_m_title
 i_m_Molecular_charge
 i_m_Spin_multiplicity
 i_m_ct_format
 :::
 Structure1
  0
  1 
  2
 m_atom[11] { 
  # First column is atom index #
  i_m_mmod_type
  r_m_x_coord
  r_m_y_coord
  r_m_z_coord
  i_m_residue_number
  i_m_color
  i_m_atomic_number
  s_m_color_rgb
  s_m_atom_name
  i_m_minimize_atom_index
  :::
  1 3 -1.593359 -1.026400 -2.195457 1 24 6 808080  "C1"  1
  2 3 -0.565628 -0.480762 -1.202775 1 2 6 A0A0A0  "C2"  2
  3 41 -1.560211 -0.493430 -3.146939 1 21 1 FFFFFF  ""  3
  4 41 -2.606114 -0.944213 -1.799475 1 21 1 FFFFFF  ""  4
  5 41 -1.413367 -2.081553 -2.404048 1 21 1 FFFFFF  ""  5
  6 41 -0.588336 -1.025736 -0.256676 1 21 1 FFFFFF  ""  6
  7 26 0.749607 -0.520029 -1.779849 1 43 7 5757FF  "N7"  7
  8 16 -0.772232 0.882337 -0.969039 1 70 8 FF2F2F  "O8"  8
  9 42 -1.525474 1.167314 -1.465735 1 21 1 FFFFFF  ""  9
  10 43 1.354087 0.054715 -1.209072 1 21 1 FFFFFF  ""  10
  11 43 1.121596 -1.457737 -1.735619 1 21 1 FFFFFF  ""  11
  :::
 } 
 m_bond[10] { 
  # First column is bond index #
  i_m_from
  i_m_to
  i_m_order
  :::
  1 1 2 1
  2 1 3 1
  3 1 4 1
  4 1 5 1
  5 2 6 1
  6 2 7 1
  7 2 8 1
  8 7 10 1
  9 7 11 1
  10 8 9 1
  :::
 } 
} 

${SCHRODINGER}/jaguar run spectroscopy.py -jobname=wflow_spectroscopy_ecd_0001 -ecd -forcefield OPLS4 -conf_search_method mixed -conf_search_steps 200 -mm_energy_window 5.0 -outconfs 100 -max_atom_dev 0.5 -qm_energy_window 5.0 -solvent None -k=dftname=b3lyp-d3 '-k=basis=lacvp**' -k=iacc=1 wflow-spectroscopy_ecd-0001.mae

{ 
 s_m_m2io_version
 :::
 2.0.0 
} 

f_m_ct { 
 s_m_title
 i_m_Molecular_charge
 i_m_Spin_multiplicity
 i_m_ct_format
 :::
 Structure1 
  0
  1 
  2
 m_atom[11] { 
  # First column is atom index #
  i_m_mmod_type
  r_m_x_coord
  r_m_y_coord
  r_m_z_coord
  i_m_residue_number
  i_m_color
  i_m_atomic_number
  s_m_color_rgb
  s_m_atom_name
  i_m_minimize_atom_index
  :::
  1 3 -1.593359 -1.026400 -2.195457 1 24 6 808080  "C1"  1
  2 3 -0.565628 -0.480762 -1.202775 1 2 6 A0A0A0  "C2"  2
  3 41 -1.560211 -0.493430 -3.146939 1 21 1 FFFFFF  ""  3
  4 41 -2.606114 -0.944213 -1.799475 1 21 1 FFFFFF  ""  4
  5 41 -1.413367 -2.081553 -2.404048 1 21 1 FFFFFF  ""  5
  6 41 -0.588336 -1.025736 -0.256676 1 21 1 FFFFFF  ""  6
  7 26 0.749607 -0.520029 -1.779849 1 43 7 5757FF  "N7"  7
  8 16 -0.772232 0.882337 -0.969039 1 70 8 FF2F2F  "O8"  8
  9 42 -1.525474 1.167314 -1.465735 1 21 1 FFFFFF  ""  9
  10 43 1.354087 0.054715 -1.209072 1 21 1 FFFFFF  ""  10
  11 43 1.121596 -1.457737 -1.735619 1 21 1 FFFFFF  ""  11
  :::
 } 
 m_bond[10] { 
  # First column is bond index #
  i_m_from
  i_m_to
  i_m_order
  :::
  1 1 2 1
  2 1 3 1
  3 1 4 1
  4 1 5 1
  5 2 6 1
  6 2 7 1
  7 2 8 1
  8 7 10 1
  9 7 11 1
  10 8 9 1
  :::
 } 
} 

${SCHRODINGER}/jaguar run spectroscopy.py -jobname=wflow_spectroscopy_nmr_0001 -nmr -nmr_types 1H 13C -only_average_nmr -forcefield OPLS4 -conf_search_method mixed -conf_search_steps 200 -mm_energy_window 5.0 -outconfs 100 -max_atom_dev 0.5 -qm_energy_window 5.0 -solvent None -k=dftname=b3lyp-d3 '-k=basis=lacvp**' -k=iacc=1 wflow-spectroscopy_nmr-0001.mae