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.
59 gmx_bool bCase = FALSE;
61 static int or_groups(atom_id nr1, atom_id *at1, atom_id nr2, atom_id *at2,
62 atom_id *nr, atom_id *at)
64 atom_id i1, i2, max = 0;
70 for (i1 = 0; i1 < nr1; i1++)
72 if ((i1 > 0) && (at1[i1] <= max))
78 for (i1 = 0; i1 < nr2; i1++)
80 if ((i1 > 0) && (at2[i1] <= max))
89 printf("One of your groups is not ascending\n");
96 while ((i1 < nr1) || (i2 < nr2))
98 if ((i2 == nr2) || ((i1 < nr1) && (at1[i1] < at2[i2])))
106 if ((i2 < nr2) && ((i1 == nr1) || (at1[i1] > at2[i2])))
115 printf("Merged two groups with OR: %d %d -> %d\n", nr1, nr2, *nr);
121 static int and_groups(atom_id nr1, atom_id *at1, atom_id nr2, atom_id *at2,
122 atom_id *nr, atom_id *at)
127 for (i1 = 0; i1 < nr1; i1++)
129 for (i2 = 0; i2 < nr2; i2++)
131 if (at1[i1] == at2[i2])
139 printf("Merged two groups with AND: %d %d -> %d\n", nr1, nr2, *nr);
144 static gmx_bool is_name_char(char c)
146 /* This string should contain all characters that can not be
147 * the first letter of a name due to the make_ndx syntax.
149 const char *spec = " !&|";
151 return (c != '\0' && strchr(spec, c) == NULL);
154 static int parse_names(char **string, int *n_names, char **names)
159 while (is_name_char((*string)[0]) || (*string)[0] == ' ')
161 if (is_name_char((*string)[0]))
163 if (*n_names >= MAXNAMES)
165 gmx_fatal(FARGS, "To many names: %d\n", *n_names+1);
168 while (is_name_char((*string)[i]))
170 names[*n_names][i] = (*string)[i];
174 printf("Name is too long, the maximum is %d characters\n", NAME_LEN);
178 names[*n_names][i] = '\0';
181 upstring(names[*n_names]);
195 static gmx_bool parse_int_char(char **string, int *nr, char *c)
202 while ((*string)[0] == ' ')
211 if (isdigit((*string)[0]))
213 *nr = (*string)[0]-'0';
215 while (isdigit((*string)[0]))
217 *nr = (*nr)*10+(*string)[0]-'0';
220 if (isalpha((*string)[0]))
225 /* Check if there is at most one non-digit character */
226 if (!isalnum((*string)[0]))
243 static gmx_bool parse_int(char **string, int *nr)
249 bRet = parse_int_char(string, nr, &c);
250 if (bRet && c != ' ')
259 static gmx_bool isquote(char c)
264 static gmx_bool parse_string(char **string, int *nr, int ngrps, char **grpname)
269 while ((*string)[0] == ' ')
275 if (isquote((*string)[0]))
279 s = strdup((*string));
283 (*string) += sp-s + 1;
285 (*nr) = find_group(s, ngrps, grpname);
289 return (*nr) != NOTSET;
292 static int select_atomnumbers(char **string, t_atoms *atoms, atom_id n1,
293 atom_id *nr, atom_id *index, char *gname)
299 while ((*string)[0] == ' ')
303 if ((*string)[0] == '-')
306 parse_int(string, &up);
307 if ((n1 < 1) || (n1 > atoms->nr) || (up < 1) || (up > atoms->nr))
309 printf("Invalid atom range\n");
313 for (i = n1-1; i <= up-1; i++)
318 printf("Found %d atom%s in range %d-%d\n", *nr, (*nr == 1) ? "" : "s", n1, up);
321 sprintf(buf, "a_%d", n1);
325 sprintf(buf, "a_%d-%d", n1, up);
336 if ((i-1 >= 0) && (i-1 < atoms->nr))
340 sprintf(buf, "_%d", i);
345 printf("Invalid atom number %d\n", i);
349 while ((*nr != 0) && (parse_int(string, &i)));
355 static int select_residuenumbers(char **string, t_atoms *atoms,
357 atom_id *nr, atom_id *index, char *gname)
364 while ((*string)[0] == ' ')
368 if ((*string)[0] == '-')
370 /* Residue number range selection */
373 printf("Error: residue insertion codes can not be used with residue range selection\n");
377 parse_int(string, &up);
379 for (i = 0; i < atoms->nr; i++)
381 ri = &atoms->resinfo[atoms->atom[i].resind];
382 for (j = n1; (j <= up); j++)
384 if (ri->nr == j && (c == ' ' || ri->ic == c))
391 printf("Found %d atom%s with res.nr. in range %d-%d\n",
392 *nr, (*nr == 1) ? "" : "s", n1, up);
395 sprintf(buf, "r_%d", n1);
399 sprintf(buf, "r_%d-%d", n1, up);
405 /* Individual residue number/insertion code selection */
410 for (i = 0; i < atoms->nr; i++)
412 ri = &atoms->resinfo[atoms->atom[i].resind];
413 if (ri->nr == j && ri->ic == c)
419 sprintf(buf, "_%d", j);
422 while (parse_int_char(string, &j, &c));
428 static int select_residueindices(char **string, t_atoms *atoms,
430 atom_id *nr, atom_id *index, char *gname)
432 /*this should be similar to select_residuenumbers except select by index (sequential numbering in file)*/
433 /*resind+1 for 1-indexing*/
439 while ((*string)[0] == ' ')
443 if ((*string)[0] == '-')
445 /* Residue number range selection */
448 printf("Error: residue insertion codes can not be used with residue range selection\n");
452 parse_int(string, &up);
454 for (i = 0; i < atoms->nr; i++)
456 ri = &atoms->resinfo[atoms->atom[i].resind];
457 for (j = n1; (j <= up); j++)
459 if (atoms->atom[i].resind+1 == j && (c == ' ' || ri->ic == c))
466 printf("Found %d atom%s with resind.+1 in range %d-%d\n",
467 *nr, (*nr == 1) ? "" : "s", n1, up);
470 sprintf(buf, "r_%d", n1);
474 sprintf(buf, "r_%d-%d", n1, up);
480 /* Individual residue number/insertion code selection */
485 for (i = 0; i < atoms->nr; i++)
487 ri = &atoms->resinfo[atoms->atom[i].resind];
488 if (atoms->atom[i].resind+1 == j && ri->ic == c)
494 sprintf(buf, "_%d", j);
497 while (parse_int_char(string, &j, &c));
504 static gmx_bool atoms_from_residuenumbers(t_atoms *atoms, int group, t_blocka *block,
505 atom_id *nr, atom_id *index, char *gname)
507 int i, j, j0, j1, resnr, nres;
509 j0 = block->index[group];
510 j1 = block->index[group+1];
512 for (j = j0; j < j1; j++)
514 if (block->a[j] >= nres)
516 printf("Index %s contains number>nres (%d>%d)\n",
517 gname, block->a[j]+1, nres);
521 for (i = 0; i < atoms->nr; i++)
523 resnr = atoms->resinfo[atoms->atom[i].resind].nr;
524 for (j = j0; j < j1; j++)
526 if (block->a[j]+1 == resnr)
534 printf("Found %d atom%s in %d residues from group %s\n",
535 *nr, (*nr == 1) ? "" : "s", j1-j0, gname);
539 static gmx_bool comp_name(char *name, char *search)
541 while (name[0] != '\0' && search[0] != '\0')
549 if (search[1] != '\0')
551 printf("WARNING: Currently '*' is only supported at the end of an expression\n");
560 /* Compare a single character */
561 if (( bCase && strncmp(name, search, 1)) ||
562 (!bCase && gmx_strncasecmp(name, search, 1)))
571 return (name[0] == '\0' && (search[0] == '\0' || search[0] == '*'));
574 static int select_chainnames(t_atoms *atoms, int n_names, char **names,
575 atom_id *nr, atom_id *index)
583 for (i = 0; i < atoms->nr; i++)
585 name[0] = atoms->resinfo[atoms->atom[i].resind].chainid;
587 while (j < n_names && !comp_name(name, names[j]))
597 printf("Found %d atom%s with chain identifier%s",
598 *nr, (*nr == 1) ? "" : "s", (n_names == 1) ? "" : "s");
599 for (j = 0; (j < n_names); j++)
601 printf(" %s", names[j]);
608 static int select_atomnames(t_atoms *atoms, int n_names, char **names,
609 atom_id *nr, atom_id *index, gmx_bool bType)
616 for (i = 0; i < atoms->nr; i++)
620 name = *(atoms->atomtype[i]);
624 name = *(atoms->atomname[i]);
627 while (j < n_names && !comp_name(name, names[j]))
637 printf("Found %d atoms with %s%s",
638 *nr, bType ? "type" : "name", (n_names == 1) ? "" : "s");
639 for (j = 0; (j < n_names); j++)
641 printf(" %s", names[j]);
648 static int select_residuenames(t_atoms *atoms, int n_names, char **names,
649 atom_id *nr, atom_id *index)
656 for (i = 0; i < atoms->nr; i++)
658 name = *(atoms->resinfo[atoms->atom[i].resind].name);
660 while (j < n_names && !comp_name(name, names[j]))
670 printf("Found %d atoms with residue name%s", *nr, (n_names == 1) ? "" : "s");
671 for (j = 0; (j < n_names); j++)
673 printf(" %s", names[j]);
680 static void copy2block(int n, atom_id *index, t_blocka *block)
687 srenew(block->index, block->nr+1);
688 block->index[block->nr] = n0+n;
689 srenew(block->a, n0+n);
690 for (i = 0; (i < n); i++)
692 block->a[n0+i] = index[i];
696 static void make_gname(int n, char **names, char *gname)
700 strcpy(gname, names[0]);
701 for (i = 1; i < n; i++)
704 strcat(gname, names[i]);
708 static void copy_group(int g, t_blocka *block, atom_id *nr, atom_id *index)
712 i0 = block->index[g];
713 *nr = block->index[g+1]-i0;
714 for (i = 0; i < *nr; i++)
716 index[i] = block->a[i0+i];
720 static void remove_group(int nr, int nr2, t_blocka *block, char ***gn)
730 for (j = 0; j <= nr2-nr; j++)
732 if ((nr < 0) || (nr >= block->nr))
734 printf("Group %d does not exist\n", nr+j);
738 shift = block->index[nr+1]-block->index[nr];
739 for (i = block->index[nr+1]; i < block->nra; i++)
741 block->a[i-shift] = block->a[i];
744 for (i = nr; i < block->nr; i++)
746 block->index[i] = block->index[i+1]-shift;
748 name = strdup((*gn)[nr]);
750 for (i = nr; i < block->nr-1; i++)
752 (*gn)[i] = (*gn)[i+1];
755 block->nra = block->index[block->nr];
756 printf("Removed group %d '%s'\n", nr+j, name);
762 static void split_group(t_atoms *atoms, int sel_nr, t_blocka *block, char ***gn,
765 char buf[STRLEN], *name;
769 printf("Splitting group %d '%s' into %s\n", sel_nr, (*gn)[sel_nr],
770 bAtom ? "atoms" : "residues");
772 n0 = block->index[sel_nr];
773 n1 = block->index[sel_nr+1];
774 srenew(block->a, block->nra+n1-n0);
775 for (i = n0; i < n1; i++)
778 resind = atoms->atom[a].resind;
779 name = *(atoms->resinfo[resind].name);
780 if (bAtom || (i == n0) || (atoms->atom[block->a[i-1]].resind != resind))
784 block->index[block->nr] = block->nra;
787 srenew(block->index, block->nr+1);
788 srenew(*gn, block->nr);
791 sprintf(buf, "%s_%s_%d", (*gn)[sel_nr], *atoms->atomname[a], a+1);
795 sprintf(buf, "%s_%s_%d", (*gn)[sel_nr], name, atoms->resinfo[resind].nr);
797 (*gn)[block->nr-1] = strdup(buf);
799 block->a[block->nra] = a;
802 block->index[block->nr] = block->nra;
805 static int split_chain(t_atoms *atoms, rvec *x,
806 int sel_nr, t_blocka *block, char ***gn)
810 atom_id i, a, natoms, *start = NULL, *end = NULL, ca_start, ca_end;
817 while (ca_start < natoms)
819 while ((ca_start < natoms) && strcmp(*atoms->atomname[ca_start], "CA"))
823 if (ca_start < natoms)
825 srenew(start, nchain+1);
826 srenew(end, nchain+1);
827 start[nchain] = ca_start;
828 while ((start[nchain] > 0) &&
829 (atoms->atom[start[nchain]-1].resind ==
830 atoms->atom[ca_start].resind))
843 while ((i < natoms) && strcmp(*atoms->atomname[i], "CA"));
846 rvec_sub(x[ca_end], x[i], vec);
849 while ((i < natoms) && (norm(vec) < 0.45));
851 end[nchain] = ca_end;
852 while ((end[nchain]+1 < natoms) &&
853 (atoms->atom[end[nchain]+1].resind == atoms->atom[ca_end].resind))
857 ca_start = end[nchain]+1;
863 printf("Found 1 chain, will not split\n");
867 printf("Found %d chains\n", nchain);
869 for (j = 0; j < nchain; j++)
871 printf("%d:%6d atoms (%d to %d)\n",
872 j+1, end[j]-start[j]+1, start[j]+1, end[j]+1);
877 srenew(block->a, block->nra+block->index[sel_nr+1]-block->index[sel_nr]);
878 for (j = 0; j < nchain; j++)
881 srenew(block->index, block->nr+1);
882 srenew(*gn, block->nr);
883 sprintf(buf, "%s_chain%d", (*gn)[sel_nr], j+1);
884 (*gn)[block->nr-1] = strdup(buf);
885 for (i = block->index[sel_nr]; i < block->index[sel_nr+1]; i++)
888 if ((a >= start[j]) && (a <= end[j]))
890 block->a[block->nra] = a;
894 block->index[block->nr] = block->nra;
895 if (block->index[block->nr-1] == block->index[block->nr])
897 remove_group(block->nr-1, NOTSET, block, gn);
907 static gmx_bool check_have_atoms(t_atoms *atoms, char *string)
911 printf("Can not process '%s' without atom info, use option -f\n", string);
920 static gmx_bool parse_entry(char **string, int natoms, t_atoms *atoms,
921 t_blocka *block, char ***gn,
922 atom_id *nr, atom_id *index, char *gname)
924 static char **names, *ostring;
925 static gmx_bool bFirst = TRUE;
926 int j, n_names, sel_nr1;
927 atom_id i, nr1, *index1;
929 gmx_bool bRet, bCompl;
934 snew(names, MAXNAMES);
935 for (i = 0; i < MAXNAMES; i++)
937 snew(names[i], NAME_LEN+1);
944 while (*string[0] == ' ')
949 if ((*string)[0] == '!')
953 while (*string[0] == ' ')
965 if (parse_int(string, &sel_nr1) ||
966 parse_string(string, &sel_nr1, block->nr, *gn))
968 if ((sel_nr1 >= 0) && (sel_nr1 < block->nr))
970 copy_group(sel_nr1, block, nr, index);
971 strcpy(gname, (*gn)[sel_nr1]);
972 printf("Copied index group %d '%s'\n", sel_nr1, (*gn)[sel_nr1]);
977 printf("Group %d does not exist\n", sel_nr1);
980 else if ((*string)[0] == 'a')
983 if (check_have_atoms(atoms, ostring))
985 if (parse_int(string, &sel_nr1))
987 bRet = select_atomnumbers(string, atoms, sel_nr1, nr, index, gname);
989 else if (parse_names(string, &n_names, names))
991 bRet = select_atomnames(atoms, n_names, names, nr, index, FALSE);
992 make_gname(n_names, names, gname);
996 else if ((*string)[0] == 't')
999 if (check_have_atoms(atoms, ostring) &&
1000 parse_names(string, &n_names, names))
1002 if (atoms->atomtype == NULL)
1004 printf("Need a run input file to select atom types\n");
1008 bRet = select_atomnames(atoms, n_names, names, nr, index, TRUE);
1009 make_gname(n_names, names, gname);
1013 else if (strncmp(*string, "res", 3) == 0)
1016 if (check_have_atoms(atoms, ostring) &&
1017 parse_int(string, &sel_nr1) &&
1018 (sel_nr1 >= 0) && (sel_nr1 < block->nr) )
1020 bRet = atoms_from_residuenumbers(atoms,
1021 sel_nr1, block, nr, index, (*gn)[sel_nr1]);
1022 sprintf(gname, "atom_%s", (*gn)[sel_nr1]);
1025 else if (strncmp(*string, "ri", 2) == 0)
1028 if (check_have_atoms(atoms, ostring) &&
1029 parse_int_char(string, &sel_nr1, &c))
1031 bRet = select_residueindices(string, atoms, sel_nr1, c, nr, index, gname);
1034 else if ((*string)[0] == 'r')
1037 if (check_have_atoms(atoms, ostring))
1039 if (parse_int_char(string, &sel_nr1, &c))
1041 bRet = select_residuenumbers(string, atoms, sel_nr1, c, nr, index, gname);
1043 else if (parse_names(string, &n_names, names))
1045 bRet = select_residuenames(atoms, n_names, names, nr, index);
1046 make_gname(n_names, names, gname);
1050 else if (strncmp(*string, "chain", 5) == 0)
1053 if (check_have_atoms(atoms, ostring) &&
1054 parse_names(string, &n_names, names))
1056 bRet = select_chainnames(atoms, n_names, names, nr, index);
1057 sprintf(gname, "ch%s", names[0]);
1058 for (i = 1; i < n_names; i++)
1060 strcat(gname, names[i]);
1066 snew(index1, natoms-*nr);
1068 for (i = 0; i < natoms; i++)
1071 while ((j < *nr) && (index[j] != i))
1077 if (nr1 >= natoms-*nr)
1079 printf("There are double atoms in your index group\n");
1087 for (i = 0; i < nr1; i++)
1089 index[i] = index1[i];
1093 for (i = strlen(gname)+1; i > 0; i--)
1095 gname[i] = gname[i-1];
1098 printf("Complemented group: %d atoms\n", *nr);
1104 static void list_residues(t_atoms *atoms)
1106 int i, j, start, end, prev_resind, resind;
1109 /* Print all the residues, assuming continuous resnr count */
1110 start = atoms->atom[0].resind;
1111 prev_resind = start;
1112 for (i = 0; i < atoms->nr; i++)
1114 resind = atoms->atom[i].resind;
1115 if ((resind != prev_resind) || (i == atoms->nr-1))
1117 if ((bDiff = strcmp(*atoms->resinfo[resind].name,
1118 *atoms->resinfo[start].name)) ||
1131 for (j = start; j <= end; j++)
1134 j+1, *(atoms->resinfo[j].name));
1139 printf(" %4d - %4d %-5s ",
1140 start+1, end+1, *(atoms->resinfo[start].name));
1145 prev_resind = resind;
1150 static void edit_index(int natoms, t_atoms *atoms, rvec *x, t_blocka *block, char ***gn, gmx_bool bVerbose)
1152 static char **atnames, *ostring;
1153 static gmx_bool bFirst = TRUE;
1154 char inp_string[STRLEN], *string;
1155 char gname[STRLEN], gname1[STRLEN], gname2[STRLEN];
1156 int i, i0, i1, sel_nr, sel_nr2, newgroup;
1157 atom_id nr, nr1, nr2, *index, *index1, *index2;
1158 gmx_bool bAnd, bOr, bPrintOnce;
1163 snew(atnames, MAXNAMES);
1164 for (i = 0; i < MAXNAMES; i++)
1166 snew(atnames[i], NAME_LEN+1);
1172 snew(index, natoms);
1173 snew(index1, natoms);
1174 snew(index2, natoms);
1181 if (bVerbose || bPrintOnce || newgroup != NOTSET)
1184 if (bVerbose || bPrintOnce || newgroup == NOTSET)
1194 for (i = i0; i < i1; i++)
1196 printf("%3d %-20s: %5d atoms\n", i, (*gn)[i],
1197 block->index[i+1]-block->index[i]);
1201 if (bVerbose || bPrintOnce)
1204 printf(" nr : group ! 'name' nr name 'splitch' nr Enter: list groups\n");
1205 printf(" 'a': atom & 'del' nr 'splitres' nr 'l': list residues\n");
1206 printf(" 't': atom type | 'keep' nr 'splitat' nr 'h': help\n");
1207 printf(" 'r': residue 'res' nr 'chain' char\n");
1208 printf(" \"name\": group 'case': case %s 'q': save and quit\n",
1209 bCase ? "insensitive" : "sensitive ");
1210 printf(" 'ri': residue index\n");
1215 if (NULL == fgets(inp_string, STRLEN, stdin))
1217 gmx_fatal(FARGS, "Error reading user input");
1219 inp_string[strlen(inp_string)-1] = 0;
1221 string = inp_string;
1222 while (string[0] == ' ')
1229 if (string[0] == 'h')
1231 printf(" nr : selects an index group by number or quoted string.\n");
1232 printf(" The string is first matched against the whole group name,\n");
1233 printf(" then against the beginning and finally against an\n");
1234 printf(" arbitrary substring. A multiple match is an error.\n");
1236 printf(" 'a' nr1 [nr2 ...] : selects atoms, atom numbering starts at 1.\n");
1237 printf(" 'a' nr1 - nr2 : selects atoms in the range from nr1 to nr2.\n");
1238 printf(" 'a' name1[*] [name2[*] ...] : selects atoms by name(s), '?' matches any char,\n");
1239 printf(" wildcard '*' allowed at the end of a name.\n");
1240 printf(" 't' type1[*] [type2[*] ...] : as 'a', but for type, run input file required.\n");
1241 printf(" 'r' nr1[ic1] [nr2[ic2] ...] : selects residues by number and insertion code.\n");
1242 printf(" 'r' nr1 - nr2 : selects residues in the range from nr1 to nr2.\n");
1243 printf(" 'r' name1[*] [name2[*] ...] : as 'a', but for residue names.\n");
1244 printf(" 'ri' nr1 - nr2 : selects residue indices, 1-indexed, (as opposed to numbers) in the range from nr1 to nr2.\n");
1245 printf(" 'chain' ch1 [ch2 ...] : selects atoms by chain identifier(s),\n");
1246 printf(" not available with a .gro file as input.\n");
1247 printf(" ! : takes the complement of a group with respect to all\n");
1248 printf(" the atoms in the input file.\n");
1249 printf(" & | : AND and OR, can be placed between any of the options\n");
1250 printf(" above, the input is processed from left to right.\n");
1251 printf(" 'name' nr name : rename group nr to name.\n");
1252 printf(" 'del' nr1 [- nr2] : deletes one group or groups in the range from nr1 to nr2.\n");
1253 printf(" 'keep' nr : deletes all groups except nr.\n");
1254 printf(" 'case' : make all name compares case (in)sensitive.\n");
1255 printf(" 'splitch' nr : split group into chains using CA distances.\n");
1256 printf(" 'splitres' nr : split group into residues.\n");
1257 printf(" 'splitat' nr : split group into atoms.\n");
1258 printf(" 'res' nr : interpret numbers in group as residue numbers\n");
1259 printf(" Enter : list the currently defined groups and commands\n");
1260 printf(" 'l' : list the residues.\n");
1261 printf(" 'h' : show this help.\n");
1262 printf(" 'q' : save and quit.\n");
1264 printf(" Examples:\n");
1265 printf(" > 2 | 4 & r 3-5\n");
1266 printf(" selects all atoms from group 2 and 4 that have residue numbers 3, 4 or 5\n");
1267 printf(" > a C* & !a C CA\n");
1268 printf(" selects all atoms starting with 'C' but not the atoms 'C' and 'CA'\n");
1269 printf(" > \"protein\" & ! \"backb\"\n");
1270 printf(" selects all atoms that are in group 'protein' and not in group 'backbone'\n");
1273 printf("\npress Enter ");
1277 else if (strncmp(string, "del", 3) == 0)
1280 if (parse_int(&string, &sel_nr))
1282 while (string[0] == ' ')
1286 if (string[0] == '-')
1289 parse_int(&string, &sel_nr2);
1295 while (string[0] == ' ')
1299 if (string[0] == '\0')
1301 remove_group(sel_nr, sel_nr2, block, gn);
1305 printf("\nSyntax error: \"%s\"\n", string);
1309 else if (strncmp(string, "keep", 4) == 0)
1312 if (parse_int(&string, &sel_nr))
1314 remove_group(sel_nr+1, block->nr-1, block, gn);
1315 remove_group(0, sel_nr-1, block, gn);
1318 else if (strncmp(string, "name", 4) == 0)
1321 if (parse_int(&string, &sel_nr))
1323 if ((sel_nr >= 0) && (sel_nr < block->nr))
1325 sscanf(string, "%s", gname);
1326 sfree((*gn)[sel_nr]);
1327 (*gn)[sel_nr] = strdup(gname);
1331 else if (strncmp(string, "case", 4) == 0)
1334 printf("Switched to case %s\n", bCase ? "sensitive" : "insensitive");
1336 else if (string[0] == 'v')
1338 bVerbose = !bVerbose;
1339 printf("Turned verbose %s\n", bVerbose ? "on" : "off");
1341 else if (string[0] == 'l')
1343 if (check_have_atoms(atoms, ostring) )
1345 list_residues(atoms);
1348 else if (strncmp(string, "splitch", 7) == 0)
1351 if (check_have_atoms(atoms, ostring) &&
1352 parse_int(&string, &sel_nr) &&
1353 (sel_nr >= 0) && (sel_nr < block->nr))
1355 split_chain(atoms, x, sel_nr, block, gn);
1358 else if (strncmp(string, "splitres", 8) == 0)
1361 if (check_have_atoms(atoms, ostring) &&
1362 parse_int(&string, &sel_nr) &&
1363 (sel_nr >= 0) && (sel_nr < block->nr))
1365 split_group(atoms, sel_nr, block, gn, FALSE);
1368 else if (strncmp(string, "splitat", 7) == 0)
1371 if (check_have_atoms(atoms, ostring) &&
1372 parse_int(&string, &sel_nr) &&
1373 (sel_nr >= 0) && (sel_nr < block->nr))
1375 split_group(atoms, sel_nr, block, gn, TRUE);
1378 else if (string[0] == '\0')
1382 else if (string[0] != 'q')
1386 if (parse_entry(&string, natoms, atoms, block, gn, &nr, index, gname))
1390 while (string[0] == ' ')
1397 if (string[0] == '&')
1401 else if (string[0] == '|')
1410 for (i = 0; i < nr; i++)
1412 index1[i] = index[i];
1414 strcpy(gname1, gname);
1415 if (parse_entry(&string, natoms, atoms, block, gn, &nr2, index2, gname2))
1419 or_groups(nr1, index1, nr2, index2, &nr, index);
1420 sprintf(gname, "%s_%s", gname1, gname2);
1424 and_groups(nr1, index1, nr2, index2, &nr, index);
1425 sprintf(gname, "%s_&_%s", gname1, gname2);
1430 while (bAnd || bOr);
1432 while (string[0] == ' ')
1438 printf("\nSyntax error: \"%s\"\n", string);
1442 copy2block(nr, index, block);
1443 srenew(*gn, block->nr);
1444 newgroup = block->nr-1;
1445 (*gn)[newgroup] = strdup(gname);
1449 printf("Group is empty\n");
1453 while (string[0] != 'q');
1460 static int block2natoms(t_blocka *block)
1465 for (i = 0; i < block->nra; i++)
1467 natoms = max(natoms, block->a[i]);
1474 void merge_blocks(t_blocka *dest, t_blocka *source)
1478 /* count groups, srenew and fill */
1481 dest->nr += source->nr;
1482 srenew(dest->index, dest->nr+1);
1483 for (i = 0; i < source->nr; i++)
1485 dest->index[i0+i] = nra0 + source->index[i];
1487 /* count atoms, srenew and fill */
1488 dest->nra += source->nra;
1489 srenew(dest->a, dest->nra);
1490 for (i = 0; i < source->nra; i++)
1492 dest->a[nra0+i] = source->a[i];
1495 /* terminate list */
1496 dest->index[dest->nr] = dest->nra;
1500 int gmx_make_ndx(int argc, char *argv[])
1502 const char *desc[] = {
1503 "Index groups are necessary for almost every gromacs program.",
1504 "All these programs can generate default index groups. You ONLY",
1505 "have to use [TT]make_ndx[tt] when you need SPECIAL index groups.",
1506 "There is a default index group for the whole system, 9 default",
1507 "index groups for proteins, and a default index group",
1508 "is generated for every other residue name.[PAR]",
1509 "When no index file is supplied, also [TT]make_ndx[tt] will generate the",
1511 "With the index editor you can select on atom, residue and chain names",
1513 "When a run input file is supplied you can also select on atom type.",
1514 "You can use NOT, AND and OR, you can split groups",
1515 "into chains, residues or atoms. You can delete and rename groups.[PAR]",
1516 "The atom numbering in the editor and the index file starts at 1."
1519 static int natoms = 0;
1520 static gmx_bool bVerbose = FALSE;
1522 { "-natoms", FALSE, etINT, {&natoms},
1523 "set number of atoms (default: read from coordinate or index file)" },
1524 { "-verbose", FALSE, etBOOL, {&bVerbose},
1525 "HIDDENVerbose output" }
1527 #define NPA asize(pa)
1532 const char *stxfile;
1534 const char *ndxoutfile;
1541 t_blocka *block, *block2;
1542 char **gnames, **gnames2;
1544 { efSTX, "-f", NULL, ffOPTRD },
1545 { efNDX, "-n", NULL, ffOPTRDMULT },
1546 { efNDX, "-o", NULL, ffWRITE }
1548 #define NFILE asize(fnm)
1550 parse_common_args(&argc, argv, 0, NFILE, fnm, NPA, pa, asize(desc), desc,
1553 stxfile = ftp2fn_null(efSTX, NFILE, fnm);
1554 if (opt2bSet("-n", NFILE, fnm))
1556 nndxin = opt2fns(&ndxinfiles, "-n", NFILE, fnm);
1562 ndxoutfile = opt2fn("-o", NFILE, fnm);
1563 bNatoms = opt2parg_bSet("-natoms", NPA, pa);
1565 if (!stxfile && !nndxin)
1567 gmx_fatal(FARGS, "No input files (structure or index)");
1573 get_stx_coordnum(stxfile, &(atoms->nr));
1574 init_t_atoms(atoms, atoms->nr, TRUE);
1577 fprintf(stderr, "\nReading structure file\n");
1578 read_stx_conf(stxfile, title, atoms, x, v, &ePBC, box);
1588 /* read input file(s) */
1589 block = new_blocka();
1591 printf("Going to read %d old index file(s)\n", nndxin);
1594 for (i = 0; i < nndxin; i++)
1596 block2 = init_index(ndxinfiles[i], &gnames2);
1597 srenew(gnames, block->nr+block2->nr);
1598 for (j = 0; j < block2->nr; j++)
1600 gnames[block->nr+j] = gnames2[j];
1603 merge_blocks(block, block2);
1605 sfree(block2->index);
1606 /* done_block(block2); */
1613 analyse(atoms, block, &gnames, FALSE, TRUE);
1618 natoms = block2natoms(block);
1619 printf("Counted atom numbers up to %d in index file\n", natoms);
1622 edit_index(natoms, atoms, x, block, &gnames, bVerbose);
1624 write_index(ndxoutfile, block, gnames);