Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add time step parameter to XTC Writers and Readers #4905

Open
pbuslaev opened this issue Feb 3, 2025 · 1 comment
Open

Add time step parameter to XTC Writers and Readers #4905

pbuslaev opened this issue Feb 3, 2025 · 1 comment

Comments

@pbuslaev
Copy link

pbuslaev commented Feb 3, 2025

Is your feature request related to a problem?

In some cases (when xtc is created with DESMOND) the timestep saved in xtc file is incorrect. It can be easier to set up the correct timestep for the loaded trajectory.

This feature is available for DCD writer and reader, and at least for the reader it should be straightforward to implement. I can work on it next week.

Describe the solution you'd like

Implement the support of timestep setting for XTC writer and reader.

Describe alternatives you've considered

My current workaround is to load xtc, then write dcd with the desired timestep and then read it back.

Additional context

@orbeckst
Copy link
Member

orbeckst commented Feb 5, 2025

Create a PR and then we can see what this would look like.

XTC and TRR files do not store the dt but include each individual time (and we just calculate dt from the first two frames) so I am not sure that the same approach as for DCD will work.

You could try a custom on-the-fly transformation https://userguide.mdanalysis.org/stable/trajectories/transformations.html where you set the Timestep.time to frame * time: ts.time = ts.frame * dt. For example:

import MDAnalysis as mda; import MDAnalysisTests.datafiles as data

u = mda.Universe(data.TPR, data.XTC)

print("Original dt", u.trajectory.dt)
print("Original times", [ts.time for ts in u.trajectory])

gives

Original dt 100.00000762939453
Original times [0.0, 100.00000762939453, 200.00001525878906, 300.0, 400.0000305175781, 500.0000305175781, 600.0, 700.0000610351562, 800.0000610351562, 900.0000610351562]

Now we are adding the transformation

def set_time(dt):
    """Set time to frame * dt"""
    def wrapped(ts):
        ts.dt = dt
        ts.time = ts.frame * dt
        return ts
    return wrapped

u.trajectory.add_transformations(set_time(33.33))

print("New dt", u.trajectory.dt)
print("New times", [ts.time for ts in u.trajectory])

and we get

New dt 33.33
New times [0.0, 33.33, 66.66, 99.99, 133.32, 166.64999999999998, 199.98, 233.31, 266.64, 299.96999999999997]

(You can then analyze the trajectory or write it out to another XTC but with correct time stamps.)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants