#!/usr/bin/python
# -*- coding: utf-8 -*-
#
#--------------------------------------------------------------------------------
#gauss2mol2.py v1.1, Copyright Bjoern Olausson
#--------------------------------------------------------------------------------
#This program is free software; you can redistribute it and/or modify
#it under the terms of the GNU General Public License as published by
#the Free Software Foundation; either version 2 of the License, or
#(at your option) any later version.
#
#This program is distributed in the hope that it will be useful,
#but WITHOUT ANY WARRANTY; without even the implied warranty of
#MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
#GNU General Public License for more details.
#
#To view the license visit
#http://www.gnu.org/licenses/old-licenses/gpl-2.0.html
#or write to
#Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
#--------------------------------------------------------------------------------
#--------------------------------------------------------------------------------
#
#This script uses antechamber to read gaussion log files, and generates mol2 files with resp charges
#Antechamber does not preserve the coordinates an correct atomtypes.
#So in this script you can choose to replace the charges in your original mol2 file with the gaussian charges
#(taken from the antechamber generated mol2 file) or to replace the coordinates and atom types in the antechamber
#generated mol2 file with the coordinates and atom types from your original mol2 file.
#
#run with --help to see how to use this script.
#
import sys, os, shutil, string, glob, subprocess, itertools
from optparse import OptionParser
parser = OptionParser()
parser.add_option(
"-f",
action="store",
type="string",
dest="F",
help="""File(s) to process. Original mol2 file(s) must be in the same place as the gaussian log file(s) and must have the same prefix (AAA.mol2, AAA.log). If you want to use wildcards, please wrap your expression in quotes (gauss2mol2.py --crd -f "*.log*" )""",
metavar="FILES"
)
parser.add_option(
"-c",
"--chg",
action="store_true",
dest="CHG",
help="Replace the charges in your mol2 with the charges from the antechamber generated mol2 file",
metavar="REPCHG"
)
parser.add_option(
"-o",
"--crd",
action="store_true",
dest="CRD",
help="Replace the coordinates and atom types in the antechamber generated mol2 file with the coordinates and atom types of your original mol2 file",
metavar="REPCRD"
)
# instruct optparse to parse the program's command line
(options, args) = parser.parse_args()
# Checking for required otions
if not options.F:
parser.error("-f must be set")
elif not options.CHG and not options.CRD:
parser.error("--chg or --crd must be set")
#Convert Options to a list
F = options.F.split(" ")
CHG = options.CHG
CRD = options.CRD
#Change to the current working dir
path = os.getcwd()
os.chdir(path)
def linefilter(lines, start, stop) :
lines = iter(lines) #state-preserving iterator
#if 'start' matches BOF, we will start to read at the beginning of the file,
#else we start at the following line that matches 'start'
if start != "BOF" :
for line in lines :
if line == start :
break
for line in lines :
#if 'stop' matches EOF, we will read to the end of file.
#else we stop reading on the line that matches 'stop'
if stop == "EOF" :
yield line
else :
if line == stop :
break
else :
yield line
#If "*" is found in a parameter, glob it.
for element in F :
if '*' in element :
F.extend(glob.glob(element))
F.remove(element)
def altercrd(gaufile, orgfile) :
#open the new file
fout = open(IN+"-CRD.mol2", "w")
#write all from beginning to '@<TRIPOS>ATOM' from the gaussian mol2 to the new file
top = linefilter(gaufile.splitlines(), 'BOF', '@<TRIPOS>ATOM')
for line in linefilter(gaufile.splitlines(), 'BOF', '@<TRIPOS>ATOM') :
fout.writelines(line+"
")
#Write the '@<TRIPOS>ATOM' match string cause i is stripped in the linefilter
fout.writelines('@<TRIPOS>ATOM'+"
")
#Now lets write the collumns 0,1 from the gaussian mol2file, the collumn 2,3,4 and 5
#from the original mol2 files and again the collumns 6,7 and 8 from the gaussian mol2 file
#So it should look like this
#
#ORIGINAL FILE EXSAMPLE:
#@<TRIPOS>ATOM
# 1 N1 -18.5878 -6.4372 8.0759 N.pl3 1 <1> -0.9900
# 2 H1 -18.4216 -6.0346 8.9973 H 1 <1> 0.3600
#@<TRIPOS>BOND
#GAUSSIAN FILE EXAMPLE:
#@<TRIPOS>ATOM
# 1 N1 5.1000 0.1990 -0.0330 N.3 1 MOL -1.0597
# 2 H1 5.9140 -0.3800 -0.1080 H 1 MOL 0.3916
#@<TRIPOS>BOND
#Merge table
#G == gaussion mol2 fil
#O == original mol2 files
#@<TRIPOS>ATOM
# G G O O O O G G G
# G G O O O O G G G
#@<TRIPOS>BOND
#MERGED FILE EXAMPLE:
#@<TRIPOS>ATOM
# 1 N1 -18.5878 -6.4372 8.0759 N.pl3 1 MOL -1.0597
# 2 H1 -18.4216 -6.0346 8.9973 H 1 MOL 0.3916
#@<TRIPOS>BOND
gau = linefilter(gaufile.splitlines(), '@<TRIPOS>ATOM', '@<TRIPOS>BOND')
org = linefilter(orgfile.splitlines(), '@<TRIPOS>ATOM', '@<TRIPOS>BOND')
for gauline, orgline in itertools.izip(gau, org) :
glist = gauline.split()
olist = orgline.split()
anbr = glist[0]
anbrF = anbr.rjust(7," ")+" "
atom = glist[1]
atomF = atom.ljust(10," ")+" "
x = olist[2]
xF = x.rjust(10," ")+" "
y = olist[3]
yF = y.rjust(10," ")+" "
z = olist[4]
zF = z.rjust(10," ")+" "
atype = olist[5]
atypeF = atype.ljust(10," ")+" "
molnbr = glist[6]
molnbrF = molnbr.rjust(10," ")+" "
molnam = glist[7]
molnamF = " "+molnam.ljust(10," ")+" "
charge = glist[8]
chargeF = charge.rjust(10," ")+" "
fout.write(anbrF+atomF+xF+yF+zF+atypeF+molnbrF+molnamF+chargeF+"
")
#Here we write the rest of the gaussian mol2 file to the new file
#Again we have to write the '@<TRIPOS>BOND' manually cause it is stripped by linefilter
fout.writelines('@<TRIPOS>BOND'+"
")
for line in linefilter(gaufile.splitlines(), '@<TRIPOS>BOND', 'EOF') :
fout.writelines(line+"
")
#finally save the new mol2 file.
fout.close()
def alterchg(gaufile, orgfile) :
sys.stdout.flush()
#open the new file
fout = open(IN+"-CHG.mol2", "w")
#write all from beginning to '@<TRIPOS>ATOM' from the gaussian mol2 to the new file
for line in linefilter(orgfile.splitlines(), 'BOF', '@<TRIPOS>ATOM') :
fout.writelines(line+"
")
#Write the '@<TRIPOS>ATOM' match string cause i is stripped in the linefilter
fout.writelines('@<TRIPOS>ATOM'+"
")
#Now lets write the collumns 0 to 7 from the original mol2file and take collumn 8
#from the gassian files
#So it should look like this
#
#ORIGINAL FILE EXSAMPLE:
#@<TRIPOS>ATOM
# 1 N1 -18.5878 -6.4372 8.0759 N.pl3 1 <1> -0.9900
# 2 H1 -18.4216 -6.0346 8.9973 H 1 <1> 0.3600
#@<TRIPOS>BOND
#GAUSSIAN FILE EXAMPLE:
#@<TRIPOS>ATOM
# 1 N1 5.1000 0.1990 -0.0330 N.3 1 MOL -1.0597
# 2 H1 5.9140 -0.3800 -0.1080 H 1 MOL 0.3916
#@<TRIPOS>BOND
#Merge table
#G == gaussion mol2 fil
#O == original mol2 files
#@<TRIPOS>ATOM
# O O O O O O O O G
# O O O O O O O O G
#@<TRIPOS>BOND
#MERGED FILE EXAMPLE:
#@<TRIPOS>ATOM
# 1 N1 -18.5878 -6.4372 8.0759 N.pl3 1 <1> -1.0597
# 2 H1 -18.4216 -6.0346 8.9973 H 1 <1> 0.3916
#@<TRIPOS>BOND
gau = linefilter(gaufile.splitlines(), '@<TRIPOS>ATOM', '@<TRIPOS>BOND')
org = linefilter(orgfile.splitlines(), '@<TRIPOS>ATOM', '@<TRIPOS>BOND')
for gauline, orgline in itertools.izip(gau, org) :
glist = gauline.split()
olist = orgline.split()
anbr = olist[0]
anbrF = anbr.rjust(7," ")+" "
atom = olist[1]
atomF = atom.ljust(10," ")+" "
x = olist[2]
xF = x.rjust(10," ")+" "
y = olist[3]
yF = y.rjust(10," ")+" "
z = olist[4]
zF = z.rjust(10," ")+" "
atype = olist[5]
atypeF = atype.ljust(10," ")+" "
molnbr = olist[6]
molnbrF = molnbr.rjust(10," ")+" "
molnam = olist[7]
molnamF = " "+molnam.ljust(10," ")+" "
charge = glist[8]
chargeF = charge.rjust(10," ")+" "
fout.write(anbrF+atomF+xF+yF+zF+atypeF+molnbrF+molnamF+chargeF+"
")
#Here we write the rest of the gaussian mol2 file to the new file
#Again we have to write the '@<TRIPOS>BOND' manually cause it is stripped by linefilter
fout.writelines('@<TRIPOS>BOND'+"
")
for line in linefilter(orgfile.splitlines(), '@<TRIPOS>BOND', 'EOF') :
fout.writelines(line+"
")
#finally save the new mol2 file.
fout.close()
#remove extension from elements, append .log and .mol2, and generate with antechamber
#sybyl conform mol2 files with resp charges calculated from the gaussian log files
for element in F :
IN = os.path.splitext(element)[0]
sys.stdout.write("\rantechamber -pf y -at sybyl -i %s.log -fi gout -o %s-g.mol2 -fo mol2 -c resp
" % (IN, IN))
sys.stdout.flush()
amberout = subprocess.Popen(["antechamber", "-i", IN+".log", "-fi", "gout", "-o", IN+"-g.mol2", "-fo", "mol2", "-c", "resp", "-pf", "y", "-at", "sybyl"], stderr=subprocess.PIPE, stdout=subprocess.PIPE)
stdout, stderr = amberout.communicate()
if stderr :
print stderr
print stdout
tmpfiles = ["esout", "punch", "qout", "QOUT"]
for remfile in tmpfiles :
os.remove(remfile)
#open the gaussian and original mol2 file
gaumol2 = open(IN+"-g.mol2", "r").read()
orgmol2 = open(IN+".mol2", "r").read()
if CRD :
altercrd(gaumol2,orgmol2)
if CHG :
alterchg(gaumol2,orgmol2)
print "If you want to make a amber prep file from your mol2 files, you can do this with"
print " # for i in *.mol2 ; do antechamber -pf y -fi mol2 -i ${i} -fo prepi -o ${i/.mol2/.prep} -at amber ; done"
print "Be shure to check the Atomnames. Numbering should start with 1"
print "If you have missing parameters, run parmcheck to generate a frcmod file with the missing parameters"
print " # for i in *.prep ; do parmchk -i ${i} -f prepi -o ${i/.prep/.frcmod} ; done"