[gmx-users] make_ndx for separating all atoms in two monolayers in a bilayer

[Ananya Debnath]

Hi,

I want to create an index file which will contain indices of all the
atoms in two monolayers separately in a bilayer. For this purpose
I have copied a python script from a post in gromacs mailing list, and
did a simple modification. But after modification, make_ndx gives
me error like "Group 1 does not exist" though output from the .py looks alright.

If someone gives me the idea why this occurs, it will be a great help.
The script I am running is as follows:
------------------------------------------------------------------------------------------

#!/usr/bin/python

import sys,string

# Usage with a pipe:
# python dppc_index.py dppc_npat.pdb |  make_ndx -f dppc_npat.pdb -o
dppc_bilayer.ndx
#

########
# main #
########
# print usage
if len(sys.argv) < 2:
    sys.exit("Usage: dppc_index.py dppc_npat.pdb")

# define reference atom
REF_ATOM = "P8"

# read data
file = open(sys.argv[1],"r")
lines = file.readlines()
file.close()

# get Z mean (see
http://www.rcsb.org/pdb/file_formats/pdb/pdbguide2.2/part_62.html for
pdb file format)
Zsum = 0.0
nP = 0
for i in lines:
    if i[:6] == "ATOM  " and string.split(i[12:16])[0] == REF_ATOM:
        Zsum += float(i[46:54])
        nP += 1
Zmean = Zsum / float(nP)

# get upper and lower leaflet
lower = []
upper = []
for i in lines:
    if i[:6] == "ATOM  " and string.split(i[12:16])[0] == REF_ATOM:
        atom_num = int(i[6:11])
        Z = float(i[46:54])
        if Z < Zmean:
           for j in range(-7,43):
                  lower.append(atom_num+j)

        if Z > Zmean:for j in range(-7,43):
                  upper.append(atom_num+j)

 # print output for make_ndx
output = "keep 0\ndel 0\n"
output += "a "
for i in range(len(lower)):
    if i < len(lower)-1:
        output += "%i | a" % lower[i]
    else:
        output += "%i" % lower[i]
output += "\nname 0 lower_leaflet\n"
output += "a "
for i in range(len(upper)):
    if i < len(upper)-1:
        output += "%i | a" % upper[i]
    else:
        output += "%i" % upper[i]
output += "\nname 1 upper_leaflet\n"
output += "q\n"
print output,

------------------------------------------------------

 Thanks and Regards,
 Ananya.