2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team,
6 * check out http://www.gromacs.org for more information.
7 * Copyright (c) 2012, by the GROMACS development team, led by
8 * David van der Spoel, Berk Hess, Erik Lindahl, and including many
9 * others, as listed in the AUTHORS file in the top-level source
10 * directory and at http://www.gromacs.org.
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
60 gmx_bool bCase = FALSE;
62 static int or_groups(atom_id nr1, atom_id *at1, atom_id nr2, atom_id *at2,
63 atom_id *nr, atom_id *at)
65 atom_id i1, i2, max = 0;
71 for (i1 = 0; i1 < nr1; i1++)
73 if ((i1 > 0) && (at1[i1] <= max))
79 for (i1 = 0; i1 < nr2; i1++)
81 if ((i1 > 0) && (at2[i1] <= max))
90 printf("One of your groups is not ascending\n");
97 while ((i1 < nr1) || (i2 < nr2))
99 if ((i2 == nr2) || ((i1 < nr1) && (at1[i1] < at2[i2])))
107 if ((i2 < nr2) && ((i1 == nr1) || (at1[i1] > at2[i2])))
116 printf("Merged two groups with OR: %u %u -> %u\n", nr1, nr2, *nr);
122 static int and_groups(atom_id nr1, atom_id *at1, atom_id nr2, atom_id *at2,
123 atom_id *nr, atom_id *at)
128 for (i1 = 0; i1 < nr1; i1++)
130 for (i2 = 0; i2 < nr2; i2++)
132 if (at1[i1] == at2[i2])
140 printf("Merged two groups with AND: %u %u -> %u\n", nr1, nr2, *nr);
145 static gmx_bool is_name_char(char c)
147 /* This string should contain all characters that can not be
148 * the first letter of a name due to the make_ndx syntax.
150 const char *spec = " !&|";
152 return (c != '\0' && strchr(spec, c) == NULL);
155 static int parse_names(char **string, int *n_names, char **names)
160 while (is_name_char((*string)[0]) || (*string)[0] == ' ')
162 if (is_name_char((*string)[0]))
164 if (*n_names >= MAXNAMES)
166 gmx_fatal(FARGS, "To many names: %d\n", *n_names+1);
169 while (is_name_char((*string)[i]))
171 names[*n_names][i] = (*string)[i];
175 printf("Name is too long, the maximum is %d characters\n", NAME_LEN);
179 names[*n_names][i] = '\0';
182 upstring(names[*n_names]);
196 static gmx_bool parse_int_char(char **string, int *nr, char *c)
203 while ((*string)[0] == ' ')
212 if (isdigit((*string)[0]))
214 *nr = (*string)[0]-'0';
216 while (isdigit((*string)[0]))
218 *nr = (*nr)*10+(*string)[0]-'0';
221 if (isalpha((*string)[0]))
226 /* Check if there is at most one non-digit character */
227 if (!isalnum((*string)[0]))
244 static gmx_bool parse_int(char **string, int *nr)
250 bRet = parse_int_char(string, nr, &c);
251 if (bRet && c != ' ')
260 static gmx_bool isquote(char c)
265 static gmx_bool parse_string(char **string, int *nr, int ngrps, char **grpname)
270 while ((*string)[0] == ' ')
276 if (isquote((*string)[0]))
280 s = strdup((*string));
284 (*string) += sp-s + 1;
286 (*nr) = find_group(s, ngrps, grpname);
290 return (*nr) != NOTSET;
293 static int select_atomnumbers(char **string, t_atoms *atoms, atom_id n1,
294 atom_id *nr, atom_id *index, char *gname)
300 while ((*string)[0] == ' ')
304 if ((*string)[0] == '-')
307 parse_int(string, &up);
308 if ((n1 < 1) || (n1 > atoms->nr) || (up < 1) || (up > atoms->nr))
310 printf("Invalid atom range\n");
314 for (i = n1-1; i <= up-1; i++)
319 printf("Found %u atom%s in range %u-%d\n", *nr, (*nr == 1) ? "" : "s", n1, up);
322 sprintf(buf, "a_%u", n1);
326 sprintf(buf, "a_%u-%d", n1, up);
337 if ((i-1 >= 0) && (i-1 < atoms->nr))
341 sprintf(buf, "_%d", i);
346 printf("Invalid atom number %d\n", i);
350 while ((*nr != 0) && (parse_int(string, &i)));
356 static int select_residuenumbers(char **string, t_atoms *atoms,
358 atom_id *nr, atom_id *index, char *gname)
365 while ((*string)[0] == ' ')
369 if ((*string)[0] == '-')
371 /* Residue number range selection */
374 printf("Error: residue insertion codes can not be used with residue range selection\n");
378 parse_int(string, &up);
380 for (i = 0; i < atoms->nr; i++)
382 ri = &atoms->resinfo[atoms->atom[i].resind];
383 for (j = n1; (j <= up); j++)
385 if (ri->nr == j && (c == ' ' || ri->ic == c))
392 printf("Found %u atom%s with res.nr. in range %u-%d\n",
393 *nr, (*nr == 1) ? "" : "s", n1, up);
396 sprintf(buf, "r_%u", n1);
400 sprintf(buf, "r_%u-%d", n1, up);
406 /* Individual residue number/insertion code selection */
411 for (i = 0; i < atoms->nr; i++)
413 ri = &atoms->resinfo[atoms->atom[i].resind];
414 if (ri->nr == j && ri->ic == c)
420 sprintf(buf, "_%d", j);
423 while (parse_int_char(string, &j, &c));
429 static int select_residueindices(char **string, t_atoms *atoms,
431 atom_id *nr, atom_id *index, char *gname)
433 /*this should be similar to select_residuenumbers except select by index (sequential numbering in file)*/
434 /*resind+1 for 1-indexing*/
440 while ((*string)[0] == ' ')
444 if ((*string)[0] == '-')
446 /* Residue number range selection */
449 printf("Error: residue insertion codes can not be used with residue range selection\n");
453 parse_int(string, &up);
455 for (i = 0; i < atoms->nr; i++)
457 ri = &atoms->resinfo[atoms->atom[i].resind];
458 for (j = n1; (j <= up); j++)
460 if (atoms->atom[i].resind+1 == j && (c == ' ' || ri->ic == c))
467 printf("Found %u atom%s with resind.+1 in range %u-%d\n",
468 *nr, (*nr == 1) ? "" : "s", n1, up);
471 sprintf(buf, "r_%u", n1);
475 sprintf(buf, "r_%u-%d", n1, up);
481 /* Individual residue number/insertion code selection */
486 for (i = 0; i < atoms->nr; i++)
488 ri = &atoms->resinfo[atoms->atom[i].resind];
489 if (atoms->atom[i].resind+1 == j && ri->ic == c)
495 sprintf(buf, "_%d", j);
498 while (parse_int_char(string, &j, &c));
505 static gmx_bool atoms_from_residuenumbers(t_atoms *atoms, int group, t_blocka *block,
506 atom_id *nr, atom_id *index, char *gname)
508 int i, j, j0, j1, resnr, nres;
510 j0 = block->index[group];
511 j1 = block->index[group+1];
513 for (j = j0; j < j1; j++)
515 if (block->a[j] >= nres)
517 printf("Index %s contains number>nres (%d>%d)\n",
518 gname, block->a[j]+1, nres);
522 for (i = 0; i < atoms->nr; i++)
524 resnr = atoms->resinfo[atoms->atom[i].resind].nr;
525 for (j = j0; j < j1; j++)
527 if (block->a[j]+1 == resnr)
535 printf("Found %u atom%s in %d residues from group %s\n",
536 *nr, (*nr == 1) ? "" : "s", j1-j0, gname);
540 static gmx_bool comp_name(char *name, char *search)
542 while (name[0] != '\0' && search[0] != '\0')
550 if (search[1] != '\0')
552 printf("WARNING: Currently '*' is only supported at the end of an expression\n");
561 /* Compare a single character */
562 if (( bCase && strncmp(name, search, 1)) ||
563 (!bCase && gmx_strncasecmp(name, search, 1)))
572 return (name[0] == '\0' && (search[0] == '\0' || search[0] == '*'));
575 static int select_chainnames(t_atoms *atoms, int n_names, char **names,
576 atom_id *nr, atom_id *index)
584 for (i = 0; i < atoms->nr; i++)
586 name[0] = atoms->resinfo[atoms->atom[i].resind].chainid;
588 while (j < n_names && !comp_name(name, names[j]))
598 printf("Found %u atom%s with chain identifier%s",
599 *nr, (*nr == 1) ? "" : "s", (n_names == 1) ? "" : "s");
600 for (j = 0; (j < n_names); j++)
602 printf(" %s", names[j]);
609 static int select_atomnames(t_atoms *atoms, int n_names, char **names,
610 atom_id *nr, atom_id *index, gmx_bool bType)
617 for (i = 0; i < atoms->nr; i++)
621 name = *(atoms->atomtype[i]);
625 name = *(atoms->atomname[i]);
628 while (j < n_names && !comp_name(name, names[j]))
638 printf("Found %u atoms with %s%s",
639 *nr, bType ? "type" : "name", (n_names == 1) ? "" : "s");
640 for (j = 0; (j < n_names); j++)
642 printf(" %s", names[j]);
649 static int select_residuenames(t_atoms *atoms, int n_names, char **names,
650 atom_id *nr, atom_id *index)
657 for (i = 0; i < atoms->nr; i++)
659 name = *(atoms->resinfo[atoms->atom[i].resind].name);
661 while (j < n_names && !comp_name(name, names[j]))
671 printf("Found %u atoms with residue name%s", *nr, (n_names == 1) ? "" : "s");
672 for (j = 0; (j < n_names); j++)
674 printf(" %s", names[j]);
681 static void copy2block(int n, atom_id *index, t_blocka *block)
688 srenew(block->index, block->nr+1);
689 block->index[block->nr] = n0+n;
690 srenew(block->a, n0+n);
691 for (i = 0; (i < n); i++)
693 block->a[n0+i] = index[i];
697 static void make_gname(int n, char **names, char *gname)
701 strcpy(gname, names[0]);
702 for (i = 1; i < n; i++)
705 strcat(gname, names[i]);
709 static void copy_group(int g, t_blocka *block, atom_id *nr, atom_id *index)
713 i0 = block->index[g];
714 *nr = block->index[g+1]-i0;
715 for (i = 0; i < *nr; i++)
717 index[i] = block->a[i0+i];
721 static void remove_group(int nr, int nr2, t_blocka *block, char ***gn)
731 for (j = 0; j <= nr2-nr; j++)
733 if ((nr < 0) || (nr >= block->nr))
735 printf("Group %d does not exist\n", nr+j);
739 shift = block->index[nr+1]-block->index[nr];
740 for (i = block->index[nr+1]; i < block->nra; i++)
742 block->a[i-shift] = block->a[i];
745 for (i = nr; i < block->nr; i++)
747 block->index[i] = block->index[i+1]-shift;
749 name = strdup((*gn)[nr]);
751 for (i = nr; i < block->nr-1; i++)
753 (*gn)[i] = (*gn)[i+1];
756 block->nra = block->index[block->nr];
757 printf("Removed group %d '%s'\n", nr+j, name);
763 static void split_group(t_atoms *atoms, int sel_nr, t_blocka *block, char ***gn,
766 char buf[STRLEN], *name;
770 printf("Splitting group %d '%s' into %s\n", sel_nr, (*gn)[sel_nr],
771 bAtom ? "atoms" : "residues");
773 n0 = block->index[sel_nr];
774 n1 = block->index[sel_nr+1];
775 srenew(block->a, block->nra+n1-n0);
776 for (i = n0; i < n1; i++)
779 resind = atoms->atom[a].resind;
780 name = *(atoms->resinfo[resind].name);
781 if (bAtom || (i == n0) || (atoms->atom[block->a[i-1]].resind != resind))
785 block->index[block->nr] = block->nra;
788 srenew(block->index, block->nr+1);
789 srenew(*gn, block->nr);
792 sprintf(buf, "%s_%s_%u", (*gn)[sel_nr], *atoms->atomname[a], a+1);
796 sprintf(buf, "%s_%s_%d", (*gn)[sel_nr], name, atoms->resinfo[resind].nr);
798 (*gn)[block->nr-1] = strdup(buf);
800 block->a[block->nra] = a;
803 block->index[block->nr] = block->nra;
806 static int split_chain(t_atoms *atoms, rvec *x,
807 int sel_nr, t_blocka *block, char ***gn)
811 atom_id i, a, natoms, *start = NULL, *end = NULL, ca_start, ca_end;
818 while (ca_start < natoms)
820 while ((ca_start < natoms) && strcmp(*atoms->atomname[ca_start], "CA"))
824 if (ca_start < natoms)
826 srenew(start, nchain+1);
827 srenew(end, nchain+1);
828 start[nchain] = ca_start;
829 while ((start[nchain] > 0) &&
830 (atoms->atom[start[nchain]-1].resind ==
831 atoms->atom[ca_start].resind))
844 while ((i < natoms) && strcmp(*atoms->atomname[i], "CA"));
847 rvec_sub(x[ca_end], x[i], vec);
850 while ((i < natoms) && (norm(vec) < 0.45));
852 end[nchain] = ca_end;
853 while ((end[nchain]+1 < natoms) &&
854 (atoms->atom[end[nchain]+1].resind == atoms->atom[ca_end].resind))
858 ca_start = end[nchain]+1;
864 printf("Found 1 chain, will not split\n");
868 printf("Found %d chains\n", nchain);
870 for (j = 0; j < nchain; j++)
872 printf("%d:%6u atoms (%u to %u)\n",
873 j+1, end[j]-start[j]+1, start[j]+1, end[j]+1);
878 srenew(block->a, block->nra+block->index[sel_nr+1]-block->index[sel_nr]);
879 for (j = 0; j < nchain; j++)
882 srenew(block->index, block->nr+1);
883 srenew(*gn, block->nr);
884 sprintf(buf, "%s_chain%d", (*gn)[sel_nr], j+1);
885 (*gn)[block->nr-1] = strdup(buf);
886 for (i = block->index[sel_nr]; i < block->index[sel_nr+1]; i++)
889 if ((a >= start[j]) && (a <= end[j]))
891 block->a[block->nra] = a;
895 block->index[block->nr] = block->nra;
896 if (block->index[block->nr-1] == block->index[block->nr])
898 remove_group(block->nr-1, NOTSET, block, gn);
908 static gmx_bool check_have_atoms(t_atoms *atoms, char *string)
912 printf("Can not process '%s' without atom info, use option -f\n", string);
921 static gmx_bool parse_entry(char **string, int natoms, t_atoms *atoms,
922 t_blocka *block, char ***gn,
923 atom_id *nr, atom_id *index, char *gname)
925 static char **names, *ostring;
926 static gmx_bool bFirst = TRUE;
927 int j, n_names, sel_nr1;
928 atom_id i, nr1, *index1;
930 gmx_bool bRet, bCompl;
935 snew(names, MAXNAMES);
936 for (i = 0; i < MAXNAMES; i++)
938 snew(names[i], NAME_LEN+1);
945 while (*string[0] == ' ')
950 if ((*string)[0] == '!')
954 while (*string[0] == ' ')
966 if (parse_int(string, &sel_nr1) ||
967 parse_string(string, &sel_nr1, block->nr, *gn))
969 if ((sel_nr1 >= 0) && (sel_nr1 < block->nr))
971 copy_group(sel_nr1, block, nr, index);
972 strcpy(gname, (*gn)[sel_nr1]);
973 printf("Copied index group %d '%s'\n", sel_nr1, (*gn)[sel_nr1]);
978 printf("Group %d does not exist\n", sel_nr1);
981 else if ((*string)[0] == 'a')
984 if (check_have_atoms(atoms, ostring))
986 if (parse_int(string, &sel_nr1))
988 bRet = select_atomnumbers(string, atoms, sel_nr1, nr, index, gname);
990 else if (parse_names(string, &n_names, names))
992 bRet = select_atomnames(atoms, n_names, names, nr, index, FALSE);
993 make_gname(n_names, names, gname);
997 else if ((*string)[0] == 't')
1000 if (check_have_atoms(atoms, ostring) &&
1001 parse_names(string, &n_names, names))
1003 if (atoms->atomtype == NULL)
1005 printf("Need a run input file to select atom types\n");
1009 bRet = select_atomnames(atoms, n_names, names, nr, index, TRUE);
1010 make_gname(n_names, names, gname);
1014 else if (strncmp(*string, "res", 3) == 0)
1017 if (check_have_atoms(atoms, ostring) &&
1018 parse_int(string, &sel_nr1) &&
1019 (sel_nr1 >= 0) && (sel_nr1 < block->nr) )
1021 bRet = atoms_from_residuenumbers(atoms,
1022 sel_nr1, block, nr, index, (*gn)[sel_nr1]);
1023 sprintf(gname, "atom_%s", (*gn)[sel_nr1]);
1026 else if (strncmp(*string, "ri", 2) == 0)
1029 if (check_have_atoms(atoms, ostring) &&
1030 parse_int_char(string, &sel_nr1, &c))
1032 bRet = select_residueindices(string, atoms, sel_nr1, c, nr, index, gname);
1035 else if ((*string)[0] == 'r')
1038 if (check_have_atoms(atoms, ostring))
1040 if (parse_int_char(string, &sel_nr1, &c))
1042 bRet = select_residuenumbers(string, atoms, sel_nr1, c, nr, index, gname);
1044 else if (parse_names(string, &n_names, names))
1046 bRet = select_residuenames(atoms, n_names, names, nr, index);
1047 make_gname(n_names, names, gname);
1051 else if (strncmp(*string, "chain", 5) == 0)
1054 if (check_have_atoms(atoms, ostring) &&
1055 parse_names(string, &n_names, names))
1057 bRet = select_chainnames(atoms, n_names, names, nr, index);
1058 sprintf(gname, "ch%s", names[0]);
1059 for (i = 1; i < n_names; i++)
1061 strcat(gname, names[i]);
1067 snew(index1, natoms-*nr);
1069 for (i = 0; i < natoms; i++)
1072 while ((j < *nr) && (index[j] != i))
1078 if (nr1 >= natoms-*nr)
1080 printf("There are double atoms in your index group\n");
1088 for (i = 0; i < nr1; i++)
1090 index[i] = index1[i];
1094 for (i = strlen(gname)+1; i > 0; i--)
1096 gname[i] = gname[i-1];
1099 printf("Complemented group: %u atoms\n", *nr);
1105 static void list_residues(t_atoms *atoms)
1107 int i, j, start, end, prev_resind, resind;
1110 /* Print all the residues, assuming continuous resnr count */
1111 start = atoms->atom[0].resind;
1112 prev_resind = start;
1113 for (i = 0; i < atoms->nr; i++)
1115 resind = atoms->atom[i].resind;
1116 if ((resind != prev_resind) || (i == atoms->nr-1))
1118 if ((bDiff = strcmp(*atoms->resinfo[resind].name,
1119 *atoms->resinfo[start].name)) ||
1132 for (j = start; j <= end; j++)
1135 j+1, *(atoms->resinfo[j].name));
1140 printf(" %4d - %4d %-5s ",
1141 start+1, end+1, *(atoms->resinfo[start].name));
1146 prev_resind = resind;
1151 static void edit_index(int natoms, t_atoms *atoms, rvec *x, t_blocka *block, char ***gn, gmx_bool bVerbose)
1153 static char **atnames, *ostring;
1154 static gmx_bool bFirst = TRUE;
1155 char inp_string[STRLEN], *string;
1156 char gname[STRLEN], gname1[STRLEN], gname2[STRLEN];
1157 int i, i0, i1, sel_nr, sel_nr2, newgroup;
1158 atom_id nr, nr1, nr2, *index, *index1, *index2;
1159 gmx_bool bAnd, bOr, bPrintOnce;
1164 snew(atnames, MAXNAMES);
1165 for (i = 0; i < MAXNAMES; i++)
1167 snew(atnames[i], NAME_LEN+1);
1173 snew(index, natoms);
1174 snew(index1, natoms);
1175 snew(index2, natoms);
1182 if (bVerbose || bPrintOnce || newgroup != NOTSET)
1185 if (bVerbose || bPrintOnce || newgroup == NOTSET)
1195 for (i = i0; i < i1; i++)
1197 printf("%3d %-20s: %5u atoms\n", i, (*gn)[i],
1198 block->index[i+1]-block->index[i]);
1202 if (bVerbose || bPrintOnce)
1205 printf(" nr : group ! 'name' nr name 'splitch' nr Enter: list groups\n");
1206 printf(" 'a': atom & 'del' nr 'splitres' nr 'l': list residues\n");
1207 printf(" 't': atom type | 'keep' nr 'splitat' nr 'h': help\n");
1208 printf(" 'r': residue 'res' nr 'chain' char\n");
1209 printf(" \"name\": group 'case': case %s 'q': save and quit\n",
1210 bCase ? "insensitive" : "sensitive ");
1211 printf(" 'ri': residue index\n");
1216 if (NULL == fgets(inp_string, STRLEN, stdin))
1218 gmx_fatal(FARGS, "Error reading user input");
1220 inp_string[strlen(inp_string)-1] = 0;
1222 string = inp_string;
1223 while (string[0] == ' ')
1230 if (string[0] == 'h')
1232 printf(" nr : selects an index group by number or quoted string.\n");
1233 printf(" The string is first matched against the whole group name,\n");
1234 printf(" then against the beginning and finally against an\n");
1235 printf(" arbitrary substring. A multiple match is an error.\n");
1237 printf(" 'a' nr1 [nr2 ...] : selects atoms, atom numbering starts at 1.\n");
1238 printf(" 'a' nr1 - nr2 : selects atoms in the range from nr1 to nr2.\n");
1239 printf(" 'a' name1[*] [name2[*] ...] : selects atoms by name(s), '?' matches any char,\n");
1240 printf(" wildcard '*' allowed at the end of a name.\n");
1241 printf(" 't' type1[*] [type2[*] ...] : as 'a', but for type, run input file required.\n");
1242 printf(" 'r' nr1[ic1] [nr2[ic2] ...] : selects residues by number and insertion code.\n");
1243 printf(" 'r' nr1 - nr2 : selects residues in the range from nr1 to nr2.\n");
1244 printf(" 'r' name1[*] [name2[*] ...] : as 'a', but for residue names.\n");
1245 printf(" 'ri' nr1 - nr2 : selects residue indices, 1-indexed, (as opposed to numbers) in the range from nr1 to nr2.\n");
1246 printf(" 'chain' ch1 [ch2 ...] : selects atoms by chain identifier(s),\n");
1247 printf(" not available with a .gro file as input.\n");
1248 printf(" ! : takes the complement of a group with respect to all\n");
1249 printf(" the atoms in the input file.\n");
1250 printf(" & | : AND and OR, can be placed between any of the options\n");
1251 printf(" above, the input is processed from left to right.\n");
1252 printf(" 'name' nr name : rename group nr to name.\n");
1253 printf(" 'del' nr1 [- nr2] : deletes one group or groups in the range from nr1 to nr2.\n");
1254 printf(" 'keep' nr : deletes all groups except nr.\n");
1255 printf(" 'case' : make all name compares case (in)sensitive.\n");
1256 printf(" 'splitch' nr : split group into chains using CA distances.\n");
1257 printf(" 'splitres' nr : split group into residues.\n");
1258 printf(" 'splitat' nr : split group into atoms.\n");
1259 printf(" 'res' nr : interpret numbers in group as residue numbers\n");
1260 printf(" Enter : list the currently defined groups and commands\n");
1261 printf(" 'l' : list the residues.\n");
1262 printf(" 'h' : show this help.\n");
1263 printf(" 'q' : save and quit.\n");
1265 printf(" Examples:\n");
1266 printf(" > 2 | 4 & r 3-5\n");
1267 printf(" selects all atoms from group 2 and 4 that have residue numbers 3, 4 or 5\n");
1268 printf(" > a C* & !a C CA\n");
1269 printf(" selects all atoms starting with 'C' but not the atoms 'C' and 'CA'\n");
1270 printf(" > \"protein\" & ! \"backb\"\n");
1271 printf(" selects all atoms that are in group 'protein' and not in group 'backbone'\n");
1274 printf("\npress Enter ");
1278 else if (strncmp(string, "del", 3) == 0)
1281 if (parse_int(&string, &sel_nr))
1283 while (string[0] == ' ')
1287 if (string[0] == '-')
1290 parse_int(&string, &sel_nr2);
1296 while (string[0] == ' ')
1300 if (string[0] == '\0')
1302 remove_group(sel_nr, sel_nr2, block, gn);
1306 printf("\nSyntax error: \"%s\"\n", string);
1310 else if (strncmp(string, "keep", 4) == 0)
1313 if (parse_int(&string, &sel_nr))
1315 remove_group(sel_nr+1, block->nr-1, block, gn);
1316 remove_group(0, sel_nr-1, block, gn);
1319 else if (strncmp(string, "name", 4) == 0)
1322 if (parse_int(&string, &sel_nr))
1324 if ((sel_nr >= 0) && (sel_nr < block->nr))
1326 sscanf(string, "%s", gname);
1327 sfree((*gn)[sel_nr]);
1328 (*gn)[sel_nr] = strdup(gname);
1332 else if (strncmp(string, "case", 4) == 0)
1335 printf("Switched to case %s\n", bCase ? "sensitive" : "insensitive");
1337 else if (string[0] == 'v')
1339 bVerbose = !bVerbose;
1340 printf("Turned verbose %s\n", bVerbose ? "on" : "off");
1342 else if (string[0] == 'l')
1344 if (check_have_atoms(atoms, ostring) )
1346 list_residues(atoms);
1349 else if (strncmp(string, "splitch", 7) == 0)
1352 if (check_have_atoms(atoms, ostring) &&
1353 parse_int(&string, &sel_nr) &&
1354 (sel_nr >= 0) && (sel_nr < block->nr))
1356 split_chain(atoms, x, sel_nr, block, gn);
1359 else if (strncmp(string, "splitres", 8) == 0)
1362 if (check_have_atoms(atoms, ostring) &&
1363 parse_int(&string, &sel_nr) &&
1364 (sel_nr >= 0) && (sel_nr < block->nr))
1366 split_group(atoms, sel_nr, block, gn, FALSE);
1369 else if (strncmp(string, "splitat", 7) == 0)
1372 if (check_have_atoms(atoms, ostring) &&
1373 parse_int(&string, &sel_nr) &&
1374 (sel_nr >= 0) && (sel_nr < block->nr))
1376 split_group(atoms, sel_nr, block, gn, TRUE);
1379 else if (string[0] == '\0')
1383 else if (string[0] != 'q')
1387 if (parse_entry(&string, natoms, atoms, block, gn, &nr, index, gname))
1391 while (string[0] == ' ')
1398 if (string[0] == '&')
1402 else if (string[0] == '|')
1411 for (i = 0; i < nr; i++)
1413 index1[i] = index[i];
1415 strcpy(gname1, gname);
1416 if (parse_entry(&string, natoms, atoms, block, gn, &nr2, index2, gname2))
1420 or_groups(nr1, index1, nr2, index2, &nr, index);
1421 sprintf(gname, "%s_%s", gname1, gname2);
1425 and_groups(nr1, index1, nr2, index2, &nr, index);
1426 sprintf(gname, "%s_&_%s", gname1, gname2);
1431 while (bAnd || bOr);
1433 while (string[0] == ' ')
1439 printf("\nSyntax error: \"%s\"\n", string);
1443 copy2block(nr, index, block);
1444 srenew(*gn, block->nr);
1445 newgroup = block->nr-1;
1446 (*gn)[newgroup] = strdup(gname);
1450 printf("Group is empty\n");
1454 while (string[0] != 'q');
1461 static int block2natoms(t_blocka *block)
1466 for (i = 0; i < block->nra; i++)
1468 natoms = max(natoms, block->a[i]);
1475 void merge_blocks(t_blocka *dest, t_blocka *source)
1479 /* count groups, srenew and fill */
1482 dest->nr += source->nr;
1483 srenew(dest->index, dest->nr+1);
1484 for (i = 0; i < source->nr; i++)
1486 dest->index[i0+i] = nra0 + source->index[i];
1488 /* count atoms, srenew and fill */
1489 dest->nra += source->nra;
1490 srenew(dest->a, dest->nra);
1491 for (i = 0; i < source->nra; i++)
1493 dest->a[nra0+i] = source->a[i];
1496 /* terminate list */
1497 dest->index[dest->nr] = dest->nra;
1501 int gmx_make_ndx(int argc, char *argv[])
1503 const char *desc[] = {
1504 "Index groups are necessary for almost every gromacs program.",
1505 "All these programs can generate default index groups. You ONLY",
1506 "have to use [TT]make_ndx[tt] when you need SPECIAL index groups.",
1507 "There is a default index group for the whole system, 9 default",
1508 "index groups for proteins, and a default index group",
1509 "is generated for every other residue name.[PAR]",
1510 "When no index file is supplied, also [TT]make_ndx[tt] will generate the",
1512 "With the index editor you can select on atom, residue and chain names",
1514 "When a run input file is supplied you can also select on atom type.",
1515 "You can use NOT, AND and OR, you can split groups",
1516 "into chains, residues or atoms. You can delete and rename groups.[PAR]",
1517 "The atom numbering in the editor and the index file starts at 1."
1520 static int natoms = 0;
1521 static gmx_bool bVerbose = FALSE;
1523 { "-natoms", FALSE, etINT, {&natoms},
1524 "set number of atoms (default: read from coordinate or index file)" },
1525 { "-verbose", FALSE, etBOOL, {&bVerbose},
1526 "HIDDENVerbose output" }
1528 #define NPA asize(pa)
1533 const char *stxfile;
1535 const char *ndxoutfile;
1542 t_blocka *block, *block2;
1543 char **gnames, **gnames2;
1545 { efSTX, "-f", NULL, ffOPTRD },
1546 { efNDX, "-n", NULL, ffOPTRDMULT },
1547 { efNDX, "-o", NULL, ffWRITE }
1549 #define NFILE asize(fnm)
1551 parse_common_args(&argc, argv, 0, NFILE, fnm, NPA, pa, asize(desc), desc,
1554 stxfile = ftp2fn_null(efSTX, NFILE, fnm);
1555 if (opt2bSet("-n", NFILE, fnm))
1557 nndxin = opt2fns(&ndxinfiles, "-n", NFILE, fnm);
1563 ndxoutfile = opt2fn("-o", NFILE, fnm);
1564 bNatoms = opt2parg_bSet("-natoms", NPA, pa);
1566 if (!stxfile && !nndxin)
1568 gmx_fatal(FARGS, "No input files (structure or index)");
1574 get_stx_coordnum(stxfile, &(atoms->nr));
1575 init_t_atoms(atoms, atoms->nr, TRUE);
1578 fprintf(stderr, "\nReading structure file\n");
1579 read_stx_conf(stxfile, title, atoms, x, v, &ePBC, box);
1589 /* read input file(s) */
1590 block = new_blocka();
1592 printf("Going to read %d old index file(s)\n", nndxin);
1595 for (i = 0; i < nndxin; i++)
1597 block2 = init_index(ndxinfiles[i], &gnames2);
1598 srenew(gnames, block->nr+block2->nr);
1599 for (j = 0; j < block2->nr; j++)
1601 gnames[block->nr+j] = gnames2[j];
1604 merge_blocks(block, block2);
1606 sfree(block2->index);
1607 /* done_block(block2); */
1614 analyse(atoms, block, &gnames, FALSE, TRUE);
1619 natoms = block2natoms(block);
1620 printf("Counted atom numbers up to %d in index file\n", natoms);
1623 edit_index(natoms, atoms, x, block, &gnames, bVerbose);
1625 write_index(ndxoutfile, block, gnames);