void print_resall(FILE *out, int nrtp, t_restp rtp[], gpp_atomtype_t atype) { int i,bt; if (nrtp == 0) { return; } print_resall_header(out, rtp); for(i=0; i<nrtp; i++) { if (rtp[i].natom > 0) { print_resatoms(out,atype,&rtp[i]); for(bt=0; bt<ebtsNR; bt++) print_resbondeds(out,bt,&rtp[i]); } } }
void read_resall(char *rrdb, int *nrtpptr, t_restp **rtp, gpp_atomtype_t atype, t_symtab *tab, gmx_bool bAllowOverrideRTP) { FILE *in; char filebase[STRLEN], line[STRLEN], header[STRLEN]; int i, nrtp, maxrtp, bt, nparam; int dum1, dum2, dum3; t_restp *rrtp, *header_settings; gmx_bool bNextResidue, bError; int firstrtp; fflib_filename_base(rrdb, filebase, STRLEN); in = fflib_open(rrdb); if (debug) { fprintf(debug, "%9s %5s", "Residue", "atoms"); for (i = 0; i < ebtsNR; i++) { fprintf(debug, " %10s", btsNames[i]); } fprintf(debug, "\n"); } snew(header_settings, 1); /* these bonded parameters will overwritten be when * * there is a [ bondedtypes ] entry in the .rtp file */ header_settings->rb[ebtsBONDS].type = 1; /* normal bonds */ header_settings->rb[ebtsANGLES].type = 1; /* normal angles */ header_settings->rb[ebtsPDIHS].type = 1; /* normal dihedrals */ header_settings->rb[ebtsIDIHS].type = 2; /* normal impropers */ header_settings->rb[ebtsEXCLS].type = 1; /* normal exclusions */ header_settings->rb[ebtsCMAP].type = 1; /* normal cmap torsions */ header_settings->bKeepAllGeneratedDihedrals = FALSE; header_settings->nrexcl = 3; header_settings->bGenerateHH14Interactions = TRUE; header_settings->bRemoveDihedralIfWithImproper = TRUE; /* Column 5 & 6 aren't really bonded types, but we include * them here to avoid introducing a new section: * Column 5 : This controls the generation of dihedrals from the bonding. * All possible dihedrals are generated automatically. A value of * 1 here means that all these are retained. A value of * 0 here requires generated dihedrals be removed if * * there are any dihedrals on the same central atoms * specified in the residue topology, or * * there are other identical generated dihedrals * sharing the same central atoms, or * * there are other generated dihedrals sharing the * same central bond that have fewer hydrogen atoms * Column 6: Number of bonded neighbors to exclude. * Column 7: Generate 1,4 interactions between two hydrogen atoms * Column 8: Remove proper dihedrals if centered on the same bond * as an improper dihedral */ get_a_line(in, line, STRLEN); if (!get_header(line, header)) { gmx_fatal(FARGS, "in .rtp file at line:\n%s\n", line); } if (gmx_strncasecmp("bondedtypes", header, 5) == 0) { get_a_line(in, line, STRLEN); if ((nparam = sscanf(line, "%d %d %d %d %d %d %d %d", &header_settings->rb[ebtsBONDS].type, &header_settings->rb[ebtsANGLES].type, &header_settings->rb[ebtsPDIHS].type, &header_settings->rb[ebtsIDIHS].type, &dum1, &header_settings->nrexcl, &dum2, &dum3)) < 4) { gmx_fatal(FARGS, "need 4 to 8 parameters in the header of .rtp file %s at line:\n%s\n", rrdb, line); } header_settings->bKeepAllGeneratedDihedrals = (dum1 != 0); header_settings->bGenerateHH14Interactions = (dum2 != 0); header_settings->bRemoveDihedralIfWithImproper = (dum3 != 0); get_a_line(in, line, STRLEN); if (nparam < 5) { fprintf(stderr, "Using default: not generating all possible dihedrals\n"); header_settings->bKeepAllGeneratedDihedrals = FALSE; } if (nparam < 6) { fprintf(stderr, "Using default: excluding 3 bonded neighbors\n"); header_settings->nrexcl = 3; } if (nparam < 7) { fprintf(stderr, "Using default: generating 1,4 H--H interactions\n"); header_settings->bGenerateHH14Interactions = TRUE; } if (nparam < 8) { fprintf(stderr, "Using default: removing proper dihedrals found on the same bond as a proper dihedral\n"); header_settings->bRemoveDihedralIfWithImproper = TRUE; } } else { fprintf(stderr, "Reading .rtp file without '[ bondedtypes ]' directive,\n" "Will proceed as if the entry was:\n"); print_resall_header(stderr, header_settings); } /* We don't know the current size of rrtp, but simply realloc immediately */ nrtp = *nrtpptr; rrtp = *rtp; maxrtp = nrtp; while (!feof(in)) { if (nrtp >= maxrtp) { maxrtp += 100; srenew(rrtp, maxrtp); } /* Initialise rtp entry structure */ rrtp[nrtp] = *header_settings; if (!get_header(line, header)) { gmx_fatal(FARGS, "in .rtp file at line:\n%s\n", line); } rrtp[nrtp].resname = gmx_strdup(header); rrtp[nrtp].filebase = gmx_strdup(filebase); get_a_line(in, line, STRLEN); bError = FALSE; bNextResidue = FALSE; do { if (!get_header(line, header)) { bError = TRUE; } else { bt = get_bt(header); if (bt != NOTSET) { /* header is an bonded directive */ bError = !read_bondeds(bt, in, line, &rrtp[nrtp]); } else if (gmx_strncasecmp("atoms", header, 5) == 0) { /* header is the atoms directive */ bError = !read_atoms(in, line, &(rrtp[nrtp]), tab, atype); } else { /* else header must be a residue name */ bNextResidue = TRUE; } } if (bError) { gmx_fatal(FARGS, "in .rtp file in residue %s at line:\n%s\n", rrtp[nrtp].resname, line); } } while (!feof(in) && !bNextResidue); if (rrtp[nrtp].natom == 0) { gmx_fatal(FARGS, "No atoms found in .rtp file in residue %s\n", rrtp[nrtp].resname); } if (debug) { fprintf(debug, "%3d %5s %5d", nrtp+1, rrtp[nrtp].resname, rrtp[nrtp].natom); for (i = 0; i < ebtsNR; i++) { fprintf(debug, " %10d", rrtp[nrtp].rb[i].nb); } fprintf(debug, "\n"); } firstrtp = -1; for (i = 0; i < nrtp; i++) { if (gmx_strcasecmp(rrtp[i].resname, rrtp[nrtp].resname) == 0) { firstrtp = i; } } if (firstrtp == -1) { nrtp++; fprintf(stderr, "\rResidue %d", nrtp); fflush(stderr); } else { if (firstrtp >= *nrtpptr) { gmx_fatal(FARGS, "Found a second entry for '%s' in '%s'", rrtp[nrtp].resname, rrdb); } if (bAllowOverrideRTP) { fprintf(stderr, "\n"); fprintf(stderr, "Found another rtp entry for '%s' in '%s', ignoring this entry and keeping the one from '%s.rtp'\n", rrtp[nrtp].resname, rrdb, rrtp[firstrtp].filebase); /* We should free all the data for this entry. * The current code gives a lot of dangling pointers. */ clear_t_restp(&rrtp[nrtp]); } else { gmx_fatal(FARGS, "Found rtp entries for '%s' in both '%s' and '%s'. If you want the first definition to override the second one, set the -rtpo option of pdb2gmx.", rrtp[nrtp].resname, rrtp[firstrtp].filebase, rrdb); } } } gmx_ffclose(in); /* give back unused memory */ srenew(rrtp, nrtp); fprintf(stderr, "\nSorting it all out...\n"); qsort(rrtp, nrtp, (size_t)sizeof(rrtp[0]), comprtp); check_rtp(nrtp, rrtp, rrdb); *nrtpptr = nrtp; *rtp = rrtp; }