1
|
|
2
|
|
3
|
|
4
|
|
5
|
|
6
|
|
7
|
|
8
|
|
9
|
|
10
|
|
11
|
|
12
|
|
13
|
|
14
|
|
15
|
|
16
|
|
17
|
|
18
|
|
19
|
|
20
|
|
21
|
|
22
|
|
23
|
|
24
|
|
25
|
|
26
|
|
27
|
|
28
|
|
29
|
|
30
|
|
31
|
|
32
|
|
33
|
|
34
|
|
35
|
#ifdef HAVE_CONFIG_H
|
36
|
#include <config.h>
|
37
|
#endif
|
38
|
|
39
|
#include <ctype.h>
|
40
|
|
41
|
#include "sysstuff.h"
|
42
|
#include "string2.h"
|
43
|
#include "vec.h"
|
44
|
#include "smalloc.h"
|
45
|
#include "typedefs.h"
|
46
|
#include "symtab.h"
|
47
|
#include "pdbio.h"
|
48
|
#include "vec.h"
|
49
|
#include "copyrite.h"
|
50
|
#include "futil.h"
|
51
|
#include "atomprop.h"
|
52
|
#include "physics.h"
|
53
|
#include "pbc.h"
|
54
|
#include "gmxfio.h"
|
55
|
|
56
|
typedef struct {
|
57
|
int ai,aj;
|
58
|
} gmx_conection_t;
|
59
|
|
60
|
typedef struct gmx_conect_t {
|
61
|
int nconect;
|
62
|
bool bSorted;
|
63
|
gmx_conection_t *conect;
|
64
|
} gmx_conect_t;
|
65
|
|
66
|
static const char *pdbtp[epdbNR]={
|
67
|
"ATOM ","HETATM", "ANISOU", "CRYST1",
|
68
|
"COMPND", "MODEL", "ENDMDL", "TER", "HEADER", "TITLE", "REMARK",
|
69
|
"CONECT"
|
70
|
};
|
71
|
|
72
|
|
73
|
|
74
|
|
75
|
static bool bWideFormat=FALSE;
|
76
|
#define REMARK_SIM_BOX "REMARK THIS IS A SIMULATION BOX"
|
77
|
|
78
|
void set_pdb_wide_format(bool bSet)
|
79
|
{
|
80
|
bWideFormat = bSet;
|
81
|
}
|
82
|
|
83
|
static void xlate_atomname_pdb2gmx(char *name)
|
84
|
{
|
85
|
int i,length;
|
86
|
char temp;
|
87
|
|
88
|
length=strlen(name);
|
89
|
if (length>3 && isdigit(name[0])) {
|
90
|
temp=name[0];
|
91
|
for(i=1; i<length; i++)
|
92
|
name[i-1]=name[i];
|
93
|
name[length-1]=temp;
|
94
|
}
|
95
|
}
|
96
|
|
97
|
static void xlate_atomname_gmx2pdb(char *name)
|
98
|
{
|
99
|
int i,length;
|
100
|
char temp;
|
101
|
|
102
|
length=strlen(name);
|
103
|
if (length>3 && isdigit(name[length-1])) {
|
104
|
temp=name[length-1];
|
105
|
for(i=length-1; i>0; --i)
|
106
|
name[i]=name[i-1];
|
107
|
name[0]=temp;
|
108
|
}
|
109
|
}
|
110
|
|
111
|
|
112
|
void gmx_write_pdb_box(FILE *out,int ePBC,matrix box)
|
113
|
{
|
114
|
real alpha,beta,gamma;
|
115
|
|
116
|
if (ePBC == -1)
|
117
|
ePBC = guess_ePBC(box);
|
118
|
|
119
|
if (ePBC == epbcNONE)
|
120
|
return;
|
121
|
|
122
|
if (norm2(box[YY])*norm2(box[ZZ])!=0)
|
123
|
alpha = RAD2DEG*acos(cos_angle_no_table(box[YY],box[ZZ]));
|
124
|
else
|
125
|
alpha = 90;
|
126
|
if (norm2(box[XX])*norm2(box[ZZ])!=0)
|
127
|
beta = RAD2DEG*acos(cos_angle_no_table(box[XX],box[ZZ]));
|
128
|
else
|
129
|
beta = 90;
|
130
|
if (norm2(box[XX])*norm2(box[YY])!=0)
|
131
|
gamma = RAD2DEG*acos(cos_angle_no_table(box[XX],box[YY]));
|
132
|
else
|
133
|
gamma = 90;
|
134
|
fprintf(out,"REMARK THIS IS A SIMULATION BOX\n");
|
135
|
if (ePBC != epbcSCREW) {
|
136
|
fprintf(out,"CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s%4d\n",
|
137
|
10*norm(box[XX]),10*norm(box[YY]),10*norm(box[ZZ]),
|
138
|
alpha,beta,gamma,"P 1",1);
|
139
|
} else {
|
140
|
|
141
|
fprintf(out,"CRYST1%9.3f%9.3f%9.3f%7.2f%7.2f%7.2f %-11s%4d\n",
|
142
|
20*norm(box[XX]),10*norm(box[YY]),10*norm(box[ZZ]),
|
143
|
alpha,beta,gamma,"P 21 1 1",1);
|
144
|
|
145
|
}
|
146
|
}
|
147
|
|
148
|
static void read_cryst1(char *line,int *ePBC,matrix box)
|
149
|
{
|
150
|
#define SG_SIZE 11
|
151
|
char sa[12],sb[12],sc[12],sg[SG_SIZE+1],ident;
|
152
|
double fa,fb,fc,alpha,beta,gamma,cosa,cosb,cosg,sing;
|
153
|
int syma,symb,symc;
|
154
|
int ePBC_file;
|
155
|
|
156
|
sscanf(line,"%*s%s%s%s%lf%lf%lf",sa,sb,sc,&alpha,&beta,&gamma);
|
157
|
|
158
|
ePBC_file = -1;
|
159
|
if (strlen(line) >= 55) {
|
160
|
strncpy(sg,line+55,SG_SIZE);
|
161
|
sg[SG_SIZE] = '\0';
|
162
|
ident = ' ';
|
163
|
syma = 0;
|
164
|
symb = 0;
|
165
|
symc = 0;
|
166
|
sscanf(sg,"%c %d %d %d",&ident,&syma,&symb,&symc);
|
167
|
if (ident == 'P' && syma == 1 && symb <= 1 && symc <= 1) {
|
168
|
fc = strtod(sc,NULL)*0.1;
|
169
|
ePBC_file = (fc > 0 ? epbcXYZ : epbcXY);
|
170
|
}
|
171
|
if (ident == 'P' && syma == 21 && symb == 1 && symc == 1) {
|
172
|
ePBC_file = epbcSCREW;
|
173
|
}
|
174
|
}
|
175
|
if (ePBC) {
|
176
|
*ePBC = ePBC_file;
|
177
|
}
|
178
|
|
179
|
if (box) {
|
180
|
fa = strtod(sa,NULL)*0.1;
|
181
|
fb = strtod(sb,NULL)*0.1;
|
182
|
fc = strtod(sc,NULL)*0.1;
|
183
|
if (ePBC_file == epbcSCREW) {
|
184
|
fa *= 0.5;
|
185
|
}
|
186
|
clear_mat(box);
|
187
|
box[XX][XX] = fa;
|
188
|
if ((alpha!=90.0) || (beta!=90.0) || (gamma!=90.0)) {
|
189
|
if (alpha != 90.0) {
|
190
|
cosa = cos(alpha*DEG2RAD);
|
191
|
} else {
|
192
|
cosa = 0;
|
193
|
}
|
194
|
if (beta != 90.0) {
|
195
|
cosb = cos(beta*DEG2RAD);
|
196
|
} else {
|
197
|
cosb = 0;
|
198
|
}
|
199
|
if (gamma != 90.0) {
|
200
|
cosg = cos(gamma*DEG2RAD);
|
201
|
sing = sin(gamma*DEG2RAD);
|
202
|
} else {
|
203
|
cosg = 0;
|
204
|
sing = 1;
|
205
|
}
|
206
|
box[YY][XX] = fb*cosg;
|
207
|
box[YY][YY] = fb*sing;
|
208
|
box[ZZ][XX] = fc*cosb;
|
209
|
box[ZZ][YY] = fc*(cosa - cosb*cosg)/sing;
|
210
|
box[ZZ][ZZ] = sqrt(fc*fc
|
211
|
- box[ZZ][XX]*box[ZZ][XX] - box[ZZ][YY]*box[ZZ][YY]);
|
212
|
} else {
|
213
|
box[YY][YY] = fb;
|
214
|
box[ZZ][ZZ] = fc;
|
215
|
}
|
216
|
}
|
217
|
}
|
218
|
|
219
|
void write_pdbfile_indexed(FILE *out,const char *title,
|
220
|
t_atoms *atoms,rvec x[],
|
221
|
int ePBC,matrix box,char chainid,
|
222
|
int model_nr, atom_id nindex, atom_id index[],
|
223
|
gmx_conect conect, bool bTerSepChains)
|
224
|
{
|
225
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
226
|
char resnm[6],nm[6],pdbform[128],pukestring[100];
|
227
|
atom_id i,ii;
|
228
|
int resind,resnr,type;
|
229
|
unsigned char resic,ch;
|
230
|
real occup,bfac;
|
231
|
bool bOccup;
|
232
|
int nlongname=0;
|
233
|
int chainnum,lastchainnum;
|
234
|
int lastresind,lastchainresind;
|
235
|
gmx_residuetype_t rt;
|
236
|
const char *p_restype;
|
237
|
const char *p_lastrestype;
|
238
|
|
239
|
gmx_residuetype_init(&rt);
|
240
|
|
241
|
bromacs(pukestring,99);
|
242
|
fprintf(out,"TITLE %s\n",(title && title[0])?title:pukestring);
|
243
|
if (bWideFormat) {
|
244
|
fprintf(out,"REMARK This file does not adhere to the PDB standard\n");
|
245
|
fprintf(out,"REMARK As a result of, some programs may not like it\n");
|
246
|
}
|
247
|
if (box && ( norm2(box[XX]) || norm2(box[YY]) || norm2(box[ZZ]) ) ) {
|
248
|
gmx_write_pdb_box(out,ePBC,box);
|
249
|
}
|
250
|
if (atoms->pdbinfo) {
|
251
|
|
252
|
|
253
|
|
254
|
bOccup = TRUE;
|
255
|
for (ii=0; (ii<nindex) && bOccup; ii++) {
|
256
|
i = index[ii];
|
257
|
bOccup = bOccup && (atoms->pdbinfo[i].occup == 0.0);
|
258
|
}
|
259
|
}
|
260
|
else
|
261
|
bOccup = FALSE;
|
262
|
|
263
|
fprintf(out,"MODEL %8d\n",model_nr>=0 ? model_nr : 1);
|
264
|
|
265
|
lastchainresind = -1;
|
266
|
lastchainnum = -1;
|
267
|
resind = -1;
|
268
|
p_restype = NULL;
|
269
|
|
270
|
for (ii=0; ii<nindex; ii++) {
|
271
|
i=index[ii];
|
272
|
lastresind = resind;
|
273
|
resind = atoms->atom[i].resind;
|
274
|
chainnum = atoms->resinfo[resind].chainnum;
|
275
|
p_lastrestype = p_restype;
|
276
|
gmx_residuetype_get_type(rt,*atoms->resinfo[resind].name,&p_restype);
|
277
|
|
278
|
|
279
|
if( bTerSepChains && ii>0 && chainnum != lastchainnum)
|
280
|
{
|
281
|
|
282
|
if(gmx_residuetype_is_protein(rt,p_lastrestype) || gmx_residuetype_is_dna(rt,p_lastrestype) || gmx_residuetype_is_rna(rt,p_lastrestype))
|
283
|
{
|
284
|
fprintf(out,"TER\n");
|
285
|
}
|
286
|
lastchainnum = chainnum;
|
287
|
lastchainresind = lastresind;
|
288
|
}
|
289
|
|
290
|
strncpy(resnm,*atoms->resinfo[resind].name,sizeof(resnm)-1);
|
291
|
strncpy(nm,*atoms->atomname[i],sizeof(nm)-1);
|
292
|
|
293
|
xlate_atomname_gmx2pdb(nm);
|
294
|
resnr = atoms->resinfo[resind].nr;
|
295
|
resic = atoms->resinfo[resind].ic;
|
296
|
if (chainid!=' ')
|
297
|
{
|
298
|
ch = chainid;
|
299
|
}
|
300
|
else
|
301
|
{
|
302
|
ch = atoms->resinfo[resind].chainid;
|
303
|
|
304
|
if (ch == 0)
|
305
|
{
|
306
|
ch = ' ';
|
307
|
}
|
308
|
}
|
309
|
if (resnr>=10000)
|
310
|
resnr = resnr % 10000;
|
311
|
if (atoms->pdbinfo) {
|
312
|
type = atoms->pdbinfo[i].type;
|
313
|
occup = bOccup ? 1.0 : atoms->pdbinfo[i].occup;
|
314
|
bfac = atoms->pdbinfo[i].bfac;
|
315
|
}
|
316
|
else {
|
317
|
type = 0;
|
318
|
occup = 1.0;
|
319
|
bfac = 0.0;
|
320
|
}
|
321
|
if (bWideFormat)
|
322
|
strcpy(pdbform,
|
323
|
"%-6s%5u %-4.4s %3.3s %c%4d%c %10.5f%10.5f%10.5f%8.4f%8.4f %2s\n");
|
324
|
else {
|
325
|
|
326
|
if ((strlen(nm)<4) && (gmx_strcasecmp(nm,atoms->atom[i].elem) != 0))
|
327
|
strcpy(pdbform,pdbformat);
|
328
|
else {
|
329
|
strcpy(pdbform,pdbformat4);
|
330
|
if (strlen(nm) > 4) {
|
331
|
int maxwln=20;
|
332
|
if (nlongname < maxwln) {
|
333
|
fprintf(stderr,"WARNING: Writing out atom name (%s) longer than 4 characters to .pdb file\n",nm);
|
334
|
} else if (nlongname == maxwln) {
|
335
|
fprintf(stderr,"WARNING: More than %d long atom names, will not write more warnings\n",maxwln);
|
336
|
}
|
337
|
nlongname++;
|
338
|
}
|
339
|
}
|
340
|
strcat(pdbform,"%6.2f%6.2f %2s\n");
|
341
|
}
|
342
|
fprintf(out,pdbform,pdbtp[type],(i+1)%100000,nm,resnm,ch,resnr,
|
343
|
(resic == '\0') ? ' ' : resic,
|
344
|
10*x[i][XX],10*x[i][YY],10*x[i][ZZ],occup,bfac,atoms->atom[i].elem);
|
345
|
if (atoms->pdbinfo && atoms->pdbinfo[i].bAnisotropic) {
|
346
|
fprintf(out,"ANISOU%5u %-4.4s%3.3s %c%4d%c %7d%7d%7d%7d%7d%7d\n",
|
347
|
(i+1)%100000,nm,resnm,ch,resnr,
|
348
|
(resic == '\0') ? ' ' : resic,
|
349
|
atoms->pdbinfo[i].uij[0],atoms->pdbinfo[i].uij[1],
|
350
|
atoms->pdbinfo[i].uij[2],atoms->pdbinfo[i].uij[3],
|
351
|
atoms->pdbinfo[i].uij[4],atoms->pdbinfo[i].uij[5]);
|
352
|
}
|
353
|
}
|
354
|
|
355
|
fprintf(out,"TER\n");
|
356
|
fprintf(out,"ENDMDL\n");
|
357
|
|
358
|
if (NULL != gc) {
|
359
|
|
360
|
for(i=0; (i<gc->nconect); i++) {
|
361
|
fprintf(out,"CONECT%5d%5d\n",gc->conect[i].ai+1,gc->conect[i].aj+1);
|
362
|
}
|
363
|
}
|
364
|
}
|
365
|
|
366
|
void write_pdbfile(FILE *out,const char *title, t_atoms *atoms,rvec x[],
|
367
|
int ePBC,matrix box,char chainid,int model_nr,gmx_conect conect,bool bTerSepChains)
|
368
|
{
|
369
|
atom_id i,*index;
|
370
|
|
371
|
snew(index,atoms->nr);
|
372
|
for(i=0; i<atoms->nr; i++)
|
373
|
index[i]=i;
|
374
|
write_pdbfile_indexed(out,title,atoms,x,ePBC,box,chainid,model_nr,
|
375
|
atoms->nr,index,conect,bTerSepChains);
|
376
|
sfree(index);
|
377
|
}
|
378
|
|
379
|
int line2type(char *line)
|
380
|
{
|
381
|
int k;
|
382
|
char type[8];
|
383
|
|
384
|
for(k=0; (k<6); k++)
|
385
|
type[k]=line[k];
|
386
|
type[k]='\0';
|
387
|
|
388
|
for(k=0; (k<epdbNR); k++)
|
389
|
if (strncmp(type,pdbtp[k],strlen(pdbtp[k])) == 0)
|
390
|
break;
|
391
|
|
392
|
return k;
|
393
|
}
|
394
|
|
395
|
static void read_anisou(char line[],int natom,t_atoms *atoms)
|
396
|
{
|
397
|
int i,j,k,atomnr;
|
398
|
char nc='\0';
|
399
|
char anr[12],anm[12];
|
400
|
|
401
|
|
402
|
j=6;
|
403
|
for(k=0; (k<5); k++,j++) anr[k]=line[j];
|
404
|
anr[k]=nc;
|
405
|
j++;
|
406
|
for(k=0; (k<4); k++,j++) anm[k]=line[j];
|
407
|
anm[k]=nc;
|
408
|
j++;
|
409
|
|
410
|
|
411
|
trim(anm);
|
412
|
|
413
|
|
414
|
atomnr = strtol(anr, NULL, 10);
|
415
|
for(i=natom-1; (i>=0); i--)
|
416
|
if ((strcmp(anm,*(atoms->atomname[i])) == 0) &&
|
417
|
(atomnr == atoms->pdbinfo[i].atomnr))
|
418
|
break;
|
419
|
if (i < 0)
|
420
|
fprintf(stderr,"Skipping ANISOU record (atom %s %d not found)\n",
|
421
|
anm,atomnr);
|
422
|
else {
|
423
|
if (sscanf(line+29,"%d%d%d%d%d%d",
|
424
|
&atoms->pdbinfo[i].uij[U11],&atoms->pdbinfo[i].uij[U22],
|
425
|
&atoms->pdbinfo[i].uij[U33],&atoms->pdbinfo[i].uij[U12],
|
426
|
&atoms->pdbinfo[i].uij[U13],&atoms->pdbinfo[i].uij[U23])
|
427
|
== 6) {
|
428
|
atoms->pdbinfo[i].bAnisotropic = TRUE;
|
429
|
}
|
430
|
else {
|
431
|
fprintf(stderr,"Invalid ANISOU record for atom %d\n",i);
|
432
|
atoms->pdbinfo[i].bAnisotropic = FALSE;
|
433
|
}
|
434
|
}
|
435
|
}
|
436
|
|
437
|
void get_pdb_atomnumber(t_atoms *atoms,gmx_atomprop_t aps)
|
438
|
{
|
439
|
int i,atomnumber;
|
440
|
size_t k;
|
441
|
char anm[6],anm_copy[6];
|
442
|
char nc='\0';
|
443
|
real eval;
|
444
|
|
445
|
if (!atoms->pdbinfo) {
|
446
|
gmx_incons("Trying to deduce atomnumbers when no pdb information is present");
|
447
|
}
|
448
|
for(i=0; (i<atoms->nr); i++) {
|
449
|
strcpy(anm,atoms->pdbinfo[i].atomnm);
|
450
|
strcpy(anm_copy,atoms->pdbinfo[i].atomnm);
|
451
|
atomnumber = NOTSET;
|
452
|
if (anm[0] != ' ') {
|
453
|
anm_copy[2] = nc;
|
454
|
if (gmx_atomprop_query(aps,epropElement,"???",anm_copy,&eval))
|
455
|
atomnumber = gmx_nint(eval);
|
456
|
else {
|
457
|
anm_copy[1] = nc;
|
458
|
if (gmx_atomprop_query(aps,epropElement,"???",anm_copy,&eval))
|
459
|
atomnumber = gmx_nint(eval);
|
460
|
}
|
461
|
}
|
462
|
if (atomnumber == NOTSET) {
|
463
|
k=0;
|
464
|
while ((k < strlen(anm)) && (isspace(anm[k]) || isdigit(anm[k])))
|
465
|
k++;
|
466
|
anm_copy[0] = anm[k];
|
467
|
anm_copy[1] = nc;
|
468
|
if (gmx_atomprop_query(aps,epropElement,"???",anm_copy,&eval))
|
469
|
atomnumber = gmx_nint(eval);
|
470
|
}
|
471
|
atoms->atom[i].atomnumber = atomnumber;
|
472
|
sprintf(atoms->atom[i].elem,"%2s",gmx_atomprop_element(aps,atomnumber));
|
473
|
if (debug)
|
474
|
fprintf(debug,"Atomnumber for atom '%s' is %d\n",anm,atomnumber);
|
475
|
}
|
476
|
}
|
477
|
|
478
|
static int read_atom(t_symtab *symtab,
|
479
|
char line[],int type,int natom,
|
480
|
t_atoms *atoms,rvec x[],int chainnum,bool bChange)
|
481
|
{
|
482
|
t_atom *atomn;
|
483
|
int j,k;
|
484
|
char nc='\0';
|
485
|
char anr[12],anm[12],anm_copy[12],altloc,resnm[12],rnr[12];
|
486
|
char xc[12],yc[12],zc[12],occup[12],bfac[12];
|
487
|
unsigned char resic;
|
488
|
char chainid;
|
489
|
int resnr,atomnumber;
|
490
|
|
491
|
if (natom>=atoms->nr)
|
492
|
gmx_fatal(FARGS,"\nFound more atoms (%d) in pdb file than expected (%d)",
|
493
|
natom+1,atoms->nr);
|
494
|
|
495
|
|
496
|
j=6;
|
497
|
for(k=0; (k<5); k++,j++) anr[k]=line[j];
|
498
|
anr[k]=nc;
|
499
|
trim(anr);
|
500
|
j++;
|
501
|
for(k=0; (k<4); k++,j++) anm[k]=line[j];
|
502
|
anm[k]=nc;
|
503
|
strcpy(anm_copy,anm);
|
504
|
atomnumber = NOTSET;
|
505
|
trim(anm);
|
506
|
altloc=line[j];
|
507
|
j++;
|
508
|
for(k=0; (k<4); k++,j++)
|
509
|
resnm[k]=line[j];
|
510
|
resnm[k]=nc;
|
511
|
trim(resnm);
|
512
|
|
513
|
chainid = line[j];
|
514
|
j++;
|
515
|
|
516
|
for(k=0; (k<4); k++,j++) {
|
517
|
rnr[k] = line[j];
|
518
|
}
|
519
|
rnr[k] = nc;
|
520
|
trim(rnr);
|
521
|
resnr = strtol(rnr, NULL, 10);
|
522
|
resic = line[j];
|
523
|
j+=4;
|
524
|
|
525
|
|
526
|
for(k=0; (k<8); k++,j++) xc[k]=line[j];
|
527
|
xc[k]=nc;
|
528
|
for(k=0; (k<8); k++,j++) yc[k]=line[j];
|
529
|
yc[k]=nc;
|
530
|
for(k=0; (k<8); k++,j++) zc[k]=line[j];
|
531
|
zc[k]=nc;
|
532
|
|
533
|
|
534
|
for(k=0; (k<6); k++,j++) occup[k]=line[j];
|
535
|
occup[k]=nc;
|
536
|
|
537
|
|
538
|
for(k=0; (k<7); k++,j++) bfac[k]=line[j];
|
539
|
bfac[k]=nc;
|
540
|
|
541
|
if (atoms->atom) {
|
542
|
atomn=&(atoms->atom[natom]);
|
543
|
if ((natom==0) ||
|
544
|
atoms->resinfo[atoms->atom[natom-1].resind].nr != resnr ||
|
545
|
atoms->resinfo[atoms->atom[natom-1].resind].ic != resic ||
|
546
|
(strcmp(*atoms->resinfo[atoms->atom[natom-1].resind].name,resnm) != 0))
|
547
|
{
|
548
|
if (natom == 0) {
|
549
|
atomn->resind = 0;
|
550
|
} else {
|
551
|
atomn->resind = atoms->atom[natom-1].resind + 1;
|
552
|
}
|
553
|
atoms->nres = atomn->resind + 1;
|
554
|
t_atoms_set_resinfo(atoms,natom,symtab,resnm,resnr,resic,chainnum,chainid);
|
555
|
}
|
556
|
else
|
557
|
{
|
558
|
atomn->resind = atoms->atom[natom-1].resind;
|
559
|
}
|
560
|
if (bChange) {
|
561
|
xlate_atomname_pdb2gmx(anm);
|
562
|
}
|
563
|
atoms->atomname[natom]=put_symtab(symtab,anm);
|
564
|
atomn->m = 0.0;
|
565
|
atomn->q = 0.0;
|
566
|
atomn->atomnumber = atomnumber;
|
567
|
atomn->elem[0] = '\0';
|
568
|
}
|
569
|
x[natom][XX]=strtod(xc,NULL)*0.1;
|
570
|
x[natom][YY]=strtod(yc,NULL)*0.1;
|
571
|
x[natom][ZZ]=strtod(zc,NULL)*0.1;
|
572
|
if (atoms->pdbinfo) {
|
573
|
atoms->pdbinfo[natom].type=type;
|
574
|
atoms->pdbinfo[natom].atomnr=strtol(anr, NULL, 10);
|
575
|
atoms->pdbinfo[natom].altloc=altloc;
|
576
|
strcpy(atoms->pdbinfo[natom].atomnm,anm_copy);
|
577
|
atoms->pdbinfo[natom].bfac=strtod(bfac,NULL);
|
578
|
atoms->pdbinfo[natom].occup=strtod(occup,NULL);
|
579
|
}
|
580
|
natom++;
|
581
|
|
582
|
return natom;
|
583
|
}
|
584
|
|
585
|
bool is_hydrogen(const char *nm)
|
586
|
{
|
587
|
char buf[30];
|
588
|
|
589
|
strcpy(buf,nm);
|
590
|
trim(buf);
|
591
|
|
592
|
if (buf[0] == 'H')
|
593
|
return TRUE;
|
594
|
else if ((isdigit(buf[0])) && (buf[1] == 'H'))
|
595
|
return TRUE;
|
596
|
return FALSE;
|
597
|
}
|
598
|
|
599
|
bool is_dummymass(const char *nm)
|
600
|
{
|
601
|
char buf[30];
|
602
|
|
603
|
strcpy(buf,nm);
|
604
|
trim(buf);
|
605
|
|
606
|
if ((buf[0] == 'M') && isdigit(buf[strlen(buf)-1]))
|
607
|
return TRUE;
|
608
|
|
609
|
return FALSE;
|
610
|
}
|
611
|
|
612
|
static void gmx_conect_addline(gmx_conect_t *con,char *line)
|
613
|
{
|
614
|
int n,ai,aj;
|
615
|
char format[32],form2[32];
|
616
|
|
617
|
sprintf(form2,"%%*s");
|
618
|
sprintf(format,"%s%%d",form2);
|
619
|
if (sscanf(line,format,&ai) == 1) {
|
620
|
do {
|
621
|
strcat(form2,"%*s");
|
622
|
sprintf(format,"%s%%d",form2);
|
623
|
n = sscanf(line,format,&aj);
|
624
|
if (n == 1) {
|
625
|
srenew(con->conect,++con->nconect);
|
626
|
con->conect[con->nconect-1].ai = ai-1;
|
627
|
con->conect[con->nconect-1].aj = aj-1;
|
628
|
}
|
629
|
} while (n == 1);
|
630
|
}
|
631
|
}
|
632
|
|
633
|
void gmx_conect_dump(FILE *fp,gmx_conect conect)
|
634
|
{
|
635
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
636
|
int i;
|
637
|
|
638
|
for(i=0; (i<gc->nconect); i++)
|
639
|
fprintf(fp,"%6s%5d%5d\n","CONECT",
|
640
|
gc->conect[i].ai+1,gc->conect[i].aj+1);
|
641
|
}
|
642
|
|
643
|
gmx_conect gmx_conect_init()
|
644
|
{
|
645
|
gmx_conect_t *gc;
|
646
|
|
647
|
snew(gc,1);
|
648
|
|
649
|
return (gmx_conect) gc;
|
650
|
}
|
651
|
|
652
|
void gmx_conect_done(gmx_conect conect)
|
653
|
{
|
654
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
655
|
|
656
|
sfree(gc->conect);
|
657
|
}
|
658
|
|
659
|
bool gmx_conect_exist(gmx_conect conect,int ai,int aj)
|
660
|
{
|
661
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
662
|
int i;
|
663
|
|
664
|
|
665
|
|
666
|
|
667
|
for(i=0; (i<gc->nconect); i++)
|
668
|
if (((gc->conect[i].ai == ai) &&
|
669
|
(gc->conect[i].aj == aj)) ||
|
670
|
((gc->conect[i].aj == ai) &&
|
671
|
(gc->conect[i].ai == aj)))
|
672
|
return TRUE;
|
673
|
return FALSE;
|
674
|
}
|
675
|
|
676
|
void gmx_conect_add(gmx_conect conect,int ai,int aj)
|
677
|
{
|
678
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
679
|
int i;
|
680
|
|
681
|
|
682
|
|
683
|
|
684
|
if (!gmx_conect_exist(conect,ai,aj)) {
|
685
|
srenew(gc->conect,++gc->nconect);
|
686
|
gc->conect[gc->nconect-1].ai = ai;
|
687
|
gc->conect[gc->nconect-1].aj = aj;
|
688
|
}
|
689
|
}
|
690
|
|
691
|
int read_pdbfile(FILE *in,char *title,int *model_nr,
|
692
|
t_atoms *atoms,rvec x[],int *ePBC,matrix box,bool bChange,
|
693
|
gmx_conect conect)
|
694
|
{
|
695
|
gmx_conect_t *gc = (gmx_conect_t *)conect;
|
696
|
t_symtab symtab;
|
697
|
bool bCOMPND;
|
698
|
bool bConnWarn = FALSE;
|
699
|
char line[STRLEN+1];
|
700
|
int line_type;
|
701
|
char *c,*d;
|
702
|
int natom,chainnum,nres_ter_prev=0;
|
703
|
char chidmax=' ';
|
704
|
bool bStop=FALSE;
|
705
|
|
706
|
if (ePBC)
|
707
|
{
|
708
|
|
709
|
*ePBC = epbcNONE;
|
710
|
}
|
711
|
if (box != NULL)
|
712
|
{
|
713
|
clear_mat(box);
|
714
|
}
|
715
|
|
716
|
open_symtab(&symtab);
|
717
|
|
718
|
bCOMPND=FALSE;
|
719
|
title[0]='\0';
|
720
|
natom=0;
|
721
|
chainnum=0;
|
722
|
while (!bStop && (fgets2(line,STRLEN,in) != NULL))
|
723
|
{
|
724
|
line_type = line2type(line);
|
725
|
|
726
|
switch(line_type)
|
727
|
{
|
728
|
case epdbATOM:
|
729
|
case epdbHETATM:
|
730
|
natom = read_atom(&symtab,line,line_type,natom,atoms,x,chainnum,bChange);
|
731
|
break;
|
732
|
|
733
|
case epdbANISOU:
|
734
|
if (atoms->pdbinfo)
|
735
|
{
|
736
|
read_anisou(line,natom,atoms);
|
737
|
}
|
738
|
break;
|
739
|
|
740
|
case epdbCRYST1:
|
741
|
read_cryst1(line,ePBC,box);
|
742
|
break;
|
743
|
|
744
|
case epdbTITLE:
|
745
|
case epdbHEADER:
|
746
|
if (strlen(line) > 6)
|
747
|
{
|
748
|
c=line+6;
|
749
|
|
750
|
while (c && (c[0]!=' ')) c++;
|
751
|
while (c && (c[0]==' ')) c++;
|
752
|
|
753
|
d=strstr(c," ");
|
754
|
if (d)
|
755
|
{
|
756
|
d[0]='\0';
|
757
|
}
|
758
|
if (strlen(c)>0)
|
759
|
{
|
760
|
strcpy(title,c);
|
761
|
}
|
762
|
}
|
763
|
break;
|
764
|
|
765
|
case epdbCOMPND:
|
766
|
if ((!strstr(line,": ")) || (strstr(line+6,"MOLECULE:")))
|
767
|
{
|
768
|
if ( !(c=strstr(line+6,"MOLECULE:")) )
|
769
|
{
|
770
|
c=line;
|
771
|
}
|
772
|
|
773
|
while (c && (c[0]!=' ')) c++;
|
774
|
while (c && (c[0]==' ')) c++;
|
775
|
|
776
|
d=strstr(c," ");
|
777
|
if (d)
|
778
|
{
|
779
|
while ( (d[-1]==';') && d>c) d--;
|
780
|
d[0]='\0';
|
781
|
}
|
782
|
if (strlen(c) > 0)
|
783
|
{
|
784
|
if (bCOMPND)
|
785
|
{
|
786
|
strcat(title,"; ");
|
787
|
strcat(title,c);
|
788
|
}
|
789
|
else
|
790
|
{
|
791
|
strcpy(title,c);
|
792
|
}
|
793
|
}
|
794
|
bCOMPND=TRUE;
|
795
|
}
|
796
|
break;
|
797
|
|
798
|
case epdbTER:
|
799
|
chainnum++;
|
800
|
break;
|
801
|
|
802
|
case epdbMODEL:
|
803
|
if(model_nr)
|
804
|
{
|
805
|
sscanf(line,"%*s%d",model_nr);
|
806
|
}
|
807
|
break;
|
808
|
|
809
|
case epdbENDMDL:
|
810
|
bStop=TRUE;
|
811
|
break;
|
812
|
case epdbCONECT:
|
813
|
if (gc)
|
814
|
{
|
815
|
gmx_conect_addline(gc,line);
|
816
|
}
|
817
|
else if (!bConnWarn)
|
818
|
{
|
819
|
fprintf(stderr,"WARNING: all CONECT records are ignored\n");
|
820
|
bConnWarn = TRUE;
|
821
|
}
|
822
|
break;
|
823
|
|
824
|
default:
|
825
|
break;
|
826
|
}
|
827
|
}
|
828
|
|
829
|
free_symtab(&symtab);
|
830
|
return natom;
|
831
|
}
|
832
|
|
833
|
void get_pdb_coordnum(FILE *in,int *natoms)
|
834
|
{
|
835
|
char line[STRLEN];
|
836
|
|
837
|
*natoms=0;
|
838
|
while (fgets2(line,STRLEN,in))
|
839
|
{
|
840
|
if ( strncmp(line,"ENDMDL",6) == 0 )
|
841
|
{
|
842
|
break;
|
843
|
}
|
844
|
if ((strncmp(line,"ATOM ",6) == 0) || (strncmp(line,"HETATM",6) == 0))
|
845
|
{
|
846
|
(*natoms)++;
|
847
|
}
|
848
|
}
|
849
|
}
|
850
|
|
851
|
void read_pdb_conf(const char *infile,char *title,
|
852
|
t_atoms *atoms,rvec x[],int *ePBC,matrix box,bool bChange,
|
853
|
gmx_conect conect)
|
854
|
{
|
855
|
FILE *in;
|
856
|
|
857
|
in = gmx_fio_fopen(infile,"r");
|
858
|
read_pdbfile(in,title,NULL,atoms,x,ePBC,box,bChange,conect);
|
859
|
gmx_fio_fclose(in);
|
860
|
}
|
861
|
|
862
|
gmx_conect gmx_conect_generate(t_topology *top)
|
863
|
{
|
864
|
int f,i;
|
865
|
gmx_conect gc;
|
866
|
|
867
|
|
868
|
gc = gmx_conect_init();
|
869
|
|
870
|
for(f=0; (f<F_NRE); f++) {
|
871
|
if (IS_CHEMBOND(f))
|
872
|
for(i=0; (i<top->idef.il[f].nr); i+=interaction_function[f].nratoms+1) {
|
873
|
gmx_conect_add(gc,top->idef.il[f].iatoms[i+1],
|
874
|
top->idef.il[f].iatoms[i+2]);
|
875
|
}
|
876
|
}
|
877
|
return gc;
|
878
|
}
|