# YASARA MACRO
# TOPIC:       3. Molecular Dynamics
# TITLE:       Convert between SIM and XTC simulation trajectories
# REQUIRES:    Dynamics 9.5.10
# AUTHOR:      Elmar Krieger
# LICENSE:     GPL
# DESCRIPTION: This macro converts an existing MD trajectory between various formats. Supported are SIM->XTC, SIM->PDB, SIM->SIM, XTC->SIM and XTC->PDB.

# Parameter section - adjust as needed
# ====================================

# The trajectory to convert must be present with a .sim or .xtc extension.
# The starting scene *_water.sce is also required.
# You can either set the target by clicking on Options > Macro > Set target,
# or by uncommenting the line below and specifying it directly.
#MacroTarget = 'c:\MyProject\1crn'

# Source format (srcformat) can be 'sim' (see SaveSim/LoadSim) or 'xtc' (see SaveXTC/LoadXTC).
# Destination format (dstformat) can be 'sim', 'xtc', 'pdb' (a series of PDB files)
# or 'pdbw' (a series of wrapped PDB files, where all atoms are inside the cell
# and potentially wrapped around periodic boundaries (i.e. broken molecules)).
srcformat='sim'
dstformat='pdb'

# Forcefield to use 
ForceField AMBER03

# In case the simulation was run without 'CorrectDrift on' and the solute diffused
# through a periodic boundary, you can keep it centered here by specifying the number
# of an atom close to the core of the solute, which will be kept at the cell center.
# If both srcformat and dstformat are 'sim', you can correct an existing trajectory.
central=0

# No changes required below this point!

# Do we have a target?
if MacroTarget==''
  RaiseError "This macro requires a target. Either edit the macro file or click Options > Macro > Set target to choose a target structure"
if srcformat!='sim' and srcformat!='xtc'
  RaiseError "The source format must be either 'sim' or 'xtc'"
if srcformat==dstformat
  if srcformat=='xtc'
    RaiseError "A conversion from 'xtc' to 'xtc' is not possible"
  elif not central
    RaiseError "A conversion from 'sim' to 'sim' only makes sense if an atom should be kept centered in the cell ('central')"
# Speed up conversion using short dummy cutoff and no longrange forces
Cutoff 2.62
Longrange None
# Load starting structure
LoadSce (MacroTarget)_water
bnd = Boundary
if dstformat=='pdb' or dstformat=='pdbw'
  # Since we can only save one object per PDB file, join all objects with atoms
  JoinObj Atom all,Atom 1
# Number of first snapshot
first=00000
# Load the first two snapshots to calculate the saving interval (steps)
for i=0 to 1
  if srcformat=='sim'
    LoadSim (MacroTarget)(first+i)
  else
    LoadXTC (MacroTarget)00000,(first+i)
  t(i) = Time
steps=0+(t1-t0)
TimeStep 1,1
if dstformat=='xtc'
  # SaveXTC does currently only extend but not overwrite an existing trajectory
  DelFile (MacroTarget)00000.xtc
if dstformat!='pdb' and dstformat!='pdbw'
  # Save the destination trajectory, since our timestep is 1fs, Steps is t (rounded)
  Save(dstformat) (MacroTarget)00000,Steps=(steps)
else
  # When saving PDB files, we can save time by pausing
  Sim Pause
Console Off
# As saving is done automatically for sim and xtc, all we need to do is load
# the snapshots from the source trajectory. This will adjust the current
# simulation time, which in turn triggers the save events. The 'Wait 1' is
# essential to suspend the macro for one cycle, so that the simulation
# continues and can be saved. (Snapshots are only saved while the simulation is running!).
i=first
do
  if srcformat=='sim'
    LoadSim (MacroTarget)(i)
    found = FileSize (MacroTarget)(i+1).sim
  else
    last = LoadXTC (MacroTarget)00000,(i+1)
    found=!last
  t = Time
  ShowMessage 'Converting (srcformat) snapshot (i) from (srcformat) to (dstformat) format, time is (0+t/1000) ps, (steps) fs/snapshot...'
  # Just in case the trajectory is not continuous, adjust the time
  Time ((i-first)*steps)
  if central
    # Keep a chosen atom at the center of the cell (Cell returns center as values 7-9) 
    _,_,_,_,_,_,cen() = Cell
    pos() = PosAtom (central)
    MoveAtom all,(cen-pos)
  if dstformat=='pdb' or dstformat=='pdbw'
    # To unwrap the soup, we need to stop the simulation
    Sim Off
    # By transfering the soup into the cell, we make sure that the current cell
    # is copied into the saved PDB file as CRYST1 (see SavePDB docs). This also
    # sets a new transformation history, so we don't have to turn off transformation.
    TransferObj Atom 1,SimCell,Local=Fix
    if dstformat=='pdb'
      # And transform back into a non-periodic cell to avoid shifts
      Boundary Wall
    Sim On
    # SavePDB expects an object selection, so 'Atom 1' selects the object with the first atom
    SavePDB Atom 1,(MacroTarget)(i)
    Sim Off
    Boundary (bnd)
    Time (t)
  # We need to proceed to save
  Wait 1
  i=i+1
while found    
HideMessage
Sim Off