From 8964d1d078a70c6c218f7bab6c445fee0a7f6d03 Mon Sep 17 00:00:00 2001
From: James Graham <J.A.Graham@soton.ac.uk>
Date: Tue, 5 Jul 2016 16:46:01 +0100
Subject: [PATCH] Fixed stopping on first frame - test added Put main() and
map_only() into pycgtool module
---
pycgtool.py | 100 ++--------------------------------------
pycgtool/interface.py | 7 ++-
pycgtool/pycgtool.py | 98 +++++++++++++++++++++++++++++++++++++++
test/data/sugar_out.xtc | Bin 0 -> 1408 bytes
test/test_pycgtool.py | 45 ++++++++++++++++++
5 files changed, 151 insertions(+), 99 deletions(-)
create mode 100755 pycgtool/pycgtool.py
create mode 100644 test/data/sugar_out.xtc
diff --git a/pycgtool.py b/pycgtool.py
index c30d715..ebd3000 100755
--- a/pycgtool.py
+++ b/pycgtool.py
@@ -2,100 +2,9 @@
import argparse
import sys
-import logging
-
-from pycgtool.frame import Frame
-from pycgtool.mapping import Mapping
-from pycgtool.bondset import BondSet
-from pycgtool.forcefield import ForceField
-from pycgtool.interface import Progress, Options
-
-logger = logging.getLogger(__name__)
-
-
-def main(args, config):
- """
- Main function of the program PyCGTOOL.
-
- Performs the complete AA->CG mapping and outputs a files dependent on given input.
-
- :param args: Arguments from argparse
- :param config: Configuration dictionary
- """
- frame = Frame(gro=args.gro, xtc=args.xtc, itp=args.itp, frame_start=args.begin, xtc_reader="mdtraj")
-
- if args.bnd:
- bonds = BondSet(args.bnd, config)
- logger.info("Bond measurements will be made")
- else:
- logger.info("Bond measurements will not be made")
-
- if args.map:
- mapping = Mapping(args.map, config, itp=args.itp)
- cgframe = mapping.apply(frame, exclude={"SOL"})
- cgframe.output(config.output_name + ".gro", format=config.output)
- logger.info("Mapping will be performed")
- else:
- logger.info("Mapping will not be performed")
-
- # Main loop - perform mapping and measurement on every frame in XTC
- def main_loop():
- nonlocal cgframe
- frame.next_frame()
- if args.map:
- cgframe = mapping.apply(frame, cgframe=cgframe, exclude={"SOL"})
- if config.output_xtc:
- cgframe.write_xtc(config.output_name + ".xtc")
- else:
- cgframe = frame
- if args.bnd:
- bonds.apply(cgframe)
-
- numframes = frame.numframes - args.begin if args.end == -1 else args.end - args.begin
- logger.info("Beginning analysis of {0} frames".format(numframes))
- Progress(numframes, postwhile=main_loop).run()
-
- if args.bnd:
- if args.map:
- logger.info("Beginning Boltzmann inversion")
- bonds.boltzmann_invert()
- if config.output_forcefield:
- logger.info("Creating GROMACS forcefield directory")
- ff = ForceField("ff" + config.output_name + ".ff")
- ff.write_rtp(config.output_name + ".rtp", mapping, bonds)
- logger.info("GROMACS forcefield directory created")
- else:
- bonds.write_itp(config.output_name + ".itp", mapping=mapping)
-
- if config.dump_measurements:
- logger.info("Dumping bond measurements to file")
- bonds.dump_values(config.dump_n_values)
-
-
-def map_only(args, config):
- """
- Perform AA->CG mapping and output coordinate file.
-
- :param args: Program arguments
- :param config: Object containing run options
- """
- frame = Frame(gro=args.gro, xtc=args.xtc)
- mapping = Mapping(args.map, config)
- cgframe = mapping.apply(frame, exclude={"SOL"})
- cgframe.output(config.output_name + ".gro", format=config.output)
-
- if args.xtc and (config.output_xtc or args.outputxtc):
- # Main loop - perform mapping and measurement on every frame in XTC
- def main_loop():
- nonlocal cgframe
- frame.next_frame()
- cgframe = mapping.apply(frame, cgframe=cgframe, exclude={"SOL"})
- cgframe.write_xtc(config.output_name + ".xtc")
-
- numframes = frame.numframes - args.begin if args.end == -1 else args.end - args.begin
- logger.info("Beginning analysis of {0} frames".format(numframes))
- Progress(numframes, postwhile=main_loop).run()
+from pycgtool.pycgtool import main, map_only
+from pycgtool.interface import Options
if __name__ == "__main__":
parser = argparse.ArgumentParser(description="Perform coarse-grain mapping of atomistic trajectory")
@@ -108,6 +17,7 @@ if __name__ == "__main__":
parser.add_argument('--interactive', default=False, action='store_true')
parser.add_argument('--outputxtc', default=False, action='store_true')
+ parser.add_argument('--quiet', default=False, action='store_true')
input_files.add_argument('--begin', type=int, default=0, help="Frame number to begin")
input_files.add_argument('--end', type=int, default=-1, help="Frame number to end")
@@ -140,10 +50,6 @@ if __name__ == "__main__":
print("Using XTC: {0}".format(args.xtc))
if config.map_only:
- logger.info("Starting: mapping only")
map_only(args, config)
else:
- logger.info("Starting: parameter measurement")
main(args, config)
-
- logger.info("Analysis complete")
diff --git a/pycgtool/interface.py b/pycgtool/interface.py
index d22d3ed..db65c37 100644
--- a/pycgtool/interface.py
+++ b/pycgtool/interface.py
@@ -263,14 +263,17 @@ class Progress:
"""
Iterate through self until stopped by maximum iterations or False condition.
"""
- collections.deque(self, maxlen=0)
+ # collections.deque(self, maxlen=0)
+ for _ in self:
+ pass
return self._its
@property
def _bar(self):
done = int(self._length * (self._its / self._maxits))
left = self._length - done
- return "{0} [".format(self._its) + done * "#" + left * "-" + "] {0}".format(self._maxits)
+ width = len(str(self._maxits))
+ return "{0:-{width}} [".format(self._its, width=width) + done * "#" + left * "-" + "] {0}".format(self._maxits)
def _stop(self):
if not self._quiet:
diff --git a/pycgtool/pycgtool.py b/pycgtool/pycgtool.py
new file mode 100755
index 0000000..c728e9b
--- /dev/null
+++ b/pycgtool/pycgtool.py
@@ -0,0 +1,98 @@
+import logging
+
+from .frame import Frame
+from .mapping import Mapping
+from .bondset import BondSet
+from .forcefield import ForceField
+from .interface import Progress
+
+logger = logging.getLogger(__name__)
+
+
+def main(args, config):
+ """
+ Main function of the program PyCGTOOL.
+
+ Performs the complete AA->CG mapping and outputs a files dependent on given input.
+
+ :param args: Arguments from argparse
+ :param config: Configuration dictionary
+ """
+ frame = Frame(gro=args.gro, xtc=args.xtc, itp=args.itp, frame_start=args.begin, xtc_reader="mdtraj")
+
+ if args.bnd:
+ bonds = BondSet(args.bnd, config)
+ logger.info("Bond measurements will be made")
+ else:
+ logger.info("Bond measurements will not be made")
+
+ if args.map:
+ mapping = Mapping(args.map, config, itp=args.itp)
+ cgframe = mapping.apply(frame, exclude={"SOL"})
+ cgframe.output(config.output_name + ".gro", format=config.output)
+ logger.info("Mapping will be performed")
+ else:
+ logger.info("Mapping will not be performed")
+
+ # Main loop - perform mapping and measurement on every frame in XTC
+ def main_loop():
+ nonlocal cgframe
+ if not frame.next_frame():
+ return False
+ if args.map:
+ cgframe = mapping.apply(frame, cgframe=cgframe, exclude={"SOL"})
+ if config.output_xtc:
+ cgframe.write_xtc(config.output_name + ".xtc")
+ else:
+ cgframe = frame
+ if args.bnd:
+ bonds.apply(cgframe)
+ return True
+
+ numframes = frame.numframes - args.begin if args.end == -1 else args.end - args.begin
+ logger.info("Beginning analysis of {0} frames".format(numframes))
+ Progress(numframes, postwhile=main_loop, quiet=args.quiet).run()
+
+ if args.bnd:
+ if args.map:
+ logger.info("Beginning Boltzmann inversion")
+ bonds.boltzmann_invert()
+ if config.output_forcefield:
+ logger.info("Creating GROMACS forcefield directory")
+ ff = ForceField("ff" + config.output_name + ".ff")
+ ff.write_rtp(config.output_name + ".rtp", mapping, bonds)
+ logger.info("GROMACS forcefield directory created")
+ else:
+ bonds.write_itp(config.output_name + ".itp", mapping=mapping)
+
+ if config.dump_measurements:
+ logger.info("Dumping bond measurements to file")
+ bonds.dump_values(config.dump_n_values)
+
+
+def map_only(args, config):
+ """
+ Perform AA->CG mapping and output coordinate file.
+
+ :param args: Program arguments
+ :param config: Object containing run options
+ """
+ frame = Frame(gro=args.gro, xtc=args.xtc)
+ mapping = Mapping(args.map, config)
+ cgframe = mapping.apply(frame, exclude={"SOL"})
+ cgframe.output(config.output_name + ".gro", format=config.output)
+
+ if args.xtc and (config.output_xtc or args.outputxtc):
+ # Main loop - perform mapping and measurement on every frame in XTC
+ def main_loop():
+ nonlocal cgframe
+ if not frame.next_frame():
+ return False
+ cgframe = mapping.apply(frame, cgframe=cgframe, exclude={"SOL"})
+ cgframe.write_xtc(config.output_name + ".xtc")
+ return True
+
+ numframes = frame.numframes - args.begin if args.end == -1 else args.end - args.begin
+ logger.info("Beginning analysis of {0} frames".format(numframes))
+ Progress(numframes, postwhile=main_loop, quiet=args.quiet).run()
+
diff --git a/test/data/sugar_out.xtc b/test/data/sugar_out.xtc
new file mode 100644
index 0000000000000000000000000000000000000000..0f9b0822191d84aa887bedd2c1c71978a93fcb20
GIT binary patch
literal 1408
zcmZQzU_Z^kz`(}9z`)4Bz`)>;e{LO!hG7i8wk9h`y@Rd)0|%SGU+h;qPjfJQ-R5BT
zvdVt<^IQj`11%1QVjcE7Z>u<HZ*z0dHa4+e{QZK1(QZBmy_@3pox9>3ECiArEX{rF
zX9T`MvY!cLr$c_@E)WgF7<{cGZD9NDe%^7g@!V~{>OhNwQJk5BqgaRiW}$osqYX+9
z#!8p$_hdCV=-f4M&@S0)zufz{gQ@UY2c0#X_C21v9jsgB9n8|w?Wdo+g=9Z7$WDj+
zBeozKhB5fsoJ+y>+dd9)uykJsvETQKgZ<o1_S?0Y91Jr`!2aJ=_sc=MRo+4K)Is}Y
z?K%#oYqmRRUz}>+b)dz;;_G1tGZsnv>E6eY>}LVl>5zX?4Mf8*24Cy;0=WHl?K<}B
zZb~>9{(j-$`1^(ZUVRw{W9MlO#>ZCJ@4fWNLHl{GgP!Xm`{n0uIhehybTFt;x9?5e
z;b2*D&B1!<6Z={7_aoU4idP1Q{5j8Hc7tdPzINzTu>Fo=9S&9xJ?z)=g*lk~{o-KH
zWoN(R>tP4uZEg+*Uys=DthnZ&FSWoyt3Sqmamik=pLO?2*>^q@aIn>vb}-NBv7a%y
z6v=*28f0+D7uW`}7ltwT+6O`IV{@>tdFNp5Jk5U1+er>4Zznl82qfF@{k+}5a95mz
zajU%j{unKA`qRm3uwTPr<Y4q?nuC53lYQ^ROb5$d@nHLBD_umgpB-eUL%v@ah=ySd
zzK+-&u>T#B{2Z)g*4l4SI0{Nj4)#(D><=g$b1-<C07{Sc`{!<Q(0b_MpjURve&z2M
z4kn*D9CSZ#x9{C6<zUUk3yKf>Ia&2c_H%&jbjbTY0Yt+v24DMlDA;~G*F_F?C$`$J
ziZgRC@!aj;xGT<npQ(d`sePS;aY~f^zOqvex|{SIbk3Z$U-i(#!DM5-gF&vfeYf~i
z2Xihv2OAC}`<W4fNcMAr>~zRC^99i`jKSBw7zwr?6#kYO#r7+<SUH&7HE^)o5oEt3
z@Qs7<^IQjmq7C4(K>OGV2hEwY>{r<qI+#dIanO7wVBdXgg@eTwD+g2iTKnl-c1ZSf
zf$Vh1>+u88FpR<1{+11~Uz-V>hC$_zN!23<`*XMK_nbNFU_7<S!7!s3oECJtT^zKJ
zt*~EFw86pjMuvl~{{#D8zAy*NnX?=$4At#t&YXp0KR3uuhul|dK{O0w@O2VE`Gn2E
zmM_e~dTNvXTK@+QhH9@IY#)c%@0!2g!SK{U2ZN`H_WQzE9JIe*aL{yL2abP`{n{6&
S*!TRXbFdWa0O$W%4?O^_<c5|2
literal 0
HcmV?d00001
diff --git a/test/test_pycgtool.py b/test/test_pycgtool.py
index 2b1edd2..eb1bf0d 100644
--- a/test/test_pycgtool.py
+++ b/test/test_pycgtool.py
@@ -2,12 +2,57 @@ import unittest
import subprocess
import os
+import numpy as np
+
+from simpletraj.trajectory import XtcTrajectory
+
+from pycgtool.interface import Options
+from pycgtool.pycgtool import main, map_only
+
+
+class Args:
+ def __init__(self, name, map=True, bnd=True):
+ self.gro = os.path.join("test/data", name+".gro")
+ self.xtc = os.path.join("test/data", name+".xtc")
+ self.map = os.path.join("test/data", name+".map") if map else None
+ self.bnd = os.path.join("test/data", name+".bnd") if bnd else None
+ self.begin = 0
+ self.end = -1
+ self.quiet = True
+
class PycgtoolTest(unittest.TestCase):
+ config = Options([("output_name", "out"),
+ ("output", "gro"),
+ ("output_xtc", True),
+ ("map_only", False),
+ ("map_center", "geom"),
+ ("constr_threshold", 100000),
+ ("dump_measurements", False),
+ ("dump_n_values", 100000),
+ ("output_forcefield", False),
+ ("temperature", 310),
+ ("angle_default_fc", True),
+ ("generate_angles", True),
+ ("generate_dihedrals", False)])
+
def test_run(self):
path = os.path.dirname(os.path.dirname(__file__))
self.assertEqual(0, subprocess.check_call([os.path.join(path, "pycgtool.py"), "-h"], stdout=subprocess.PIPE))
+ def test_map_only(self):
+ map_only(Args("sugar"), self.config)
+
+ xtc = XtcTrajectory("out.xtc")
+ xtc_ref = XtcTrajectory("test/data/sugar_out.xtc")
+ self.assertEqual(xtc_ref.numframes, xtc.numframes)
+
+ for i in range(xtc_ref.numframes):
+ xtc.get_frame(i)
+ xtc_ref.get_frame(i)
+ np.testing.assert_array_almost_equal(xtc_ref.box, xtc.box, decimal=3)
+ np.testing.assert_array_almost_equal(xtc_ref.x, xtc.x, decimal=3)
+
# TODO more tests
--
GitLab