GIT repositories md_scripts / master python / init_umbrella.py
master

Tree @master (Download .tar.gz)

init_umbrella.py @masterraw · history · blame

#!/usr/bin/python

import argparse
import os
import numpy
import shutil
import math



# defaults
#--------------------------------
traj = "../pull/run/traj.trr"
index = "index.ndx"
struct = "../pull/run/topol.tpr"
top = "../topology/topol.top"
mdp = "umbrella.mdp"
window = 0.1
do_g_dist = 1



# commandline parser
#--------------------------------
parser = argparse.ArgumentParser()
parser.add_argument("-s", help="struct")
parser.add_argument("-f", help="traj")
parser.add_argument("-n", help="index")
parser.add_argument("-p", help="top")
parser.add_argument("-w", help="window width")
parser.add_argument("--skip_g_dist", help='Dont calculate distances.', action='store_true')
args = parser.parse_args()

if args.s:
  struct = args.s
if args.f:
  traj = args.f
if args.n:
  index = args.n
if args.p:
  top = args.p
if args.skip_g_dist:
  do_g_dist = 0

# absolute paths:
traj=os.path.abspath(traj)
struct=os.path.abspath(struct)
index=os.path.abspath(index)
mdp=os.path.abspath(mdp)
top=os.path.abspath(top)



# g_dist
#--------------------------------
if do_g_dist > 0:
  g_dist_command = ('g_dist' 
      + ' -f ' + traj
      + ' -s ' + struct 
      + ' -n ' + index
      )

  os.system(g_dist_command)



# find frames
#--------------------------------
frame_list = ([])
dists = numpy.loadtxt('dist.xvg',comments='@',skiprows=8,usecols=(0,4))

# Maximum/minimum distance
maxdist = dists[:,1].max()
mindist = dists[:,1].min()

itr = math.ceil(mindist)
while itr < maxdist:
  diff = abs( dists[:,1] - itr )
  indice = diff.argmin()
  time = dists[indice,0]
  dist = dists[indice,1]
  frame_list.append([time,dist])
  print(str(itr) + '   ' + str(dist))
  itr = itr + window


# extract .gro-files and create run input files
#--------------------------------
for itr in frame_list:

  time = int(round(itr[0]))
  dist = round(itr[1],4)
  foldername = str(dist)
  filename = 'confin.gro'
  gromplog = 'grompp.log'

  # backup old directory 
  backup_num = 0
  backup_folder = foldername
  while os.path.exists(backup_folder):
    print(backup_folder + ' exists')
    backup_num = backup_num + 1
    backup_folder = '#' + foldername + '.' + str(backup_num) + '#'
  
  if backup_num > 0:
    shutil.move(foldername, backup_folder)


  # create new directory
  os.mkdir(foldername)
  os.chdir(foldername)


  # dump gro with trjconv
  trjconv_command = ('trjconv' 
      + ' -dump ' + str(time)
      + ' -o ' + filename
      + ' -b ' + str( max( [time - 1, 0] ) )
      + ' -f ' + traj
      + ' -s ' + struct 
      )
  os.system('echo 0 | ' + trjconv_command)


  # grompp
  grompp_command = ('grompp' 
      + ' -c ' + filename
      + ' -p ' + top
      + ' -f ' + mdp
      + ' -n ' + index
      + ' 2> grompp.log'
      )
  os.system(grompp_command)


  # back to main directory
  os.chdir('../')