if(!-var msf) getvar msf msg="No MSF file specified"

if(!-var mds) getvar mds msg="No MDS file specified"

if(!-var trj) getvar trj msg="No Trj file specified"

if(!-var pdb) getvar pdb msg="No ref PDB file specified"

#load force field

#include "loadff.inc"

#load md system

loadmsf file="$msf"

#load a reference for RMSD calculation

loadpdb file="$pdb" copy2ref1=true

#load MD core and setup

loadsetup file=$mds

#feed in trj

trjin file="$trj"

#remove files if already exists

if(-f rmsd.txt) exec rm -rf rmsd.txt

if(-f bbrmsd.txt) exec rm -rf bbrmsd.txt

#define mainpro

defcls main select="segn PRO* DNA* RNA* && not hydrogen"

defcls backbone select="(segn DNA* RNA* && name C3' C4' C5' P) \

                     || (segn PRO* && name C CA N O)"

ifrm=0

while($ifrm<$ntotframes){

  trjreadnextframe

  mdtime=$curmdtime

  if($mdtime<0) goto done #double check whether trjreadnextframe successfully or not

  #print progress for every 10ps

  tmod=$mdtime%10

  if($tmod==0) echo -hs "current MD time: $mdtime ps"

  #calculate RMSD for all heavy atoms comparing to initial structure

  calc rmsd fitgroup=main calcgroup=main

  echo "$mdtime $value" >> rmsd.txt

  #calculate RMSD for main chain atoms

  calc rmsd fitgroup="backbone" calcgroup="backbone"

  echo "$mdtime $value" >> bbrmsd.txt

  ifrm++

}

echo -hs "output files:"

echo -hs "  rmsd.txt    the RMSD fluctuation of whole structure refer to the initial structure"

echo -hs "  bbrmsd.txt  the RMSD fluctuation of backbone refer to the initial structure"

stop

:done

echo -hs "something went wrong when reading trjactory data"

abort