Skip to content

Commit 9b16f8f

Browse files
committed
Improve input flexibility and show example of mapping option.
1 parent b6d5f07 commit 9b16f8f

File tree

3 files changed

+106
-59
lines changed

3 files changed

+106
-59
lines changed

README.md

Lines changed: 13 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -48,14 +48,25 @@ are the corresponding atom names in the perturbed state. This allows the mapping
4848
be generated programmatically by any suitable external tool.
4949

5050
```bash
51-
ghostly reference.prm7 reference.rst7 perturbed.prm7 perturbed.rst7 --log-level debug
51+
ghostly --reference reference.* --perturbed.* --mapping '{0: 0, 1: 4, 2: 3, 3: 2, 4: 1}' --log-level debug
5252
```
5353

54+
Alternatively, you can pass a stream file containing a perturbable [BioSimSpace](https://biosimspace.openbiosim.org)
55+
system, which already contains the merged end states, using the `--system` option.
56+
57+
```bash
58+
ghostly --system system.bss --log-level debug
59+
```
60+
61+
> [!NOTE]
62+
> The `--system` option takes precedence over the `--reference` and `--perturbed` options.
63+
5464
When finished, the program will output a [BioSimSpace](https://biosimspace.openbiosim.org)
5565
stream file for the perturbable molecule, or AMBER or GROMACS files for the two end states.
5666
The format can be specified using the the `--output-format` option. If you require input
5767
for a free-energy perturbation simulation, e.g. a hybrid GROMACS toplogy file, the you can
5868
use the stream file with [BioSimSpace](https://biosimspace.openbiosim.org) to generate the
59-
required input files.
69+
required input files. Note that when a stream file is used as input, then the ouput will
70+
always be a stream file, regardless of the `--output-format` option.
6071

6172
Please run `ghostly --help` for more details of other configuration options.

src/ghostly/_cli.py

Lines changed: 93 additions & 55 deletions
Original file line numberDiff line numberDiff line change
@@ -48,27 +48,30 @@ def run():
4848
)
4949

5050
parser.add_argument(
51-
"reference_topology",
51+
"--system",
5252
type=str,
53-
help="Topology file for the reference molecule.",
53+
help=(
54+
"Stream file for the perturbable system. This takes precedence over ",
55+
"the reference and perturbed molecule files.",
56+
),
57+
nargs=1,
58+
required=False,
5459
)
5560

5661
parser.add_argument(
57-
"reference_coordinates",
62+
"--reference",
5863
type=str,
59-
help="Coordinate file for the reference molecule.",
64+
help="Topology and coordinate file for the reference molecule.",
65+
nargs=2,
66+
required=False,
6067
)
6168

6269
parser.add_argument(
63-
"perturbed_topology",
70+
"--perturbed",
6471
type=str,
6572
help="Topology file for the perturbed molecule.",
66-
)
67-
68-
parser.add_argument(
69-
"perturbed_coordinates",
70-
type=str,
71-
help="Coordinate file for the perturbed molecule.",
73+
nargs=2,
74+
required=False,
7275
)
7376

7477
parser.add_argument(
@@ -151,23 +154,49 @@ def run():
151154
logger.remove()
152155
logger.add(sys.stderr, level=args.log_level.upper(), enqueue=True)
153156

154-
# Try to load the reference.
155-
try:
156-
reference = BSS.IO.readMolecules(
157-
[args.reference_topology, args.reference_coordinates]
158-
)[0]
159-
except Exception as e:
160-
logger.error(f"An error occurred while loading the reference molecule: {e}")
161-
sys.exit(1)
157+
# Make sure that we've got a steam file, or files for the
158+
# reference and perturbed molecules.
159+
if args.system is None:
160+
if args.reference is None or args.perturbed is None:
161+
logger.error(
162+
"You must provide either a stream file or both reference and perturbed "
163+
"molecule files."
164+
)
165+
sys.exit(1)
162166

163-
# Try to load the perturbed molecule.
164-
try:
165-
perturbed = BSS.IO.readMolecules(
166-
[args.perturbed_topology, args.perturbed_coordinates]
167-
)[0]
168-
except Exception as e:
169-
logger.error(f"An error occurred while loading the perturbed molecule: {e}")
170-
sys.exit(1)
167+
# Try to load the reference.
168+
try:
169+
reference = BSS.IO.readMolecules(args.reference)[0]
170+
except Exception as e:
171+
logger.error(f"An error occurred while loading the reference molecule: {e}")
172+
sys.exit(1)
173+
174+
# Try to load the perturbed molecule.
175+
try:
176+
perturbed = BSS.IO.readMolecules(args.perturbed)[0]
177+
except Exception as e:
178+
logger.error(f"An error occurred while loading the perturbed molecule: {e}")
179+
sys.exit(1)
180+
else:
181+
if args.reference is not None or args.perturbed is not None:
182+
logger.warning(
183+
"Stream file and molecule files provided. Stream file takes precedence."
184+
)
185+
186+
# Try to load the system from the stream file.
187+
try:
188+
system = sr.stream.load(args.system[0])
189+
except Exception as e:
190+
logger.error(
191+
f"An error occurred while loading the system from the stream file: {e}"
192+
)
193+
sys.exit(1)
194+
195+
if not isinstance(system, (sr.system.System, sr.legacy.System.System)):
196+
raise TypeError("The provided stream file does not contain a valid system.")
197+
198+
# The molecule is already merged, so set the mapping to None.
199+
args.mapping = None
171200

172201
# Try to parse the mapping.
173202
if args.mapping is not None:
@@ -212,20 +241,22 @@ def run():
212241
sys.exit(1)
213242

214243
# Try to merge the reference and perturbed molecules.
215-
try:
216-
merged = BSS.Align.merge(
217-
reference,
218-
perturbed,
219-
mapping=mapping,
220-
force=True,
221-
)
222-
except Exception as e:
223-
logger.error(f"An error occurred while merging the molecules: {e}")
224-
sys.exit(1)
244+
if args.system is None:
245+
try:
246+
merged = BSS.Align.merge(
247+
reference,
248+
perturbed,
249+
mapping=mapping,
250+
force=True,
251+
)
252+
except Exception as e:
253+
logger.error(f"An error occurred while merging the molecules: {e}")
254+
sys.exit(1)
225255

226256
# Try to apply the modifications.
227257
try:
228-
system = sr.system.System(merged.toSystem()._sire_object)
258+
if args.system is None:
259+
system = merged.toSystem()._sire_object
229260
system = modify(system)
230261
except Exception as e:
231262
logger.error(
@@ -234,21 +265,28 @@ def run():
234265
sys.exit(1)
235266

236267
# Try to save the system.
237-
try:
238-
if args.output_format == "bss":
239-
sr.stream.save(system, f"{args.output_prefix}.bss")
240-
else:
241-
if args.output_format == "amber":
242-
formats = ["prm7", "rst7"]
243-
elif args.output_format == "gromacs":
244-
formats = ["grotop", "gro87"]
245-
246-
reference = merged._toRegularMolecule()
247-
perturbed = merged._toRegularMolecule(is_lambda1=True)
248-
249-
BSS.IO.saveMolecules(args.output_prefix + "_reference", reference, formats)
250-
BSS.IO.saveMolecules(args.output_prefix + "_perturbed", perturbed, formats)
268+
if args.system is None:
269+
try:
270+
if args.output_format == "bss":
271+
sr.stream.save(system, f"{args.output_prefix}.bss")
272+
else:
273+
if args.output_format == "amber":
274+
formats = ["prm7", "rst7"]
275+
elif args.output_format == "gromacs":
276+
formats = ["grotop", "gro87"]
277+
278+
reference = merged._toRegularMolecule()
279+
perturbed = merged._toRegularMolecule(is_lambda1=True)
280+
281+
BSS.IO.saveMolecules(
282+
args.output_prefix + "_reference", reference, formats
283+
)
284+
BSS.IO.saveMolecules(
285+
args.output_prefix + "_perturbed", perturbed, formats
286+
)
251287

252-
except Exception as e:
253-
logger.error(f"An error occurred while saving the modified end states: {e}")
254-
sys.exit(1)
288+
except Exception as e:
289+
logger.error(f"An error occurred while saving the modified end states: {e}")
290+
sys.exit(1)
291+
else:
292+
sr.stream.save(system, f"{args.output_prefix}.bss")

src/ghostly/_ghostly.py

Lines changed: 0 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -80,8 +80,6 @@ def modify(system, k_hard=100, k_soft=5, optimise_angles=True):
8080
https://pubs.acs.org/doi/10.1021/acs.jctc.0c01328
8181
"""
8282

83-
_logger.info(f"Applying Boresch modifications to ghost atom bonded terms")
84-
8583
# Check the system is a Sire system.
8684
if not isinstance(system, (_System, _LegacySystem)):
8785
raise TypeError(

0 commit comments

Comments
 (0)