Skip to content

Commit fd23412

Browse files
committed
Since we moved the command line input parsing to the main function, arg_dict is not needed. Removing arg_dict and replacing uses of it with args.variable_name.
1 parent 04cb247 commit fd23412

1 file changed

Lines changed: 72 additions & 97 deletions

File tree

CodeEntropy/main_mcc.py

Lines changed: 72 additions & 97 deletions
Original file line numberDiff line numberDiff line change
@@ -31,36 +31,36 @@ def main():
3131
*.csv = results from different calculateion,
3232
*.pkl - Pickled reduced universe for further analysis,
3333
*.out - detailed output such as matrix and spectra""")
34-
35-
34+
35+
3636
parser.add_argument('-f', '--top_traj_file',
3737
required=True,
3838
dest="filePath",
3939
action='store',
4040
nargs='+',
4141
help="Path to Structure/topology file (AMBER PRMTOP, GROMACS TPR which contains topology and dihedral information) followed by Trajectory file(s) (AMBER NETCDF or GROMACS TRR) you will need to output the coordinates and forces to the same file. Required.")
42-
parser.add_argument('-l', '--selectString',
42+
parser.add_argument('-l', '--selectString',
4343
action='store',
44-
dest="selectionString",
44+
dest="selection_string",
4545
type=str,
4646
default='all',
4747
help='Selection string for CodeEntropy such as protein or resid, refer to MDAnalysis.select_atoms for more information.')
48-
parser.add_argument('-b', '--begin',
49-
action="store",
50-
dest="start",
51-
help="Start analysing the trajectory from this frame index. Defaults to 0",
52-
default=0,
48+
parser.add_argument('-b', '--begin',
49+
action="store",
50+
dest="start",
51+
help="Start analysing the trajectory from this frame index. Defaults to 0",
52+
default=0,
5353
type= int)
54-
parser.add_argument('-e', '--end',
55-
action="store",
56-
dest="end",
57-
help="Stop analysing the trajectory at this frame index. Defaults to -1 (end of trajectory file)",
54+
parser.add_argument('-e', '--end',
55+
action="store",
56+
dest="end",
57+
help="Stop analysing the trajectory at this frame index. Defaults to -1 (end of trajectory file)",
5858
default=-1,
5959
type=int)
60-
parser.add_argument('-d', '--step',
61-
action="store",
62-
dest="step",
63-
help="interval between two consecutive frames to be read index. Defaults to 1",
60+
parser.add_argument('-d', '--step',
61+
action="store",
62+
dest="step",
63+
help="interval between two consecutive frames to be read index. Defaults to 1",
6464
default=1,
6565
type=int)
6666
parser.add_argument("-n","--bin_width",
@@ -69,69 +69,69 @@ def main():
6969
default=30,
7070
type=int,
7171
help="Bin width in degrees for making the histogram of the dihedral angles for the conformational entropy. Default: 30")
72-
parser.add_argument('-k', '--tempra',
73-
action="store",
74-
dest="temp",
75-
help="Temperature for entropy calculation (K). Default to 298.0 K",
76-
default=298.0,
72+
parser.add_argument('-k', '--tempra',
73+
action="store",
74+
dest="temp",
75+
help="Temperature for entropy calculation (K). Default to 298.0 K",
76+
default=298.0,
7777
type=float)
7878
parser.add_argument("-v","--verbose",
7979
action="store",
8080
dest="verbose",
8181
default=False,
8282
type=bool,
8383
help="True/False flag for noisy or quiet output. Default: False")
84-
parser.add_argument('-t', '--thread',
85-
action="store",
86-
dest="thread",
87-
help="How many multiprocess to use. Default 1 for single core execution.",
88-
default=1,
84+
parser.add_argument('-t', '--thread',
85+
action="store",
86+
dest="thread",
87+
help="How many multiprocess to use. Default 1 for single core execution.",
88+
default=1,
8989
type=int)
9090
parser.add_argument("-o","--out",
9191
action="store",
92-
dest ="outFile",
92+
dest ="outfile",
9393
default="outfile.out",
9494
help ="Name of the file where the output will be written. Default: outfile.out")
9595
parser.add_argument("-r","--resout",
9696
action="store",
97-
dest ="resOutFile",
97+
dest ="resfile",
9898
default="res_outfile.out",
9999
help ="Name of the file where the residue entropy output will be written. Default: res_outfile.out")
100100
parser.add_argument("-m", "--mout",
101101
action="store",
102-
dest ="moutFile",
102+
dest ="moutfile",
103103
default=None,
104104
help ="Name of the file where certain matrices will be written (default: None).")
105105

106-
parser.add_argument('-c', '--cutShell',
106+
parser.add_argument('-c', '--cutShell',
107107
action='store',
108108
dest="cutShell",
109-
default=None,
109+
default=None,
110110
type=float,
111111
help='include cutoff shell analysis, add cutoff distance in angstrom Default None will ust the RAD Algorithm')
112-
parser.add_argument('-p', '--pureAtomNum',
112+
parser.add_argument('-p', '--pureAtomNum',
113113
action='store',
114-
dest="puteAtomNum",
115-
default=1,
114+
dest="puteAtomNum",
115+
default=1,
116116
type=int,
117117
help='Reference molecule resid for system of pure liquid. Default to 1')
118-
parser.add_argument('-x', '--excludedResnames',
118+
parser.add_argument('-x', '--excludedResnames',
119119
dest="excludedResnames",
120-
action='store',
120+
action='store',
121121
nargs='+',
122-
default=None,
122+
default=None,
123123
help='exclude a list of molecule names from nearest non-like analysis. Default: None. Multiples are gathered into list.')
124124
parser.add_argument('-w', '--water',
125125
dest = "waterResnames",
126-
action='store',
126+
action='store',
127127
default='WAT',
128-
nargs='+',
128+
nargs='+',
129129
help='resname for water molecules. Default: WAT. Multiples are gathered into list.')
130130
parser.add_argument('-s', '--solvent',
131-
dest="solventResnames",
132-
action='store',
131+
dest="solventResnames",
132+
action='store',
133133
nargs='+',
134-
default=None,
134+
default=None,
135135
help='include resname of solvent molecules (case-sensitive) Default: None. Multiples are gathered into list.')
136136
parser.add_argument("--solContact",
137137
action="store_true",
@@ -144,55 +144,30 @@ def main():
144144
print('Command line arguments are ill-defined, please check the arguments')
145145
raise
146146

147-
147+
148148
############## REPLACE INPUTS ##############
149149
print("printing all input")
150150
for arg in vars(args):
151151
print(' {} {}'.format(arg, getattr(args, arg) or ''))
152152

153153
startTime = datetime.now()
154154

155-
# create dictonary of inputs to be passed to the main function
156-
arg_dict = {}
157-
158-
arg_dict['outfile'] = args.outFile
159-
arg_dict['temper'] = args.temp
160-
arg_dict['selection_string'] = args.selectionString
161-
arg_dict['start'] = args.start
162-
arg_dict['end'] = args.end
163-
arg_dict['step'] = args.step
164-
arg_dict['bin_width'] = args.bin_width
165-
arg_dict['verbose'] = args.verbose
166-
arg_dict['thread'] = args.thread
167-
arg_dict['moutFile'] = args.moutFile
168-
arg_dict['resfile'] = args.resOutFile
169-
170-
arg_dict['cutShell'] = args.cutShell
171-
arg_dict['puteAtomNum'] = args.puteAtomNum
172-
arg_dict['excludedResnames'] = args.excludedResnames
173-
arg_dict['waterResnames'] = args.waterResnames
174-
arg_dict['solventResnames'] = args.solventResnames
175-
176155
# Get topology and trajectory file names and make universe
177156
tprfile = args.filePath[0]
178157
trrfile = args.filePath[1:]
179158
u = mda.Universe(tprfile, trrfile)
180-
arg_dict['universe'] = u
181-
182-
# Define MDAnalysis Universe from inputs
183-
u = arg_dict['universe']
184159

185160
# Define bin_width for histogram from inputs
186-
bin_width = arg_dict['bin_width']
161+
bin_width = args.bin_width
187162

188163
# Define trajectory slicing from inputs
189-
start = arg_dict['start']
164+
start = args.start
190165
if start is None:
191166
start = 0
192-
end = arg_dict['end']
167+
end = args.end
193168
if end is None:
194169
end = -1
195-
step = arg_dict['step']
170+
step = args.step
196171
if step is None:
197172
step = 1
198173
# Count number of frames, easy if not slicing
@@ -210,22 +185,22 @@ def main():
210185
residue_results_df = pd.DataFrame(columns=['Molecule ID', 'Residue','Type', 'Result'])
211186

212187
# printing headings for output files
213-
with open(arg_dict['outfile'], "a") as out:
188+
with open(args.outfile, "a") as out:
214189
print("Molecule\tLevel\tType\tResult (J/mol/K)\n", file=out)
215190

216-
with open(arg_dict['resfile'], "a") as res:
191+
with open(args.resfile, "a") as res:
217192
print("Molecule\tResidue\tType\tResult (J/mol/K)\n", file=res)
218193

219194
# Reduce number of atoms in MDA universe to selection_string arg (default all atoms included)
220-
if arg_dict['selection_string'] == 'all':
195+
if args.selection_string == 'all':
221196
reduced_atom = u
222197
else:
223-
reduced_atom = MDAHelper.new_U_select_atom(u, arg_dict['selection_string'])
198+
reduced_atom = MDAHelper.new_U_select_atom(u, args.selection_string)
224199
reduced_atom_name = f"{len(reduced_atom.trajectory)}_frame_dump_atom_selection"
225200
MDAHelper.write_universe(reduced_atom, reduced_atom_name)
226201

227202
# Scan system for molecules and select levels (united atom, residue, polymer) for each
228-
number_molecules, levels = LF.select_levels(reduced_atom, arg_dict['verbose'])
203+
number_molecules, levels = LF.select_levels(reduced_atom, args.verbose)
229204

230205
# Loop over molecules
231206
for molecule in range(number_molecules):
@@ -264,28 +239,28 @@ def main():
264239

265240
## Vibrational entropy at every level
266241
# Get the force and torque matrices for the beads at the relevant level
267-
force_matrix, torque_matrix = LF.get_matrices(residue_container, level, arg_dict['verbose'], start, end, step, number_frames, highest_level)
242+
force_matrix, torque_matrix = LF.get_matrices(residue_container, level, args.verbose, start, end, step, number_frames, highest_level)
268243

269244
# Calculate the entropy from the diagonalisation of the matrices
270-
S_trans_residue = EF.vibrational_entropy(force_matrix, "force", arg_dict['temper'],highest_level)
245+
S_trans_residue = EF.vibrational_entropy(force_matrix, "force", args.temp,highest_level)
271246
S_trans += S_trans_residue
272247
print(f"S_trans_{level}_{residue} = {S_trans_residue}")
273248
new_row = pd.DataFrame({'Molecule ID': [molecule], 'Residue': [residue],
274249
'Type':['Transvibrational (J/mol/K)'],
275250
'Result': [S_trans_residue],})
276251
residue_results_df = pd.concat([residue_results_df, new_row], ignore_index=True)
277-
with open(arg_dict['resfile'], "a") as res:
252+
with open(args.resfile, "a") as res:
278253
print(molecule,"\t",residue,"\tTransvibration\t",S_trans_residue, file=res)
279254

280255

281-
S_rot_residue = EF.vibrational_entropy(torque_matrix, "torque", arg_dict['temper'], highest_level)
256+
S_rot_residue = EF.vibrational_entropy(torque_matrix, "torque", args.temp, highest_level)
282257
S_rot += S_rot_residue
283258
print(f"S_rot_{level}_{residue} = {S_rot_residue}")
284259
new_row = pd.DataFrame({'Molecule ID': [molecule], 'Residue': [residue],
285260
'Type':['Rovibrational (J/mol/K)'],
286261
'Result': [S_rot_residue],})
287262
residue_results_df = pd.concat([residue_results_df, new_row], ignore_index=True)
288-
with open(arg_dict['resfile'], "a") as res:
263+
with open(args.resfile, "a") as res:
289264
# print(new_row, file=res)
290265
print(molecule,"\t",residue,"\tRovibrational \t",S_rot_residue, file=res)
291266

@@ -304,7 +279,7 @@ def main():
304279
'Type':['Conformational (J/mol/K)'],
305280
'Result': [S_conf_residue],})
306281
residue_results_df = pd.concat([residue_results_df, new_row], ignore_index=True)
307-
with open(arg_dict['resfile'], "a") as res:
282+
with open(args.resfile, "a") as res:
308283
print(molecule,"\t",residue,"\tConformational\t",S_conf_residue, file=res)
309284

310285

@@ -313,7 +288,7 @@ def main():
313288
new_row = pd.DataFrame({'Molecule ID': [molecule], 'Level': [level],
314289
'Type':['Transvibrational (J/mol/K)'],
315290
'Result': [S_trans],})
316-
with open(arg_dict['outfile'], "a") as out:
291+
with open(args.outfile, "a") as out:
317292
print(molecule,"\t",level,"\tTransvibration\t",S_trans, file=out)
318293

319294
results_df = pd.concat([results_df, new_row], ignore_index=True)
@@ -323,7 +298,7 @@ def main():
323298
'Type':['Rovibrational (J/mol/K)'],
324299
'Result': [S_rot],})
325300
results_df = pd.concat([results_df, new_row], ignore_index=True)
326-
with open(arg_dict['outfile'], "a") as out:
301+
with open(args.outfile, "a") as out:
327302
print(molecule,"\t",level,"\tRovibrational \t",S_rot, file=out)
328303

329304

@@ -332,7 +307,7 @@ def main():
332307
'Type':['Conformational (J/mol/K)'],
333308
'Result': [S_conf],})
334309
results_df = pd.concat([results_df, new_row], ignore_index=True)
335-
with open(arg_dict['outfile'], "a") as out:
310+
with open(args.outfile, "a") as out:
336311
print(molecule,"\t",level,"\tConformational\t",S_conf, file=out)
337312

338313

@@ -341,27 +316,27 @@ def main():
341316
if level in ('polymer', 'residue'):
342317
## Vibrational entropy at every level
343318
# Get the force and torque matrices for the beads at the relevant level
344-
force_matrix, torque_matrix = LF.get_matrices(molecule_container, level, arg_dict['verbose'], start, end, step, number_frames, highest_level)
319+
force_matrix, torque_matrix = LF.get_matrices(molecule_container, level, args.verbose, start, end, step, number_frames, highest_level)
345320

346321
# Calculate the entropy from the diagonalisation of the matrices
347-
S_trans = EF.vibrational_entropy(force_matrix, "force", arg_dict['temper'],highest_level)
322+
S_trans = EF.vibrational_entropy(force_matrix, "force", args.temp,highest_level)
348323
print(f"S_trans_{level} = {S_trans}")
349324
new_row = pd.DataFrame({'Molecule ID': [molecule], 'Level': [level],
350325
'Type':['Transvibrational (J/mol/K)'],
351326
'Result': [S_trans],})
352327
results_df = pd.concat([results_df, new_row], ignore_index=True)
353-
with open(arg_dict['outfile'], "a") as out:
328+
with open(args.outfile, "a") as out:
354329
print(molecule,"\t",level,"\tTransvibrational\t",S_trans, file=out)
355330

356331

357332

358-
S_rot = EF.vibrational_entropy(torque_matrix, "torque", arg_dict['temper'], highest_level)
333+
S_rot = EF.vibrational_entropy(torque_matrix, "torque", args.temp, highest_level)
359334
print(f"S_rot_{level} = {S_rot}")
360335
new_row = pd.DataFrame({'Molecule ID': [molecule], 'Level': [level],
361336
'Type':['Rovibrational (J/mol/K)'],
362337
'Result': [S_rot],})
363338
results_df = pd.concat([results_df, new_row], ignore_index=True)
364-
with open(arg_dict['outfile'], "a") as out:
339+
with open(args.outfile, "a") as out:
365340
print(molecule,"\t",level,"\tRovibrational \t",S_rot, file=out)
366341

367342

@@ -381,7 +356,7 @@ def main():
381356
'Type':['Conformational (J/mol/K)'],
382357
'Result': [S_conf],})
383358
results_df = pd.concat([results_df, new_row], ignore_index=True)
384-
with open(arg_dict['outfile'], "a") as out:
359+
with open(args.outfile, "a") as out:
385360
print(molecule,"\t",level,"\tConformational\t",S_conf, file=out)
386361

387362

@@ -395,7 +370,7 @@ def main():
395370
# 'Type':['Orientational (J/mol/K)'],
396371
# 'Result': [S_orient],})
397372
# results_df = pd.concat([results_df, new_row], ignore_index=True)
398-
# with open(arg_dict['outfile'], "a") as out:
373+
# with open(args.outfile, "a") as out:
399374
# print(molecule,"\t",level,"\tOrientational\t",S_orient, file=out)
400375

401376
# Report total entropy for the molecule
@@ -405,12 +380,12 @@ def main():
405380
'Type':['Molecule Total Entropy '],
406381
'Result': [S_molecule],})
407382
results_df = pd.concat([results_df, new_row], ignore_index=True)
408-
with open(arg_dict['outfile'], "a") as out:
383+
with open(args.outfile, "a") as out:
409384
print(molecule,"\t Molecule\tTotal Entropy\t",S_molecule, file=out)
410385

411386

412387
# END main function
413388

414389
if __name__ == "__main__":
415-
416-
main()
390+
391+
main()

0 commit comments

Comments
 (0)