Skip to content

Commit

Permalink
NMR shift calculation and experimental NMR description reading works
Browse files Browse the repository at this point in the history
  • Loading branch information
Kristaps Ermanis committed Jun 20, 2019
1 parent 05b84b3 commit c744c93
Show file tree
Hide file tree
Showing 2 changed files with 31 additions and 11 deletions.
29 changes: 19 additions & 10 deletions NMR.py
Original file line number Diff line number Diff line change
Expand Up @@ -73,8 +73,8 @@ def ExpNMRFromDesc(self):

ExpNMR_file.close()

self.Clabels, self.Cshifts = self.ParseExp(self, Cexp)
self.Hlabels, self.Hshifts = self.ParseExp(self, Hexp)
self.Clabels, self.Cshifts = self.ParseExp(Cexp)
self.Hlabels, self.Hshifts = self.ParseExp(Hexp)
self.Equivalents = equivalents
self.Omits = omits

Expand All @@ -95,18 +95,19 @@ def ParseExp(self, exp):
return expLabels, expShifts


def CalcBoltzmannWeightedShifts(Isomers):
def CalcBoltzmannWeightedShieldings(Isomers):

energies = []

for i, iso in enumerate(Isomers):

# Calculate rel. energies in kJ/mol
minE = min(iso.DFTEnergies)

relEs = []

for e in energies:
relEs.append(((e - minE) * hartreeEnergy))
for e in iso.DFTEnergies:
relEs.append((e - minE) * hartreeEnergy)

Isomers[i].Energies = relEs

Expand All @@ -118,16 +119,18 @@ def CalcBoltzmannWeightedShifts(Isomers):

q = sum(populations)

for i in range(0, len(populations)):
populations[i] = populations[i] / q
for p in range(0, len(populations)):
populations[p] = populations[p] / q

Isomers[i].Populations = populations

print(Isomers[i].Populations)

# Calculate Boltzmann weighed shielding constants
# by summing the shifts multiplied by the isomers population
BoltzmannShieldings = []

for atom in range(iso.Atoms):
for atom in range(len(iso.Atoms)):
shielding = 0
for population, shieldings in zip(iso.Populations, iso.ConformerShieldings):
shielding = shielding + shieldings[atom] * population
Expand Down Expand Up @@ -181,12 +184,12 @@ def NMRDataValid(Isomers):
return True


def GetCalcShifts(Isomers, settings):
def CalcNMRShifts(Isomers, settings):

print('WARNING: NMR shift calculation currently ignores the instruction to exclude atoms from analysis')
for i, iso in enumerate(Isomers):

BShieldings = iso.IsomerShieldings
BShieldings = iso.BoltzmannShieldings
Cvalues = []
Hvalues = []

Expand All @@ -203,6 +206,12 @@ def GetCalcShifts(Isomers, settings):
Isomers[i].Cshifts = Cvalues
Isomers[i].Hshifts = Hvalues

print('C shifts for isomer ' + str(i) + ": ")
print(', '.join(['{0:.3f}'.format(x) for x in Isomers[i].Cshifts]))

print('H shifts for isomer ' + str(i) + ": ")
print(', '.join(['{0:.3f}'.format(x) for x in Isomers[i].Hshifts]))

return Isomers


Expand Down
13 changes: 12 additions & 1 deletion PyDP4.py
Original file line number Diff line number Diff line change
Expand Up @@ -158,7 +158,7 @@ def __init__(self, InputFile, Charge=-100):
self.NMROutputFiles = [] # list of DFT NMR output file names
self.ShieldingLabels = [] # A list of atom labels corresponding to the shielding values
self.ConformerShieldings = [] # list of calculated NMR shielding constant lists for every conformer
self.IsomerShieldings = [] # Boltzmann weighted NMR shielding constant list for the isomer
self.BolztmannShieldings = [] # Boltzmann weighted NMR shielding constant list for the isomer
self.Cshifts = [] # Calculated C NMR shifts
self.Hshifts = [] # Calculated H NMR shifts

Expand Down Expand Up @@ -324,6 +324,17 @@ def main(settings):
print('\nReading experimental NMR data...')
NMRData = NMR.NMRData(settings.NMRsource)

print('Clabels: ' + str(NMRData.Clabels))
print('Cshifts: ' + str(NMRData.Cshifts))

print('Hlabels: ' + str(NMRData.Hlabels))
print('Hshifts: ' + str(NMRData.Hshifts))

print('Equivalents: ' + str(NMRData.Equivalents))
print('Omits: ' + str(NMRData.Omits))

quit()

print("Conformation data:")
NMR.PrintConformationData(Isomers)

Expand Down

0 comments on commit c744c93

Please sign in to comment.