getting the trajectory from tracing a walker

[Ibrahim Mohamed]

i am using this code to obtain a complete trajectory of nacl association 

from westpa.analysis import BasicMDTrajectory

from westpa.analysis import Run

run = Run.open('west.h5')

iter_num = int(input('iteration number: ')) # 32

walker_num = int(input('segment number: ')) # 18

iteration = run.iteration(iter_num)

walker = iteration.walker(walker_num)

trace = walker.trace()

trajectory = BasicMDTrajectory(traj_ext = '.xtc', top='step4.2_equilibration_nowater.gro', state_ext='.gro')

traj = trajectory(walker)

traj.save_xtc(f'traj_iter_{iter_num}_walker_{walker_num}')

run.close()

i have a file called step4.2_equilibration_nowater.gro in the common_files folder and the trajectory is .xtc (i am using GROMACS) then i tried to open one segment with the .gro file in VMD just to make sure that the files are correct and it opened.
However, this code produces an error:

Traceback (most recent call last):
  File "/lfs01/workdirs/cairo029u1/1PHD/WESTPA/trial/nacl_gromacs/gromacs/westpa_files/get_traj.py", line 14, in <module>
    traj = trajectory(walker)
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/westpa/analysis/trajectories.py", line 74, in __call__
    value = self.fget(obj, include_initpoint=include_initpoint, **kwargs)
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/westpa/analysis/trajectories.py", line 273, in fget
    frame = mdtraj.load(path, top=traj.top)
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/mdtraj/core/trajectory.py", line 430, in load
    t = loader(tmp_file, **kwargs)
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/mdtraj/formats/gro.py", line 101, in load_gro
    return f.read_as_traj(n_frames=n_frames, stride=stride,
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/mdtraj/formats/gro.py", line 222, in read_as_traj
    coordinates, time, unitcell_vectors = self.read(stride=stride, atom_indices=atom_indices)
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/mdtraj/formats/gro.py", line 279, in read
    frame_xyz, frame_box, frame_time = self._read_frame()
  File "/home/cairo029u1/.local/lib/python3.9/site-packages/mdtraj/formats/gro.py", line 368, in _read_frame
    assert self.n_atoms == int(line.strip())
AssertionError

what is the problem here?

Note: i modified the post_iter.sh to make it remove the pbc and the water (keep only two atoms Na and Cl) and saved the file as seg.xtc

[Ibrahim Mohamed]

i tried to add a print command for the assertion and this was the output from a single run of the code:

self.n_atoms =  2 int(line.strip()) = 2

self.n_atoms =  2 int(line.strip()) = 3368

i think this means that for the first segment the topology (self.n_atoms) and .gro file ( int(line.strip()) ) are the same but in the second segment the used .gro file is the one from the simulation (inside the segment folder).

does this mean that i need to modify the .gro file in each segment folder to have the same atom number as the top file (step4.2_equilibration_nowater.gro)? if yes is not there a more easier way to do it?

[Ibrahim Mohamed]

i tried to write a python code to concatenate the .xtc files from the file generated from w_trace command 

import mdtraj

trace_file = input('enter the name of the file: ')

bstate = './common_files/step4.2_equilibration_nowater.gro'

base_traj = './traj_segs/'

trajs = []

with open(trace_file+'.txt') as f:

    for l in f:

        if not l.startswith('#'):

            elems = l.split()

            iter_number = elems[0]

            seg_number = elems[1]

            if int(iter_number) != 0:

                iter_folder = ('0' * (6-len(iter_number)) )  + iter_number

                seg_folder =  ('0' * (6-len(seg_number)) )  + seg_number

                traj_file = base_traj + iter_folder + '/' + seg_folder + '/seg.xtc'

                trajs.append(traj_file)

#print(*trajs, sep='\n')

whole_traj = mdtraj.load(trajs, discard_overlapping_frames=True, top = bstate)

whole_traj.save_xtc(trace_file+'.xtc')

the code save a xtc file of all segments in the same order as in the txt file generated. However when i tried to see the trajectory on VMD i noticed that the ions were jumping at some frames, though the gmx trjconv was used after each iteration to remove the PBC. is this jumping is caused by the PBC or other reasons?

attached are the gro and the concatenated xtc file

[Leung, Jeremy]

Hi Ibrahim,

Let me answer in order.

1. Yes. The whole problem with loading with westpa.analysis is that you have two different system setups (one with water, one without) so you need two gro/topology files. Most molecular programs require trajectories to share the same topology to concatenate. I would make copies of your .gro restart files s.t. it would only contain the 2 atoms just like in your xtc.

2. The jumping is likely due to PBC and imaging. Since you removed the PBC information, I would try aligning Na+ all frames in vmd to the first frame's Na+ to see if that looked ok. Without the PBC information, it'll be impossible to correct the the imaging problem though.

-- JL

---
Jeremy M. G. Leung
PhD Candidate, Chemistry
Graduate Student Researcher, Chemistry (Chong Lab)
University of Pittsburgh | 219 Parkman Avenue, Pittsburgh, PA 15260
jml...@pitt.edu | [He, Him, His]

--
You received this message because you are subscribed to the Google Groups "westpa-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to westpa-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/westpa-users/78c25594-ab79-4090-b9b4-c5b943436fb7n%40googlegroups.com.
<traj_32_18_trace.xtc><step4.2_equilibration_nowater.gro>