1 /**
2 * @file valist.c
3 * @author Nathan Baker
4 * @brief Class Valist methods
5 * @version $Id$
6 * @attention
7 * @verbatim
8 *
9 * APBS -- Adaptive Poisson-Boltzmann Solver
10 *
11 * Nathan A. Baker (nathan.baker@pnnl.gov)
12 * Pacific Northwest National Laboratory
13 *
14 * Additional contributing authors listed in the code documentation.
15 *
16 * Copyright (c) 2010-2014 Battelle Memorial Institute. Developed at the
17 * Pacific Northwest National Laboratory, operated by Battelle Memorial
18 * Institute, Pacific Northwest Division for the U.S. Department of Energy.
19 *
20 * Portions Copyright (c) 2002-2010, Washington University in St. Louis.
21 * Portions Copyright (c) 2002-2010, Nathan A. Baker.
22 * Portions Copyright (c) 1999-2002, The Regents of the University of
23 * California.
24 * Portions Copyright (c) 1995, Michael Holst.
25 * All rights reserved.
26 *
27 * Redistribution and use in source and binary forms, with or without
28 * modification, are permitted provided that the following conditions are met:
29 *
30 * Redistributions of source code must retain the above copyright notice, this
31 * list of conditions and the following disclaimer.
32 *
33 * Redistributions in binary form must reproduce the above copyright notice,
34 * this list of conditions and the following disclaimer in the documentation
35 * and/or other materials provided with the distribution.
36 *
37 * Neither the name of the developer nor the names of its contributors may be
38 * used to endorse or promote products derived from this software without
39 * specific prior written permission.
40 *
41 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
42 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
43 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
44 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
45 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
46 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
47 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
48 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
49 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
50 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
51 * THE POSSIBILITY OF SUCH DAMAGE.
52 *
53 * @endverbatim
54 */
55
56 #include "valist.h"
57
58 VEMBED(rcsid="$Id$")
59
60 VPRIVATE char *Valist_whiteChars = " \t\r\n";
61 VPRIVATE char *Valist_commChars = "#%";
62 VPRIVATE char *Valist_xmlwhiteChars = " \t\r\n<>";
63
64 #if !defined(VINLINE_VATOM)
65
Valist_getCenterX(Valist * thee)66 VPUBLIC double Valist_getCenterX(Valist *thee) {
67
68 if (thee == NULL) {
69 Vnm_print(2, "Valist_getCenterX: Found null pointer when getting the center of X coordinate!\n");
70 VASSERT(0);
71 }
72 return thee->center[0];
73
74 }
75
Valist_getCenterY(Valist * thee)76 VPUBLIC double Valist_getCenterY(Valist *thee) {
77
78 if (thee == NULL) {
79 Vnm_print(2, "Valist_getCenterY: Found null pointer when getting the center of Y coordinate!\n");
80 VASSERT(0);
81 }
82 return thee->center[1];
83
84 }
Valist_getCenterZ(Valist * thee)85 VPUBLIC double Valist_getCenterZ(Valist *thee) {
86
87 if (thee == NULL) {
88 Vnm_print(2, "Valist_getCenterZ: Found null pointer when getting the center of Z coordinate!\n");
89 VASSERT(0);
90 }
91 return thee->center[2];
92
93 }
94
Valist_getAtomList(Valist * thee)95 VPUBLIC Vatom* Valist_getAtomList(Valist *thee) {
96
97 if (thee == NULL) {
98 Vnm_print(2, "Valist_getAtomList: Found null pointer when getting the atom list!\n");
99 VASSERT(0);
100 }
101 return thee->atoms;
102
103 }
104
Valist_getNumberAtoms(Valist * thee)105 VPUBLIC int Valist_getNumberAtoms(Valist *thee) {
106
107 if (thee == NULL) {
108 Vnm_print(2, "Valist_getNumberAtoms: Found null pointer when getting the number of atoms!\n");
109 VASSERT(0);
110 }
111 return thee->number;
112
113 }
114
Valist_getAtom(Valist * thee,int i)115 VPUBLIC Vatom* Valist_getAtom(Valist *thee, int i) {
116
117 if (thee == NULL) {
118 Vnm_print(2, "Valist_getAtom: Found null pointer when getting atoms!\n");
119 VASSERT(0);
120 }
121 if (i >= thee->number) {
122 Vnm_print(2, "Valist_getAtom: Requested atom number (%d) outside of atom list range (%d)!\n", i, thee->number);
123 VASSERT(0);
124 }
125 return &(thee->atoms[i]);
126
127 }
128
Valist_memChk(Valist * thee)129 VPUBLIC unsigned long int Valist_memChk(Valist *thee) {
130
131 if (thee == NULL) return 0;
132 return Vmem_bytes(thee->vmem);
133
134 }
135
136 #endif /* if !defined(VINLINE_VATOM) */
137
Valist_ctor()138 VPUBLIC Valist* Valist_ctor() {
139
140 /* Set up the structure */
141 Valist *thee = VNULL;
142 thee = (Valist*)Vmem_malloc(VNULL, 1, sizeof(Valist));
143 if ( thee == VNULL) {
144 Vnm_print(2, "Valist_ctor: Got NULL pointer when constructing the atom list object!\n");
145 VASSERT(0);
146 }
147 if ( Valist_ctor2(thee) != VRC_SUCCESS) {
148 Vnm_print(2, "Valist_ctor: Error in constructing the atom list object!\n");
149 VASSERT(0);
150 }
151
152 return thee;
153 }
154
Valist_ctor2(Valist * thee)155 VPUBLIC Vrc_Codes Valist_ctor2(Valist *thee) {
156
157 thee->atoms = VNULL;
158 thee->number = 0;
159
160 /* Initialize the memory management object */
161 thee->vmem = Vmem_ctor("APBS:VALIST");
162
163 return VRC_SUCCESS;
164
165 }
166
Valist_dtor(Valist ** thee)167 VPUBLIC void Valist_dtor(Valist **thee)
168 {
169 if ((*thee) != VNULL) {
170 Valist_dtor2(*thee);
171 Vmem_free(VNULL, 1, sizeof(Valist), (void **)thee);
172 (*thee) = VNULL;
173 }
174 }
175
Valist_dtor2(Valist * thee)176 VPUBLIC void Valist_dtor2(Valist *thee) {
177
178 Vmem_free(thee->vmem, thee->number, sizeof(Vatom), (void **)&(thee->atoms));
179 thee->atoms = VNULL;
180 thee->number = 0;
181
182 Vmem_dtor(&(thee->vmem));
183 }
184
185 /* Read serial number from PDB ATOM/HETATM field */
Valist_readPDBSerial(Valist * thee,Vio * sock,int * serial)186 VPRIVATE Vrc_Codes Valist_readPDBSerial(Valist *thee, Vio *sock, int *serial) {
187
188 char tok[VMAX_BUFSIZE];
189 int ti = 0;
190
191 if (Vio_scanf(sock, "%s", tok) != 1) {
192 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing serial!\n");
193 return VRC_FAILURE;
194 }
195 if (sscanf(tok, "%d", &ti) != 1) {
196 Vnm_print(2, "Valist_readPDB: Unable to parse serial token (%s) as int!\n",
197 tok);
198 return VRC_FAILURE;
199 }
200 *serial = ti;
201
202 return VRC_SUCCESS;
203 }
204
205 /* Read atom name from PDB ATOM/HETATM field */
Valist_readPDBAtomName(Valist * thee,Vio * sock,char atomName[VMAX_ARGLEN])206 VPRIVATE Vrc_Codes Valist_readPDBAtomName(Valist *thee, Vio *sock,
207 char atomName[VMAX_ARGLEN]) {
208
209 char tok[VMAX_BUFSIZE];
210
211 if (Vio_scanf(sock, "%s", tok) != 1) {
212 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing atom name!\n");
213 return VRC_FAILURE;
214 }
215 if (strlen(tok) < VMAX_ARGLEN) strcpy(atomName, tok);
216 else {
217 Vnm_print(2, "Valist_readPDB: Atom name (%s) too long!\n", tok);
218 return VRC_FAILURE;
219 }
220 return VRC_SUCCESS;
221 }
222
223 /* Read residue name from PDB ATOM/HETATM field */
Valist_readPDBResidueName(Valist * thee,Vio * sock,char resName[VMAX_ARGLEN])224 VPRIVATE Vrc_Codes Valist_readPDBResidueName(Valist *thee, Vio *sock,
225 char resName[VMAX_ARGLEN]) {
226
227 char tok[VMAX_BUFSIZE];
228
229 if (Vio_scanf(sock, "%s", tok) != 1) {
230 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing residue name!\n");
231 return VRC_FAILURE;
232 }
233 if (strlen(tok) < VMAX_ARGLEN) strcpy(resName, tok);
234 else {
235 Vnm_print(2, "Valist_readPDB: Residue name (%s) too long!\n", tok);
236 return VRC_FAILURE;
237 }
238 return VRC_SUCCESS;
239 }
240
241 /* Read residue number from PDB ATOM/HETATM field */
Valist_readPDBResidueNumber(Valist * thee,Vio * sock,int * resSeq)242 VPRIVATE Vrc_Codes Valist_readPDBResidueNumber(
243 Valist *thee, Vio *sock, int *resSeq) {
244
245 char tok[VMAX_BUFSIZE];
246 char *resstring;
247 int ti = 0;
248
249 if (Vio_scanf(sock, "%s", tok) != 1) {
250 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing resSeq!\n");
251 return VRC_FAILURE;
252 }
253 if (sscanf(tok, "%d", &ti) != 1) {
254
255 /* One of three things can happen here:
256 1) There is a chainID in the line: THR A 1
257 2) The chainID is merged with resSeq: THR A1001
258 3) An actual error: THR foo
259
260 */
261
262 if (strlen(tok) == 1) {
263 /* Case 1: Chain ID Present
264 Read the next field and hope its a float */
265
266 if (Vio_scanf(sock, "%s", tok) != 1) {
267 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing resSeq!\n");
268 return VRC_FAILURE;
269 }
270 if (sscanf(tok, "%d", &ti) != 1) {
271 Vnm_print(2, "Valist_readPDB: Unable to parse resSeq token (%s) as int!\n",
272 tok);
273 return VRC_FAILURE;
274 }
275
276 } else {
277 /* Case 2: Chain ID, merged string.
278 Move pointer forward past the chainID and check
279 */
280 //strcpy(resstring, tok);
281 resstring = tok;
282 resstring++;
283
284 if (sscanf(resstring, "%d", &ti) != 1) {
285 /* Case 3: More than one non-numeral char is present. Error.*/
286 Vnm_print(2, "Valist_readPDB: Unable to parse resSeq token (%s) as int!\n",
287 resstring);
288 return VRC_FAILURE;
289 }
290 }
291 }
292 *resSeq = ti;
293
294 return VRC_SUCCESS;
295 }
296
297 /* Read atom coordinate from PDB ATOM/HETATM field */
Valist_readPDBAtomCoord(Valist * thee,Vio * sock,double * coord)298 VPRIVATE Vrc_Codes Valist_readPDBAtomCoord(Valist *thee, Vio *sock, double *coord) {
299
300 char tok[VMAX_BUFSIZE];
301 double tf = 0;
302
303 if (Vio_scanf(sock, "%s", tok) != 1) {
304 Vnm_print(2, "Valist_readPDB: Ran out of tokens while parsing atom coordinate!\n");
305 return VRC_FAILURE;
306 }
307 if (sscanf(tok, "%lf", &tf) != 1) {
308 return VRC_FAILURE;
309 }
310 *coord = tf;
311
312 return VRC_SUCCESS;
313 }
314
315 /* Read charge and radius from PQR ATOM/HETATM field */
Valist_readPDBChargeRadius(Valist * thee,Vio * sock,double * charge,double * radius)316 VPRIVATE Vrc_Codes Valist_readPDBChargeRadius(Valist *thee, Vio *sock,
317 double *charge, double *radius) {
318
319 char tok[VMAX_BUFSIZE];
320 double tf = 0;
321
322 if (Vio_scanf(sock, "%s", tok) != 1) {
323 Vnm_print(2, "Valist_readPQR: Ran out of tokens while parsing charge!\n");
324 return VRC_FAILURE;
325 }
326 if (sscanf(tok, "%lf", &tf) != 1) {
327 return VRC_FAILURE;
328 }
329 *charge = tf;
330
331 if (Vio_scanf(sock, "%s", tok) != 1) {
332 Vnm_print(2, "Valist_readPQR: Ran out of tokens while parsing radius!\n");
333 return VRC_FAILURE;
334 }
335 if (sscanf(tok, "%lf", &tf) != 1) {
336 return VRC_FAILURE;
337 }
338 *radius = tf;
339
340 return VRC_SUCCESS;
341 }
342
343 /* Read ATOM/HETATM field of PDB through the X/Y/Z fields */
Valist_readPDB_throughXYZ(Valist * thee,Vio * sock,int * serial,char atomName[VMAX_ARGLEN],char resName[VMAX_ARGLEN],int * resSeq,double * x,double * y,double * z)344 VPRIVATE Vrc_Codes Valist_readPDB_throughXYZ(
345 Valist *thee,
346 Vio *sock, /* Socket ready for reading */
347 int *serial, /* Set to atom number */
348 char atomName[VMAX_ARGLEN], /* Set to atom name */
349 char resName[VMAX_ARGLEN], /* Set to residue name */
350 int *resSeq, /* Set to residue number */
351 double *x, /* Set to x-coordinate */
352 double *y, /* Set to y-coordinate */
353 double *z /* Set to z-coordinate */
354 ) {
355
356
357 int i, njunk, gotit;
358
359 /* Grab serial */
360 if (Valist_readPDBSerial(thee, sock, serial) == VRC_FAILURE) {
361 Vnm_print(2, "Valist_readPDB: Error while parsing serial!\n");
362 }
363
364 /* Grab atom name */
365 if (Valist_readPDBAtomName(thee, sock, atomName) == VRC_FAILURE) {
366 Vnm_print(2, "Valist_readPDB: Error while parsing atom name!\n");
367 return VRC_FAILURE;
368 }
369
370 /* Grab residue name */
371 if (Valist_readPDBResidueName(thee, sock, resName) == VRC_FAILURE) {
372 Vnm_print(2, "Valist_readPDB: Error while parsing residue name!\n");
373 return VRC_FAILURE;
374 }
375
376
377 /* Grab residue number */
378 if (Valist_readPDBResidueNumber(thee, sock, resSeq) == VRC_FAILURE) {
379 Vnm_print(2, "Valist_readPDB: Error while parsing residue number!\n");
380 return VRC_FAILURE;
381 }
382
383
384 /* Read tokens until we find one that can be parsed as an atom
385 * x-coordinate. We will allow njunk=1 intervening field that
386 * cannot be parsed as a coordinate */
387 njunk = 1;
388 gotit = 0;
389 for (i=0; i<(njunk+1); i++) {
390 if (Valist_readPDBAtomCoord(thee, sock, x) == VRC_SUCCESS) {
391 gotit = 1;
392 break;
393 }
394 }
395 if (!gotit) {
396 Vnm_print(2, "Valist_readPDB: Can't find x!\n");
397 return VRC_FAILURE;
398 }
399 /* Read y-coordinate */
400 if (Valist_readPDBAtomCoord(thee, sock, y) == VRC_FAILURE) {
401 Vnm_print(2, "Valist_readPDB: Can't find y!\n");
402 return VRC_FAILURE;
403 }
404 /* Read z-coordinate */
405 if (Valist_readPDBAtomCoord(thee, sock, z) == VRC_FAILURE) {
406 Vnm_print(2, "Valist_readPDB: Can't find z!\n");
407 return VRC_FAILURE;
408 }
409
410 #if 0 /* Set to 1 if you want to debug */
411 Vnm_print(1, "Valist_readPDB: serial = %d\n", *serial);
412 Vnm_print(1, "Valist_readPDB: atomName = %s\n", atomName);
413 Vnm_print(1, "Valist_readPDB: resName = %s\n", resName);
414 Vnm_print(1, "Valist_readPDB: resSeq = %d\n", *resSeq);
415 Vnm_print(1, "Valist_readPDB: pos = (%g, %g, %g)\n",
416 *x, *y, *z);
417 #endif
418
419 return VRC_SUCCESS;
420 }
421
422 /* Get a the next available atom storage location, increasing the storage
423 * space if necessary. Return VNULL if something goes wrong. */
Valist_getAtomStorage(Valist * thee,Vatom ** plist,int * pnlist,int * pnatoms)424 VPRIVATE Vatom* Valist_getAtomStorage(
425 Valist *thee,
426 Vatom **plist, /* Pointer to existing list of atoms */
427 int *pnlist, /* Size of existing list, may be changed */
428 int *pnatoms /* Existing number of atoms in list; incremented
429 before exit */
430 ) {
431
432 Vatom *oldList, *newList, *theList;
433 Vatom *oldAtom, *newAtom;
434 int iatom, inext, oldLength, newLength, natoms;
435
436 newList = VNULL;
437
438 /* See if we need more space */
439 if (*pnatoms >= *pnlist) {
440
441 /* Double the storage space */
442 oldLength = *pnlist;
443 newLength = 2*oldLength;
444 newList = (Vatom*)Vmem_malloc(thee->vmem, newLength, sizeof(Vatom));
445 oldList = *plist;
446
447 /* Check the allocation */
448 if (newList == VNULL) {
449 Vnm_print(2, "Valist_readPDB: failed to allocate space for %d (Vatom)s!\n", newLength);
450 return VNULL;
451 }
452
453 /* Copy the atoms over */
454 natoms = *pnatoms;
455 for (iatom=0; iatom<natoms; iatom++) {
456 oldAtom = &(oldList[iatom]);
457 newAtom = &(newList[iatom]);
458 Vatom_copyTo(oldAtom, newAtom);
459 Vatom_dtor2(oldAtom);
460 }
461
462 /* Free the old list */
463 Vmem_free(thee->vmem, oldLength, sizeof(Vatom), (void **)plist);
464
465 /* Copy new list to plist */
466 *plist = newList;
467 *pnlist = newLength;
468 }
469
470 theList = *plist;
471 inext = *pnatoms;
472
473 /* Get the next available spot and increment counters */
474 newAtom = &(theList[inext]);
475 *pnatoms = inext + 1;
476
477 return newAtom;
478 }
479
Valist_setAtomArray(Valist * thee,Vatom ** plist,int nlist,int natoms)480 VPRIVATE Vrc_Codes Valist_setAtomArray(Valist *thee,
481 Vatom **plist, /* Pointer to list of atoms to store */
482 int nlist, /* Length of list */
483 int natoms /* Number of real atom entries in list */
484 ) {
485
486 Vatom *list, *newAtom, *oldAtom;
487 int i;
488
489 list = *plist;
490
491 /* Allocate necessary space */
492 thee->number = 0;
493 thee->atoms = (Vatom*)Vmem_malloc(thee->vmem, natoms, sizeof(Vatom));
494 if (thee->atoms == VNULL) {
495 Vnm_print(2, "Valist_readPDB: Unable to allocate space for %d (Vatom)s!\n",
496 natoms);
497 return VRC_FAILURE;
498 }
499 thee->number = natoms;
500
501 /* Copy over data */
502 for (i=0; i<thee->number; i++) {
503 newAtom = &(thee->atoms[i]);
504 oldAtom = &(list[i]);
505 Vatom_copyTo(oldAtom, newAtom);
506 Vatom_dtor2(oldAtom);
507 }
508
509 /* Free old array */
510 Vmem_free(thee->vmem, nlist, sizeof(Vatom), (void **)plist);
511
512 return VRC_SUCCESS;
513 }
514
Valist_readPDB(Valist * thee,Vparam * param,Vio * sock)515 VPUBLIC Vrc_Codes Valist_readPDB(Valist *thee, Vparam *param, Vio *sock) {
516
517 /* WE DO NOT DIRECTLY CONFORM TO PDB STANDARDS -- TO ALLOW LARGER FILES, WE
518 * REQUIRE ALL FIELDS TO BE WHITESPACE DELIMITED */
519
520 Vatom *atoms = VNULL;
521 Vatom *nextAtom = VNULL;
522 Vparam_AtomData *atomData = VNULL;
523
524 char tok[VMAX_BUFSIZE];
525 char atomName[VMAX_ARGLEN], resName[VMAX_ARGLEN];
526
527 int nlist, natoms, serial, resSeq;
528
529 double x, y, z, charge, radius, epsilon;
530 double pos[3];
531
532 if (thee == VNULL) {
533 Vnm_print(2, "Valist_readPDB: Got NULL pointer when reading PDB file!\n");
534 VASSERT(0);
535 }
536 thee->number = 0;
537
538 Vio_setWhiteChars(sock, Valist_whiteChars);
539 Vio_setCommChars(sock, Valist_commChars);
540
541 /* Allocate some initial space for the atoms */
542 nlist = 200;
543 atoms = (Vatom*)Vmem_malloc(thee->vmem, nlist, sizeof(Vatom));
544
545 natoms = 0;
546 /* Read until we run out of lines */
547 while (Vio_scanf(sock, "%s", tok) == 1) {
548
549 /* Parse only ATOM/HETATOM fields */
550 if ((Vstring_strcasecmp(tok, "ATOM") == 0) ||
551 (Vstring_strcasecmp(tok, "HETATM") == 0)) {
552
553 /* Read ATOM/HETATM field of PDB through the X/Y/Z fields */
554 if (Valist_readPDB_throughXYZ(thee, sock, &serial, atomName,
555 resName, &resSeq, &x, &y, &z) == VRC_FAILURE) {
556 Vnm_print(2, "Valist_readPDB: Error parsing atom %d!\n",
557 serial);
558 return VRC_FAILURE;
559 }
560
561 /* Try to find the parameters. */
562 atomData = Vparam_getAtomData(param, resName, atomName);
563 if (atomData == VNULL) {
564 Vnm_print(2, "Valist_readPDB: Couldn't find parameters for \
565 atom = %s, residue = %s\n", atomName, resName);
566 return VRC_FAILURE;
567 }
568 charge = atomData->charge;
569 radius = atomData->radius;
570 epsilon = atomData->epsilon;
571
572 /* Get pointer to next available atom position */
573 nextAtom = Valist_getAtomStorage(thee, &atoms, &nlist, &natoms);
574 if (nextAtom == VNULL) {
575 Vnm_print(2, "Valist_readPDB: Error in allocating spacing for atoms!\n");
576 return VRC_FAILURE;
577 }
578
579 /* Store the information */
580 pos[0] = x; pos[1] = y; pos[2] = z;
581 Vatom_setPosition(nextAtom, pos);
582 Vatom_setCharge(nextAtom, charge);
583 Vatom_setRadius(nextAtom, radius);
584 Vatom_setEpsilon(nextAtom, epsilon);
585 Vatom_setAtomID(nextAtom, natoms-1);
586 Vatom_setResName(nextAtom, resName);
587 Vatom_setAtomName(nextAtom, atomName);
588
589 } /* if ATOM or HETATM */
590 } /* while we haven't run out of tokens */
591
592 Vnm_print(0, "Valist_readPDB: Counted %d atoms\n", natoms);
593 fflush(stdout);
594
595 /* Store atoms internally */
596 if (Valist_setAtomArray(thee, &atoms, nlist, natoms) == VRC_FAILURE) {
597 Vnm_print(2, "Valist_readPDB: unable to store atoms!\n");
598 return VRC_FAILURE;
599 }
600
601 return Valist_getStatistics(thee);
602
603
604 }
605
Valist_readPQR(Valist * thee,Vparam * params,Vio * sock)606 VPUBLIC Vrc_Codes Valist_readPQR(Valist *thee, Vparam *params, Vio *sock) {
607
608 /* WE DO NOT DIRECTLY CONFORM TO PDB STANDARDS -- TO ALLOW LARGER FILES, WE
609 * REQUIRE ALL FIELDS TO BE WHITESPACE DELIMITED */
610
611
612 Vatom *atoms = VNULL;
613 Vatom *nextAtom = VNULL;
614 Vparam_AtomData *atomData = VNULL;
615
616 char tok[VMAX_BUFSIZE];
617 char atomName[VMAX_ARGLEN], resName[VMAX_ARGLEN];
618 char chs[VMAX_BUFSIZE];
619
620 int use_params = 0;
621 int nlist, natoms, serial, resSeq;
622
623 double x, y, z, charge, radius, epsilon;
624 double pos[3];
625
626 epsilon = 0.0;
627
628 if (thee == VNULL) {
629 Vnm_print(2, "Valist_readPQR: Got NULL pointer when reading PQR file!\n");
630 VASSERT(0);
631 }
632 thee->number = 0;
633
634 Vio_setWhiteChars(sock, Valist_whiteChars);
635 Vio_setCommChars(sock, Valist_commChars);
636
637 /* Allocate some initial space for the atoms */
638 nlist = 200;
639 atoms = (Vatom*)Vmem_malloc(thee->vmem, nlist, sizeof(Vatom));
640
641 /* Check if we are using a parameter file or not */
642 if(params != VNULL) use_params = 1;
643
644 natoms = 0;
645 /* Read until we run out of lines */
646 while (Vio_scanf(sock, "%s", tok) == 1) {
647
648 /* Parse only ATOM/HETATOM fields */
649 if ((Vstring_strcasecmp(tok, "ATOM") == 0) ||
650 (Vstring_strcasecmp(tok, "HETATM") == 0)) {
651
652 /* Read ATOM/HETATM field of PDB through the X/Y/Z fields */
653 if (Valist_readPDB_throughXYZ(thee, sock, &serial, atomName,
654 resName, &resSeq, &x, &y, &z) == VRC_FAILURE) {
655 Vnm_print(2, "Valist_readPQR: Error parsing atom %d!\n",serial);
656 Vnm_print(2, "Please double check this atom in the pqr file, e.g., make sure there are no concatenated fields.\n");
657 return VRC_FAILURE;
658 }
659
660 /* Read Q/R fields */
661 if (Valist_readPDBChargeRadius(thee, sock, &charge, &radius) == VRC_FAILURE) {
662 Vnm_print(2, "Valist_readPQR: Error parsing atom %d!\n",
663 serial);
664 Vnm_print(2, "Please double check this atom in the pqr file, e.g., make sure there are no concatenated fields.\n");
665 return VRC_FAILURE;
666 }
667
668 if(use_params){
669 /* Try to find the parameters. */
670 atomData = Vparam_getAtomData(params, resName, atomName);
671 if (atomData == VNULL) {
672 Vnm_print(2, "Valist_readPDB: Couldn't find parameters for \
673 atom = %s, residue = %s\n", atomName, resName);
674 return VRC_FAILURE;
675 }
676 charge = atomData->charge;
677 radius = atomData->radius;
678 epsilon = atomData->epsilon;
679 }
680
681 /* Get pointer to next available atom position */
682 nextAtom = Valist_getAtomStorage(thee, &atoms, &nlist, &natoms);
683 if (nextAtom == VNULL) {
684 Vnm_print(2, "Valist_readPQR: Error in allocating spacing for atoms!\n");
685 return VRC_FAILURE;
686 }
687
688 /* Store the information */
689 pos[0] = x; pos[1] = y; pos[2] = z;
690 Vatom_setPosition(nextAtom, pos);
691 Vatom_setCharge(nextAtom, charge);
692 Vatom_setRadius(nextAtom, radius);
693 Vatom_setEpsilon(nextAtom, epsilon);
694 Vatom_setAtomID(nextAtom, natoms-1);
695 Vatom_setResName(nextAtom, resName);
696 Vatom_setAtomName(nextAtom, atomName);
697
698 } /* if ATOM or HETATM */
699 else {
700 /*
701 * nop
702 * Note that if we find a line that starts with something that's not
703 * ATOM or HETATM we'll just keep parsing strings until we find one
704 * of the acceptable keywords.
705 * Extraordinary measures are not necessary, and only add to the
706 * befuddlement.
707 */
708 }
709 } /* while we haven't run out of tokens */
710
711 Vnm_print(0, "Valist_readPQR: Counted %d atoms\n", natoms);
712 fflush(stdout);
713
714 /* Store atoms internally */
715 if (Valist_setAtomArray(thee, &atoms, nlist, natoms) == VRC_FAILURE) {
716 Vnm_print(2, "Valist_readPDB: unable to store atoms!\n");
717 return VRC_FAILURE;
718 }
719
720 return Valist_getStatistics(thee);
721
722
723 }
724
Valist_readXML(Valist * thee,Vparam * params,Vio * sock)725 VPUBLIC Vrc_Codes Valist_readXML(Valist *thee, Vparam *params, Vio *sock) {
726
727 Vatom *atoms = VNULL;
728 Vatom *nextAtom = VNULL;
729
730 char tok[VMAX_BUFSIZE];
731 char endtag[VMAX_BUFSIZE];
732
733 int nlist, natoms;
734 int xset, yset, zset, chgset, radset;
735
736 double x, y, z, charge, radius, dtmp;
737 double pos[3];
738
739 if (thee == VNULL) {
740 Vnm_print(2, "Valist_readXML: Got NULL pointer when reading XML file!\n");
741 VASSERT(0);
742 }
743 thee->number = 0;
744
745 Vio_setWhiteChars(sock, Valist_xmlwhiteChars);
746 Vio_setCommChars(sock, Valist_commChars);
747
748 /* Allocate some initial space for the atoms */
749 nlist = 200;
750 atoms = (Vatom*)Vmem_malloc(thee->vmem, nlist, sizeof(Vatom));
751
752 /* Initialize some variables */
753 natoms = 0;
754 xset = 0;
755 yset = 0;
756 zset = 0;
757 chgset = 0;
758 radset = 0;
759 strcpy(endtag,"/");
760
761 if(params == VNULL){
762 Vnm_print(1,"\nValist_readXML: Warning Warning Warning Warning Warning\n");
763 Vnm_print(1,"Valist_readXML: The use of XML input files with parameter\n");
764 Vnm_print(1,"Valist_readXML: files is currently not supported.\n");
765 Vnm_print(1,"Valist_readXML: Warning Warning Warning Warning Warning\n\n");
766 }
767
768 /* Read until we run out of lines */
769 while (Vio_scanf(sock, "%s", tok) == 1) {
770
771 /* The first tag taken is the start tag - save it to detect end */
772 if (Vstring_strcasecmp(endtag, "/") == 0) strcat(endtag, tok);
773
774 if (Vstring_strcasecmp(tok, "x") == 0) {
775 Vio_scanf(sock, "%s", tok);
776 if (sscanf(tok, "%lf", &dtmp) != 1) {
777 Vnm_print(2, "Valist_readXML: Unexpected token (%s) while \
778 reading x!\n", tok);
779 return VRC_FAILURE;
780 }
781 x = dtmp;
782 xset = 1;
783 } else if (Vstring_strcasecmp(tok, "y") == 0) {
784 Vio_scanf(sock, "%s", tok);
785 if (sscanf(tok, "%lf", &dtmp) != 1) {
786 Vnm_print(2, "Valist_readXML: Unexpected token (%s) while \
787 reading y!\n", tok);
788 return VRC_FAILURE;
789 }
790 y = dtmp;
791 yset = 1;
792 } else if (Vstring_strcasecmp(tok, "z") == 0) {
793 Vio_scanf(sock, "%s", tok);
794 if (sscanf(tok, "%lf", &dtmp) != 1) {
795 Vnm_print(2, "Valist_readXML: Unexpected token (%s) while \
796 reading z!\n", tok);
797 return VRC_FAILURE;
798 }
799 z = dtmp;
800 zset = 1;
801 } else if (Vstring_strcasecmp(tok, "charge") == 0) {
802 Vio_scanf(sock, "%s", tok);
803 if (sscanf(tok, "%lf", &dtmp) != 1) {
804 Vnm_print(2, "Valist_readXML: Unexpected token (%s) while \
805 reading charge!\n", tok);
806 return VRC_FAILURE;
807 }
808 charge = dtmp;
809 chgset = 1;
810 } else if (Vstring_strcasecmp(tok, "radius") == 0) {
811 Vio_scanf(sock, "%s", tok);
812 if (sscanf(tok, "%lf", &dtmp) != 1) {
813 Vnm_print(2, "Valist_readXML: Unexpected token (%s) while \
814 reading radius!\n", tok);
815 return VRC_FAILURE;
816 }
817 radius = dtmp;
818 radset = 1;
819 } else if (Vstring_strcasecmp(tok, "/atom") == 0) {
820
821 /* Get pointer to next available atom position */
822 nextAtom = Valist_getAtomStorage(thee, &atoms, &nlist, &natoms);
823 if (nextAtom == VNULL) {
824 Vnm_print(2, "Valist_readXML: Error in allocating spacing for atoms!\n");
825 return VRC_FAILURE;
826 }
827
828 if (xset && yset && zset && chgset && radset){
829
830 /* Store the information */
831 pos[0] = x; pos[1] = y; pos[2] = z;
832 Vatom_setPosition(nextAtom, pos);
833 Vatom_setCharge(nextAtom, charge);
834 Vatom_setRadius(nextAtom, radius);
835 Vatom_setAtomID(nextAtom, natoms-1);
836
837 /* Reset the necessary flags */
838 xset = 0;
839 yset = 0;
840 zset = 0;
841 chgset = 0;
842 radset = 0;
843 } else {
844 Vnm_print(2, "Valist_readXML: Missing field(s) in atom tag:\n");
845 if (!xset) Vnm_print(2,"\tx value not set!\n");
846 if (!yset) Vnm_print(2,"\ty value not set!\n");
847 if (!zset) Vnm_print(2,"\tz value not set!\n");
848 if (!chgset) Vnm_print(2,"\tcharge value not set!\n");
849 if (!radset) Vnm_print(2,"\tradius value not set!\n");
850 return VRC_FAILURE;
851 }
852 } else if (Vstring_strcasecmp(tok, endtag) == 0) break;
853 }
854
855 Vnm_print(0, "Valist_readXML: Counted %d atoms\n", natoms);
856 fflush(stdout);
857
858 /* Store atoms internally */
859 if (Valist_setAtomArray(thee, &atoms, nlist, natoms) == VRC_FAILURE) {
860 Vnm_print(2, "Valist_readXML: unable to store atoms!\n");
861 return VRC_FAILURE;
862 }
863
864 return Valist_getStatistics(thee);
865
866 }
867
868 /* Load up Valist with various statistics */
Valist_getStatistics(Valist * thee)869 VPUBLIC Vrc_Codes Valist_getStatistics(Valist *thee) {
870
871 Vatom *atom;
872 int i, j;
873
874 if (thee == VNULL) {
875 Vnm_print(2, "Valist_getStatistics: Got NULL pointer when loading up Valist with various statistics!\n");
876 VASSERT(0);
877 }
878
879 thee->center[0] = 0.;
880 thee->center[1] = 0.;
881 thee->center[2] = 0.;
882 thee->maxrad = 0.;
883 thee->charge = 0.;
884
885 if (thee->number == 0) return VRC_FAILURE;
886
887 /* Reset stat variables */
888 atom = &(thee->atoms[0]);
889 for (i=0; i<3; i++) {
890 thee->maxcrd[i] = thee->mincrd[i] = atom->position[i];
891 }
892 thee->maxrad = atom->radius;
893 thee->charge = 0.0;
894
895 for (i=0; i<thee->number; i++) {
896
897 atom = &(thee->atoms[i]);
898 for (j=0; j<3; j++) {
899 if (atom->position[j] < thee->mincrd[j])
900 thee->mincrd[j] = atom->position[j];
901 if (atom->position[j] > thee->maxcrd[j])
902 thee->maxcrd[j] = atom->position[j];
903 }
904 if (atom->radius > thee->maxrad) thee->maxrad = atom->radius;
905 thee->charge = thee->charge + atom->charge;
906 }
907
908 thee->center[0] = 0.5*(thee->maxcrd[0] + thee->mincrd[0]);
909 thee->center[1] = 0.5*(thee->maxcrd[1] + thee->mincrd[1]);
910 thee->center[2] = 0.5*(thee->maxcrd[2] + thee->mincrd[2]);
911
912 Vnm_print(0, "Valist_getStatistics: Max atom coordinate: (%g, %g, %g)\n",
913 thee->maxcrd[0], thee->maxcrd[1], thee->maxcrd[2]);
914 Vnm_print(0, "Valist_getStatistics: Min atom coordinate: (%g, %g, %g)\n",
915 thee->mincrd[0], thee->mincrd[1], thee->mincrd[2]);
916 Vnm_print(0, "Valist_getStatistics: Molecule center: (%g, %g, %g)\n",
917 thee->center[0], thee->center[1], thee->center[2]);
918
919 return VRC_SUCCESS;
920 }
921