diff --git a/pycgtool.py b/pycgtool.py
index c30d7150ff40327cfd2beda3e22a2fa349e67a4f..ebd3000544ead28601add8eb21422de170358ef4 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 d22d3ed2d8b6cf88446870cb41d2fb1d68fecce5..db65c372b3a14fac46ea135873892fc9371883a2 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 0000000000000000000000000000000000000000..c728e9b00310b32446263a0cbe8c3fef61148bf7
--- /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
Binary files /dev/null and b/test/data/sugar_out.xtc differ
diff --git a/test/test_pycgtool.py b/test/test_pycgtool.py
index 2b1edd275a4bf73771e3cd4c3dc60a797487e987..eb1bf0db8b7f6766b781a46c6509b9caa55d5640 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