00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111 #include "largefiles.h"
00112
00113 #include <stdio.h>
00114 #include <stdlib.h>
00115 #include <string.h>
00116 #include <ctype.h>
00117 #include <time.h>
00118 #include <math.h>
00119 #include <errno.h>
00120 #include "molfile_plugin.h"
00121
00122 #include "periodic_table.h"
00123 #define THISPLUGIN plugin
00124 #include "vmdconio.h"
00125
00126
00127
00128
00129
00130
00131 #include <expat.h>
00132
00138 #define HOOMD_FORMAT_MAJV 1
00139 #define HOOMD_FORMAT_MINV 4
00140
00141 #define SMALL 1.0e-7f
00143
00144 #define HOOMD_NONE 0
00145 #define HOOMD_XML 1
00146 #define HOOMD_CONFIG 2
00147 #define HOOMD_BOX 3
00148 #define HOOMD_POS 4
00149 #define HOOMD_IMAGE 5
00150 #define HOOMD_VEL 6
00151 #define HOOMD_TYPE 7
00152 #define HOOMD_BODY 8
00153 #define HOOMD_MASS 9
00154 #define HOOMD_CHARGE 10
00155 #define HOOMD_DIAMETER 11
00156 #define HOOMD_BOND 12
00157 #define HOOMD_ANGLE 13
00158 #define HOOMD_DIHEDRAL 14
00159 #define HOOMD_IMPROPER 15
00160 #define HOOMD_WALL 16
00161 #define HOOMD_COORD 17
00162 #define HOOMD_ACCEL 18
00163 #define HOOMD_ORIENT 19
00164 #define HOOMD_UNKNOWN 20
00165 #define HOOMD_NUMTAGS 21
00167
00168 static const char *hoomd_tag_name[HOOMD_NUMTAGS] = {
00169 "(none)", "hoomd_xml", "configuration", "box",
00170 "position", "image", "velocity",
00171 "type", "body", "mass", "charge", "diameter",
00172 "bond", "angle", "dihedral", "improper",
00173 "wall", "coord", "acceleration", "orientation",
00174 "(error)"
00175 };
00176
00178 #define HOOMD_MAXDEPTH 5
00179
00180 typedef struct {
00181 FILE *fp;
00182 XML_Parser p;
00183 void *buffer;
00184 int parsedepth;
00185 int parse_error;
00186 int currtag[HOOMD_MAXDEPTH+1];
00187 int counter;
00188 int majv, minv;
00189 int optflags;
00190 int numframe;
00191 int doneframe;
00192 int numdims;
00193 int numatoms;
00194 int numtypes;
00195 int numbonds;
00196 int numangles;
00197 int numdihedrals;
00198 int numimpropers;
00199 int numbondtypes;
00200 int numangletypes;
00201 int numdihedraltypes;
00202 int numimpropertypes;
00203 int *from;
00204 int *to;
00205 int *bondtype;
00206 int *angle;
00207 int *dihedral;
00208 int *improper;
00209 int *angletype;
00210 int *dihedraltype;
00211 int *impropertype;
00212 char **bondtypename;
00213 char **angletypename;
00214 char **dihedraltypename;
00215 char **impropertypename;
00216 char *filename;
00217 float mysigma;
00218 int *imagecounts;
00219 molfile_atom_t *atomlist;
00220 #if vmdplugin_ABIVERSION > 10
00221 molfile_timestep_metadata_t ts_meta;
00222 #endif
00223 molfile_timestep_t ts;
00224 } hoomd_data_t;
00225
00226
00227 static molfile_plugin_t plugin;
00228
00229
00230
00231
00232 #ifndef HOOMD_DEBUG
00233 #define HOOMD_DEBUG 0
00234 #endif
00235 #if HOOMD_DEBUG
00236 #define PRINTERR vmdcon_printf(VMDCON_ERROR, \
00237 "\n In file %s, line %d: \n %s \n \n", \
00238 __FILE__, __LINE__, strerror(errno))
00239 #else
00240 #define PRINTERR (void)(0)
00241 #endif
00242
00243
00244 #define SAFE_FREE(ptr) free(ptr); ptr=NULL
00245
00246 #define SAFE_CALLOC(ptr,type,count) \
00247 ptr = (type *)calloc(count,sizeof(type)); \
00248 if (ptr == NULL) { \
00249 \
00250 PRINTERR; \
00251 return MOLFILE_ERROR; \
00252 }
00253
00255 static void make_lower(char *s)
00256 {
00257 while (*s) {
00258 *s = tolower(*s);
00259 ++s;
00260 }
00261 }
00262
00265 static int addstring(char ***list, const char *key)
00266 {
00267 int i=0;
00268 char **p = *list;
00269
00270 for (; p[i] != NULL;) {
00271 if (strcmp(p[i], key) == 0)
00272 return i;
00273 ++i;
00274 }
00275
00276 *list = realloc(*list, (i+2)*sizeof(char *));
00277 (*list)[i] = strdup(key);
00278 (*list)[i+1] = NULL;
00279
00280 return i;
00281 }
00282
00283
00285 #define XML_ATTR_COPY(a,n,k,v) \
00286 strncpy(k, a[n], sizeof(k)-1); make_lower(k); \
00287 strncpy(v, a[n+1], sizeof(v)-1); make_lower(value)
00288
00290 #define XML_UNK_ATTR(e,k,v) \
00291 vmdcon_printf(VMDCON_WARN, "hoomdplugin) ignoring unknown " \
00292 " attribute %s='%s' for element %s\n", k, v, e);
00293
00294
00307 static void xml_new_tag(void *data,
00308 const XML_Char *elmt,
00309 const XML_Char **attr)
00310 {
00311 int i, mytag;
00312 hoomd_data_t *d=data;
00313 char element[MOLFILE_BUFSIZ];
00314 char key[MOLFILE_BUFSIZ];
00315 char value[MOLFILE_BUFSIZ];
00316
00317
00318 if (d == NULL) return;
00319
00320 d->parsedepth ++;
00321
00322 #if HOOMD_DEBUG
00323 for (i=0; i < d->parsedepth; ++i) {
00324 printf(">");
00325 }
00326 printf("%s",elmt);
00327 for (i=0;attr[i]; i+=2) {
00328 printf(" %s='%s'",attr[i],attr[i+1]);
00329 }
00330 printf("\n");
00331 #endif
00332
00333 strncpy(element, elmt, sizeof(element)-1);
00334 make_lower(element);
00335
00336
00337 mytag = HOOMD_UNKNOWN;
00338 for (i=0; i < HOOMD_NUMTAGS; ++i) {
00339 if (strcmp(elmt, hoomd_tag_name[i]) == 0) {
00340 mytag = i;
00341 break;
00342 }
00343 }
00344 d->currtag[d->parsedepth] = mytag;
00345 d->counter = 0;
00346
00347 switch (mytag) {
00348
00349 case HOOMD_XML:
00350 for (i=0; attr[i]; i+=2) {
00351 XML_ATTR_COPY(attr, i, key, value);
00352 if(strcmp(key,"version") == 0) {
00353 d->majv = atoi(strtok(value, "."));
00354 d->minv = atoi(strtok(NULL, "."));
00355 } else {
00356 XML_UNK_ATTR(element, key, value);
00357 }
00358 }
00359 if (d->majv < 0) {
00360 vmdcon_printf(VMDCON_INFO, "hoomdplugin) No version attribute found "
00361 "on <%s> element assuming version=\"1.0\"\n",
00362 hoomd_tag_name[HOOMD_XML] );
00363 d->majv = 1;
00364 d->minv = 0;
00365 }
00366 break;
00367
00368 case HOOMD_CONFIG:
00369 for (i=0; attr[i]; i+=2) {
00370 XML_ATTR_COPY(attr, i, key, value);
00371 if(strcmp(key,"time_step") == 0) {
00372 #if vmdplugin_ABIVERSION > 10
00373 d->ts.physical_time = atof(value);
00374 #else
00375 ;
00376 #endif
00377 } else if(strcmp(key,"natoms") == 0) {
00378 d->numatoms = atoi(value);
00379 } else if(strcmp(key,"ntypes") == 0) {
00380 d->numtypes = atoi(value);
00381 } else if(strcmp(key,"dimensions") == 0) {
00382 if (d->majv == 1 && d->minv < 2) {
00383 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found dimensions "
00384 "attribute in a pre-1.2 version format file. "
00385 "Ignoring it...\n");
00386 } else {
00387 d->numdims = atoi(value);
00388 }
00389 } else if(strcmp(key,"vizsigma") == 0) {
00390 if (d->majv == 1 && d->minv < 4) {
00391 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found vizsigma "
00392 "attribute in a pre-1.4 version format file. "
00393 "Ignoring it...\n");
00394 } else {
00395 d->mysigma = atof(value);
00396 }
00397 } else {
00398 XML_UNK_ATTR(element, key, value);
00399 }
00400 }
00401 break;
00402
00403 case HOOMD_BOX:
00404
00405 d->ts.A = d->ts.B = d->ts.C = 0.0f;
00406 d->ts.alpha = d->ts.beta = d->ts.gamma = 90.0f;
00407 for (i=0; attr[i]; i+=2) {
00408 XML_ATTR_COPY(attr, i, key, value);
00409 if (strcmp(key,"units") == 0) {
00410 } else if(strcmp(key,"lx") == 0) {
00411 d->ts.A = atof(value);
00412 } else if(strcmp(key,"ly") == 0) {
00413 d->ts.B = atof(value);
00414 } else if(strcmp(key,"lz") == 0) {
00415 d->ts.C = atof(value);
00416 } else {
00417 XML_UNK_ATTR(element, key, value);
00418 }
00419 }
00420 break;
00421
00422 case HOOMD_POS:
00423 for (i=0; attr[i]; i+=2) {
00424 XML_ATTR_COPY(attr, i, key, value);
00425 if (strcmp(key,"units") == 0) {
00426 ;
00427 } else if (strcmp(key,"num") == 0) {
00428 if (d->numatoms < 1) {
00429 d->numatoms = atoi(value);
00430 }
00431 } else {
00432 XML_UNK_ATTR(element, key, value);
00433 }
00434 }
00435 break;
00436
00437 case HOOMD_VEL:
00438 #if vmdplugin_ABIVERSION > 10
00439 d->ts_meta.has_velocities = 1;
00440 for (i=0; attr[i]; i+=2) {
00441 XML_ATTR_COPY(attr, i, key, value);
00442 if (strcmp(key, "units") == 0) {
00443 ;
00444 } else if (strcmp(key,"num") == 0) {
00445 ;
00446 } else {
00447 XML_UNK_ATTR(element, key, value);
00448 }
00449 }
00450 #endif
00451 break;
00452
00453 case HOOMD_BODY:
00454 if (d->majv == 1 && d->minv < 4) {
00455 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found <%s> section in "
00456 "a pre-1.4 version format file. Ignoring it...\n",
00457 hoomd_tag_name[mytag]);
00458 } else {
00459 for (i=0; attr[i]; i+=2) {
00460 XML_ATTR_COPY(attr, i, key, value);
00461 if (strcmp(key,"num") == 0) {
00462 ;
00463 } else {
00464 XML_UNK_ATTR(element, key, value);
00465 }
00466 }
00467 }
00468 break;
00469
00470 case HOOMD_MASS:
00471 d->optflags |= MOLFILE_MASS;
00472 for (i=0; attr[i]; i+=2) {
00473 XML_ATTR_COPY(attr, i, key, value);
00474 if (strcmp(key,"units") == 0) {
00475 ;
00476 } else if (strcmp(key,"num") == 0) {
00477 ;
00478 } else {
00479 XML_UNK_ATTR(element, key, value);
00480 }
00481 }
00482 break;
00483
00484 case HOOMD_CHARGE:
00485 if (d->majv == 1 && d->minv < 3) {
00486 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found charge "
00487 "section in a pre-1.3 version format file. "
00488 "Ignoring it...\n");
00489 } else {
00490 d->optflags |= MOLFILE_CHARGE;
00491 for (i=0; attr[i]; i+=2) {
00492 XML_ATTR_COPY(attr, i, key, value);
00493 if (strcmp(key,"units") == 0) {
00494 ;
00495 } else if (strcmp(key,"num") == 0) {
00496 ;
00497 } else {
00498 XML_UNK_ATTR(element, key, value);
00499 }
00500 }
00501 }
00502 break;
00503
00504 case HOOMD_DIAMETER:
00505 d->optflags |= MOLFILE_RADIUS;
00506 for (i=0; attr[i]; i+=2) {
00507 XML_ATTR_COPY(attr, i, key, value);
00508 if (strcmp(key,"units") == 0) {
00509 ;
00510 } else if (strcmp(key,"num") == 0) {
00511 ;
00512 } else {
00513 XML_UNK_ATTR(element, key, value);
00514 }
00515 }
00516 break;
00517
00518 case HOOMD_IMAGE:
00519 case HOOMD_TYPE:
00520 case HOOMD_BOND:
00521 for (i=0; attr[i]; i+=2) {
00522 XML_ATTR_COPY(attr, i, key, value);
00523 if (strcmp(key,"num") == 0) {
00524 ;
00525 } else {
00526 XML_UNK_ATTR(element, key, value);
00527 }
00528 }
00529 break;
00530
00531
00532 case HOOMD_ANGLE:
00533 case HOOMD_DIHEDRAL:
00534 case HOOMD_IMPROPER:
00535 if (d->majv == 1 && d->minv < 1) {
00536 if ( ((mytag == HOOMD_ANGLE) && (d->numangles < -1)) ||
00537 ((mytag == HOOMD_DIHEDRAL) && (d->numdihedrals < -1)) ||
00538 ((mytag == HOOMD_IMPROPER) && (d->numimpropers < -1)) ) {
00539 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found <%s> section in "
00540 "a pre-1.1 version format file. Ignoring it...\n",
00541 hoomd_tag_name[mytag]);
00542 } else {
00543 for (i=0; attr[i]; i+=2) {
00544 XML_ATTR_COPY(attr, i, key, value);
00545 if (strcmp(key,"num") == 0) {
00546 ;
00547 } else {
00548 XML_UNK_ATTR(element, key, value);
00549 }
00550 }
00551 }
00552 }
00553 break;
00554
00555 case HOOMD_WALL:
00556 for (i=0; attr[i]; i+=2) {
00557 XML_ATTR_COPY(attr, i, key, value);
00558 XML_UNK_ATTR(element, key, value);
00559 }
00560 break;
00561
00562 case HOOMD_COORD:
00563 for (i=0; attr[i]; i+=2) {
00564 XML_ATTR_COPY(attr, i, key, value);
00565 if (strcmp(key,"units") == 0) {
00566 ;
00567 } else if(strcmp(key,"ox") == 0) {
00568 ;
00569 } else if(strcmp(key,"oy") == 0) {
00570 ;
00571 } else if(strcmp(key,"oz") == 0) {
00572 ;
00573 } else if(strcmp(key,"nx") == 0) {
00574 ;
00575 } else if(strcmp(key,"ny") == 0) {
00576 ;
00577 } else if(strcmp(key,"nz") == 0) {
00578 ;
00579 } else {
00580 XML_UNK_ATTR(element, key, value);
00581 }
00582 }
00583 break;
00584
00585 case HOOMD_ACCEL:
00586 if (d->majv == 1 && d->minv < 2) {
00587 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found <%s> section in "
00588 "a pre-1.2 version format file. Ignoring it...\n",
00589 hoomd_tag_name[mytag]);
00590 } else {
00591 for (i=0; attr[i]; i+=2) {
00592 XML_ATTR_COPY(attr, i, key, value);
00593 if (strcmp(key,"num") == 0) {
00594 ;
00595 } else {
00596 XML_UNK_ATTR(element, key, value);
00597 }
00598 XML_UNK_ATTR(element, key, value);
00599 }
00600 }
00601 break;
00602
00603 case HOOMD_ORIENT:
00604 if (d->majv == 1 && d->minv < 4) {
00605 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Found <%s> section in "
00606 "a pre-1.4 version format file. Ignoring it...\n",
00607 hoomd_tag_name[mytag]);
00608 } else {
00609 for (i=0; attr[i]; i+=2) {
00610 XML_ATTR_COPY(attr, i, key, value);
00611 if (strcmp(key,"num") == 0) {
00612 ;
00613 } else {
00614 XML_UNK_ATTR(element, key, value);
00615 }
00616 }
00617 }
00618 break;
00619
00620 default:
00621 d->currtag[d->parsedepth] = HOOMD_UNKNOWN;
00622 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Ignoring unknown XML "
00623 "element: %s.\n", element);
00624 break;
00625 }
00626 }
00627
00629 static void xml_end_tag(void *data, const XML_Char *elmt)
00630 {
00631 hoomd_data_t *d=data;
00632 int mytag;
00633
00634 if (d == NULL) return;
00635
00636 mytag = d->currtag[d->parsedepth];
00637 switch (mytag) {
00638
00639 case HOOMD_CONFIG:
00640 d->doneframe=1;
00641 break;
00642
00643 case HOOMD_TYPE:
00644
00645
00646
00647
00648 if (d->numatoms < 0) {
00649 d->numatoms = d->counter;
00650 }
00651 break;
00652
00653 case HOOMD_XML:
00654 case HOOMD_BOX:
00655 case HOOMD_IMAGE:
00656 case HOOMD_WALL:
00657 case HOOMD_COORD:
00658 case HOOMD_ACCEL:
00659 case HOOMD_ORIENT:
00660 ;
00661 break;
00662
00663 case HOOMD_POS:
00664
00665 if ( (d->counter > 0) && (d->numatoms < 0) ) {
00666 d->numatoms = d->counter/3;
00667 }
00668
00669 case HOOMD_VEL:
00670 if ( (d->counter > 0) && (d->numatoms > 0) && (d->counter != (3*d->numatoms)) ) {
00671 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Inconsistent %s data. Expected %d, but "
00672 "got %d items.\n", hoomd_tag_name[mytag], 3*d->numatoms, d->counter);
00673 d->parse_error = 1;
00674 }
00675 break;
00676
00677 case HOOMD_MASS:
00678 case HOOMD_CHARGE:
00679 case HOOMD_BODY:
00680 case HOOMD_DIAMETER:
00681 if ( (d->counter > 0) && (d->numatoms > 0) && (d->counter != d->numatoms) ) {
00682 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Inconsistent %s data. Expected %d, but "
00683 "got %d items.\n", hoomd_tag_name[mytag], d->numatoms, d->counter);
00684 d->parse_error = 1;
00685 }
00686
00687 break;
00688
00689 case HOOMD_BOND:
00690
00691
00692
00693 if (d->numbonds < 0) {
00694 d->numbonds = d->counter/3;
00695 }
00696 break;
00697
00698 case HOOMD_ANGLE:
00699
00700
00701
00702 if (d->numangles < 0) {
00703 if (d->majv == 1 && d->minv > 0) {
00704 d->numangles = d->counter/4;
00705 } else {
00706 d->numangles = 0;
00707 }
00708 }
00709 break;
00710
00711 case HOOMD_DIHEDRAL:
00712
00713
00714
00715 if (d->numdihedrals < 0) {
00716 if (d->majv == 1 && d->minv > 0) {
00717 d->numdihedrals = d->counter/5;
00718 } else {
00719 d->numdihedrals = 0;
00720 }
00721 }
00722 break;
00723
00724 case HOOMD_IMPROPER:
00725
00726
00727
00728 if (d->numimpropers < 0) {
00729 if (d->majv == 1 && d->minv > 0) {
00730 d->numimpropers = d->counter/5;
00731 } else {
00732 d->numimpropers = 0;
00733 }
00734 }
00735 break;
00736
00737 default:
00738 if (mytag < HOOMD_NUMTAGS) {
00739 vmdcon_printf(VMDCON_WARN, "hoomdplugin) No end handler for HOOMD tag '%s'.\n",
00740 hoomd_tag_name[mytag]);
00741 } else {
00742 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Unknown HOOMD tag id: '%d'.\n", mytag);
00743 }
00744
00745 break;
00746 }
00747
00748 d->currtag[d->parsedepth] = HOOMD_NONE;
00749 d->parsedepth--;
00750 d->counter=0;
00751 #if HOOMD_DEBUG
00752 printf("end of tag %s. parsedepth=%d\n",elmt,d->parsedepth);
00753 #endif
00754 }
00755
00757 static void xml_data_block(void *data, const XML_Char *s, int len)
00758 {
00759 hoomd_data_t *d =NULL;
00760 int i,lmax, mytag;
00761
00762 d = (hoomd_data_t *)data;
00763 if (d == NULL) return;
00764 if (d->parsedepth < 1) return;
00765 if (d->parsedepth > HOOMD_MAXDEPTH) return;
00766 if (len < 1) return;
00767
00768 mytag = d->currtag[d->parsedepth];
00769 switch (mytag) {
00770
00771 case HOOMD_TYPE:
00772
00773
00774 if (d->numatoms < 0) {
00775 char buffer[1024];
00776 char *p;
00777 lmax=1023;
00778
00779 if (len < lmax) lmax=len;
00780 memcpy(buffer,s,lmax);
00781 buffer[lmax]='\0';
00782
00783 p=strtok(buffer," \t\n");
00784 while (p) {
00785 d->counter ++;
00786 p=strtok(NULL," \t\n");
00787 }
00788 } else {
00789
00790 if (d->numframe < 2) {
00791 char buffer[1024];
00792 char *p;
00793 molfile_atom_t *atom;
00794
00795 if (len==1 && ((*s == ' ') || (*s == '\n') || (*s == '\t'))) return;
00796
00797 lmax=1023;
00798 if (len < lmax) lmax=len;
00799 memcpy(buffer,s,lmax);
00800 buffer[lmax]='\0';
00801
00802 p=strtok(buffer," \t\n");
00803 while (p) {
00804 atom=d->atomlist + d->counter;
00805 if (atom == NULL) return;
00806
00807 strncpy(atom->name, p, sizeof(atom->name));
00808 strncpy(atom->type, p, sizeof(atom->type));
00809 atom->atomicnumber = get_pte_idx(atom->type);
00810 d->counter ++;
00811 p=strtok(NULL," \t\n");
00812 }
00813 }
00814 }
00815 break;
00816
00817 case HOOMD_BOND:
00818
00819 if ((d->numbonds < 0) || (d->numframe < 2)) {
00820 char buffer[1024];
00821 char *p;
00822 int idx, num, n;
00823 lmax=1023;
00824
00825 if (len < lmax) lmax=len;
00826 memcpy(buffer,s,lmax);
00827 buffer[lmax]='\0';
00828
00829 p=strtok(buffer," \t\n");
00830 while (p) {
00831 num = d->counter / 3;
00832 if (d->numbonds > 0) {
00833 i=d->counter % 3;
00834 if (i == 0) {
00835 n = d->numbondtypes;
00836 idx=addstring(&(d->bondtypename), p);
00837 if (idx < 0)
00838 d->parse_error = 1;
00839 if (idx >= n)
00840 d->numbondtypes=idx+1;
00841 d->bondtype[num] = idx;
00842
00843
00844 } else if (i == 1) {
00845 d->from[num] = atoi(p) + 1;
00846 } else if (i == 2) {
00847 d->to[num] = atoi(p) + 1;
00848 }
00849 }
00850 d->counter ++;
00851 p=strtok(NULL," \t\n");
00852 }
00853 }
00854 break;
00855
00856 case HOOMD_ANGLE:
00857 if (d->majv == 1 && d->minv > 0) {
00858
00859 if ((d->numangles < 0) || (d->numframe < 2)) {
00860 char buffer[1024];
00861 char *p;
00862 int idx, num, n;
00863 lmax=1023;
00864
00865 if (len < lmax) lmax=len;
00866 memcpy(buffer,s,lmax);
00867 buffer[lmax]='\0';
00868
00869 p=strtok(buffer," \t\n");
00870 while (p) {
00871 num = d->counter / 4;
00872 if (d->numangles > 0) {
00873 i=d->counter % 4;
00874 if (i == 0) {
00875 n = d->numangletypes;
00876 idx=addstring(&(d->angletypename), p);
00877 if (idx < 0)
00878 d->parse_error = 1;
00879 if (idx >= n)
00880 d->numangletypes=idx+1;
00881 d->angletype[num] = idx;
00882
00883
00884 } else if (i == 1) {
00885 d->angle[3*num ] = atoi(p) + 1;
00886 } else if (i == 2) {
00887 d->angle[3*num+1] = atoi(p) + 1;
00888 } else if (i == 3) {
00889 d->angle[3*num+2] = atoi(p) + 1;
00890 }
00891 }
00892 d->counter ++;
00893 p=strtok(NULL," \t\n");
00894 }
00895 }
00896 }
00897 break;
00898
00899 case HOOMD_DIHEDRAL:
00900 if (d->majv == 1 && d->minv > 0) {
00901
00902 if ((d->numdihedrals < 0) || (d->numframe < 2)) {
00903 char buffer[1024];
00904 char *p;
00905 int idx, num, n;
00906 lmax=1023;
00907
00908 if (len < lmax) lmax=len;
00909 memcpy(buffer,s,lmax);
00910 buffer[lmax]='\0';
00911
00912 p=strtok(buffer," \t\n");
00913 while (p) {
00914 num = d->counter / 5;
00915 if (d->numdihedrals > 0) {
00916 i=d->counter % 5;
00917 if (i == 0) {
00918 n = d->numdihedraltypes;
00919 idx=addstring(&(d->dihedraltypename), p);
00920 if (idx < 0)
00921 d->parse_error = 1;
00922 if (idx >= n)
00923 d->numdihedraltypes=idx+1;
00924 d->dihedraltype[num] = idx;
00925
00926
00927 } else if (i == 1) {
00928 d->dihedral[4*num ] = atoi(p) + 1;
00929 } else if (i == 2) {
00930 d->dihedral[4*num+1] = atoi(p) + 1;
00931 } else if (i == 3) {
00932 d->dihedral[4*num+2] = atoi(p) + 1;
00933 } else if (i == 4) {
00934 d->dihedral[4*num+3] = atoi(p) + 1;
00935 }
00936 }
00937 d->counter ++;
00938 p=strtok(NULL," \t\n");
00939 }
00940 }
00941 }
00942 break;
00943
00944 case HOOMD_IMPROPER:
00945 if (d->majv == 1 && d->minv > 0) {
00946
00947 if ((d->numimpropers < 0) || (d->numframe < 2)) {
00948 char buffer[1024];
00949 char *p;
00950 int idx, num, n;
00951 lmax=1023;
00952
00953 if (len < lmax) lmax=len;
00954 memcpy(buffer,s,lmax);
00955 buffer[lmax]='\0';
00956
00957 p=strtok(buffer," \t\n");
00958 while (p) {
00959 num = d->counter / 5;
00960 if (d->numimpropers > 0) {
00961 i=d->counter % 5;
00962 if (i == 0) {
00963 n = d->numimpropertypes;
00964 idx=addstring(&(d->impropertypename), p);
00965 if (idx < 0)
00966 d->parse_error = 1;
00967 if (idx >= n)
00968 d->numimpropertypes=idx+1;
00969 d->impropertype[num] = idx;
00970
00971
00972 } else if (i == 1) {
00973 d->improper[4*num ] = atoi(p) + 1;
00974 } else if (i == 2) {
00975 d->improper[4*num+1] = atoi(p) + 1;
00976 } else if (i == 3) {
00977 d->improper[4*num+2] = atoi(p) + 1;
00978 } else if (i == 4) {
00979 d->improper[4*num+3] = atoi(p) + 1;
00980 }
00981 }
00982 d->counter ++;
00983 p=strtok(NULL," \t\n");
00984 }
00985 }
00986 }
00987 break;
00988
00989 case HOOMD_DIAMETER:
00990
00991 if (d->numatoms > 0) {
00992 char buffer[1024];
00993 char *p;
00994 molfile_atom_t *atom;
00995
00996 if (len==1 && ((*s == ' ') || (*s == '\n') || (*s == '\t'))) return;
00997
00998 lmax=1023;
00999 if (len < lmax) lmax=len;
01000 memcpy(buffer,s,lmax);
01001 buffer[lmax]='\0';
01002
01003 p=strtok(buffer," \t\n");
01004 while (p) {
01005 atom=d->atomlist + d->counter;
01006 if (atom == NULL) return;
01007
01008 atom->radius = 0.5 * atof(p);
01009 d->counter ++;
01010 p=strtok(NULL," \t\n");
01011 }
01012 }
01013 break;
01014
01015 case HOOMD_MASS:
01016
01017 if (d->numatoms > 0) {
01018 char buffer[1024];
01019 char *p;
01020 molfile_atom_t *atom;
01021
01022 if (len==1 && ((*s == ' ') || (*s == '\n') || (*s == '\t'))) return;
01023
01024 lmax=1023;
01025 if (len < lmax) lmax=len;
01026 memcpy(buffer,s,lmax);
01027 buffer[lmax]='\0';
01028
01029 p=strtok(buffer," \t\n");
01030 while (p) {
01031 atom=d->atomlist + d->counter;
01032 if (atom == NULL) return;
01033
01034 atom->mass = atof(p);
01035 d->counter ++;
01036 p=strtok(NULL," \t\n");
01037 }
01038 }
01039 break;
01040
01041 case HOOMD_BODY:
01042
01043 if (d->numatoms > 0) {
01044 char buffer[1024];
01045 char *p;
01046 molfile_atom_t *atom;
01047
01048 if (len==1 && ((*s == ' ') || (*s == '\n') || (*s == '\t'))) return;
01049
01050 lmax=1023;
01051 if (len < lmax) lmax=len;
01052 memcpy(buffer,s,lmax);
01053 buffer[lmax]='\0';
01054
01055 p=strtok(buffer," \t\n");
01056 while (p) {
01057 atom=d->atomlist + d->counter;
01058 if (atom == NULL) return;
01059
01060 atom->resid = atoi(p);
01061 d->counter ++;
01062 p=strtok(NULL," \t\n");
01063 }
01064 }
01065 break;
01066
01067 case HOOMD_CHARGE:
01068
01069 if (d->majv == 1 && d->minv > 2) {
01070 if (d->numatoms > 0) {
01071 char buffer[1024];
01072 char *p;
01073 molfile_atom_t *atom;
01074
01075 if (len==1 && ((*s == ' ') || (*s == '\n') || (*s == '\t'))) return;
01076
01077 lmax=1023;
01078 if (len < lmax) lmax=len;
01079 memcpy(buffer,s,lmax);
01080 buffer[lmax]='\0';
01081
01082 p=strtok(buffer," \t\n");
01083 while (p) {
01084 atom=d->atomlist + d->counter;
01085 if (atom == NULL) return;
01086
01087 atom->charge = atof(p);
01088 d->counter ++;
01089 p=strtok(NULL," \t\n");
01090 }
01091 }
01092 }
01093 break;
01094
01095 case HOOMD_POS:
01096
01097 if (d->ts.coords) {
01098 char buffer[1024];
01099 int lmax=1023;
01100 char *p;
01101
01102 if (len < lmax) lmax=len;
01103 memcpy(buffer,s,lmax);
01104 buffer[lmax]='\0';
01105
01106 p=strtok(buffer," \t\n");
01107 while (p) {
01108 if (d->counter < (3 * d->numatoms))
01109 d->ts.coords[d->counter] = atof(p);
01110 d->counter ++;
01111 p=strtok(NULL," \t\n");
01112 }
01113
01114 }
01115 break;
01116
01117 #if vmdplugin_ABIVERSION > 10
01118 case HOOMD_VEL:
01119 if ((d->numatoms > 0) && (d->ts.velocities != NULL)) {
01120 char buffer[1024];
01121 int lmax=1023;
01122 char *p;
01123
01124 if (len < lmax) lmax=len;
01125 memcpy(buffer,s,lmax);
01126 buffer[lmax]='\0';
01127
01128 p=strtok(buffer," \t\n");
01129 while (p) {
01130 if (d->counter < (3 * d->numatoms))
01131 d->ts.velocities[d->counter] = atof(p);
01132 d->counter ++;
01133 p=strtok(NULL," \t\n");
01134 }
01135 }
01136 break;
01137 #endif
01138
01139 case HOOMD_IMAGE:
01140 if ((d->numatoms > 0) && (d->imagecounts != NULL)) {
01141 char buffer[1024];
01142 int lmax=1023;
01143 char *p;
01144
01145 if (len < lmax) lmax=len;
01146 memcpy(buffer,s,lmax);
01147 buffer[lmax]='\0';
01148
01149 p=strtok(buffer," \t\n");
01150 while (p) {
01151 d->imagecounts[d->counter] = atoi(p);
01152 d->counter ++;
01153 p=strtok(NULL," \t\n");
01154 }
01155 }
01156 break;
01157
01158 default:
01159 break;
01160 }
01161 }
01162
01164 static void xml_comment(void *data, const XML_Char *s)
01165 {
01166 hoomd_data_t *d=data;
01167 if (d==NULL) return;
01168 #if HOOMD_DEBUG
01169 printf("COMMENT: %s. parsedepth=%d\n", (char *)s, d->parsedepth);
01170 #endif
01171 }
01172
01174 static int hoomd_parse_line(hoomd_data_t *data)
01175 {
01176 int done, len;
01177
01178 data->buffer = XML_GetBuffer(data->p, MOLFILE_BIGBUFSIZ);
01179 done = (NULL == fgets(data->buffer,MOLFILE_BIGBUFSIZ,data->fp));
01180
01181 if (!done)
01182 len=strlen(data->buffer);
01183 else
01184 len=0;
01185
01186 if (ferror(data->fp)) {
01187 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) problem reading HOOMD"
01188 " data file '%s'\n", data->filename);
01189 return MOLFILE_ERROR;
01190 }
01191
01192 if (! XML_ParseBuffer(data->p, len, done)) {
01193 vmdcon_printf(VMDCON_ERROR,
01194 "hoomdplugin) XML syntax error at line %d:\n%s\n",
01195 XML_GetCurrentLineNumber(data->p),
01196 XML_ErrorString(XML_GetErrorCode(data->p)));
01197 return MOLFILE_ERROR;
01198 }
01199
01200 if (data->parse_error > 0) {
01201 vmdcon_printf(VMDCON_ERROR,
01202 "hoomdplugin) XML data parse error at line %d.\n",
01203 XML_GetCurrentLineNumber(data->p));
01204 return MOLFILE_ERROR;
01205 }
01206
01207 return MOLFILE_SUCCESS;
01208 }
01209
01210
01215 static void *open_hoomd_read(const char *filename, const char *filetype,
01216 int *natoms) {
01217 FILE *fp;
01218 XML_Parser p;
01219 hoomd_data_t *data;
01220
01221 fp = fopen(filename, "rb");
01222 if (!fp) return NULL;
01223
01224 data = (hoomd_data_t *)calloc(1,sizeof(hoomd_data_t));
01225 if (data) {
01226
01227 data->counter = 0;
01228 data->numatoms = -1;
01229 data->numtypes = -1;
01230 data->numbonds = -1;
01231 data->numangles = -1;
01232 data->numdihedrals = -1;
01233 data->numimpropers = -1;
01234 data->numframe = -1;
01235 data->numbondtypes = -1;
01236 data->numangletypes = -1;
01237 data->numdihedraltypes = -1;
01238 data->numimpropertypes = -1;
01239 data->bondtype = NULL;
01240 data->bondtypename = NULL;
01241 data->angle = NULL;
01242 data->angletype = NULL;
01243 data->angletypename = NULL;
01244 data->dihedral = NULL;
01245 data->dihedraltype = NULL;
01246 data->dihedraltypename = NULL;
01247 data->dihedral = NULL;
01248 data->dihedraltype = NULL;
01249 data->dihedraltypename = NULL;
01250 data->parse_error = 0;
01251
01252 data->majv = -1;
01253 data->minv = -1;
01254
01255 data->optflags = MOLFILE_NOOPTIONS;
01256
01257 data->mysigma = 1.0f;
01258
01259 p = XML_ParserCreate(NULL);
01260 if (!p) {
01261 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Could not create XML"
01262 " parser for HOOMD-blue data file '%s'\n", filename);
01263 SAFE_FREE(data);
01264 fclose(fp);
01265 return NULL;
01266 }
01267
01268 XML_SetElementHandler(p, xml_new_tag, xml_end_tag);
01269 XML_SetCommentHandler(p, xml_comment);
01270 XML_SetCharacterDataHandler(p, xml_data_block);
01271 XML_SetUserData(p,data);
01272 data->p = p;
01273 data->fp = fp;
01274 data->filename = strdup(filename);
01275
01276
01277 do {
01278 if (MOLFILE_ERROR == hoomd_parse_line(data)) {
01279 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) XML Parse error "
01280 "while reading HOOMD-blue data file '%s'\n", filename);
01281 XML_ParserFree(data->p);
01282 data->p = NULL;
01283 SAFE_FREE(data->filename);
01284 fclose(data->fp);
01285 SAFE_FREE(data);
01286 return NULL;
01287 }
01288 } while (!feof(fp) && !data->doneframe);
01289
01290 if ( data->majv > HOOMD_FORMAT_MAJV ) {
01291 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Encountered incompatible "
01292 "HOOMD-blue data file format version '%d.%d.\n",
01293 data->majv, data->minv);
01294 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) This plugin supports only "
01295 "HOOMD-blue data files up to version '%d.%d'.\n",
01296 HOOMD_FORMAT_MAJV, HOOMD_FORMAT_MINV);
01297 XML_ParserFree(data->p);
01298 data->p = NULL;
01299 SAFE_FREE(data->filename);
01300 fclose(data->fp);
01301 SAFE_FREE(data);
01302 return NULL;
01303 } else {
01304 if ( (data->majv == HOOMD_FORMAT_MAJV) &&
01305 (data->minv > HOOMD_FORMAT_MINV) ) {
01306 vmdcon_printf(VMDCON_WARN, "hoomdplugin) Encountered newer HOOMD-blue "
01307 "data file format version '%d.%d'.\n",
01308 data->majv, data->minv);
01309 vmdcon_printf(VMDCON_WARN, "hoomdplugin) This plugin supports HOOMD-blue "
01310 "data files up to version '%d.%d'. Continuing...\n",
01311 HOOMD_FORMAT_MAJV, HOOMD_FORMAT_MINV);
01312 }
01313 }
01314
01315 if (data->numatoms < 0) {
01316 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Could not determine "
01317 "number of atoms in HOOMD-blue data file '%s'\n", filename);
01318 XML_ParserFree(data->p);
01319 data->p = NULL;
01320 SAFE_FREE(data->filename);
01321 fclose(data->fp);
01322 SAFE_FREE(data);
01323 return NULL;
01324 }
01325
01326
01327 XML_ParserFree(p);
01328 data->p = NULL;
01329 rewind(fp);
01330
01331 data->counter = 0;
01332 data->numframe = 0;
01333 *natoms=data->numatoms;
01334 }
01335
01336 return data;
01337 }
01338
01339 static int read_hoomd_structure(void *mydata, int *optflags,
01340 molfile_atom_t *atoms) {
01341 molfile_atom_t *a;
01342 XML_Parser p;
01343 const char *envvar;
01344 int i;
01345 hoomd_data_t *data = (hoomd_data_t *)mydata;
01346
01347 data->parsedepth = 0;
01348 data->counter = 0;
01349 data->numframe = 0;
01350 data->doneframe = 0;
01351 data->atomlist = atoms;
01352 SAFE_CALLOC(data->imagecounts, int, 3*data->numatoms);
01353 SAFE_CALLOC(data->bondtypename, char*, 1);
01354 SAFE_CALLOC(data->angletypename, char*, 1);
01355 SAFE_CALLOC(data->dihedraltypename, char*, 1);
01356 SAFE_CALLOC(data->impropertypename, char*, 1);
01357 if (data->numbonds > 0) {
01358 SAFE_CALLOC(data->from, int, data->numbonds);
01359 SAFE_CALLOC(data->to, int, data->numbonds);
01360 SAFE_CALLOC(data->bondtype, int, data->numbonds);
01361 }
01362 if (data->numangles > 0) {
01363 SAFE_CALLOC(data->angle, int, data->numangles * 3);
01364 SAFE_CALLOC(data->angletype, int, data->numangles);
01365 }
01366 if (data->numdihedrals > 0) {
01367 SAFE_CALLOC(data->dihedral, int, data->numdihedrals * 4);
01368 SAFE_CALLOC(data->dihedraltype, int, data->numdihedrals);
01369 }
01370 if (data->numimpropers > 0) {
01371 SAFE_CALLOC(data->improper, int, data->numimpropers * 4);
01372 SAFE_CALLOC(data->impropertype, int, data->numimpropers);
01373 }
01374 SAFE_CALLOC(data->ts.coords,float,3*data->numatoms);
01375 data->ts.A = data->ts.B = data->ts.C = 0.0f;
01376 data->ts.alpha = data->ts.beta = data->ts.gamma = 90.0f;
01377
01378 #if vmdplugin_ABIVERSION > 10
01379 data->ts_meta.count = -1;
01380 data->ts_meta.has_velocities = 0;
01381 SAFE_CALLOC(data->ts.velocities,float,3*data->numatoms);
01382 #endif
01383 p = XML_ParserCreate(NULL);
01384 if (!p) {
01385 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Could not create XML"
01386 " parser for HOOMD-blue data file '%s'\n", data->filename);
01387 return MOLFILE_ERROR;
01388 }
01389
01390 XML_SetElementHandler(p, xml_new_tag, xml_end_tag);
01391 XML_SetCommentHandler(p, xml_comment);
01392 XML_SetCharacterDataHandler(p,xml_data_block);
01393 XML_SetUserData(p,data);
01394 data->p = p;
01395
01396
01397 for (i=0, a=atoms; i < data->numatoms; ++i, ++a) {
01398 a->radius = 1.0f;
01399 a->mass = 1.0f;
01400 a->resname[0]='\0';
01401 a->resid=0;
01402 a->segid[0]='\0';
01403 a->chain[0]='\0';
01404 }
01405
01406
01407 do {
01408 if (MOLFILE_ERROR == hoomd_parse_line(data)) {
01409 XML_ParserFree(data->p);
01410 data->p = NULL;
01411 return MOLFILE_ERROR;
01412 }
01413 } while (!feof(data->fp) && !data->doneframe);
01414
01415
01416 envvar = getenv("VMDHOOMDSIGMA");
01417 if (envvar) data->mysigma = atof(envvar);
01418
01419
01420 for (i=0, a=atoms; i < data->numatoms; ++i, ++a) {
01421
01422 if (!(data->optflags & MOLFILE_RADIUS)) {
01423
01424
01425 a->radius = data->mysigma * get_pte_vdw_radius(a->atomicnumber);
01426 }
01427
01428 if (!(data->optflags & MOLFILE_MASS)) {
01429
01430 a->mass = get_pte_mass(a->atomicnumber);
01431 }
01432 }
01433 *optflags = data->optflags | MOLFILE_RADIUS | MOLFILE_MASS;
01434
01435 data->numframe=1;
01436
01437 return MOLFILE_SUCCESS;
01438 }
01439
01440 #if vmdplugin_ABIVERSION > 10
01441
01442 static int read_timestep_metadata(void *mydata,
01443 molfile_timestep_metadata_t *meta) {
01444 hoomd_data_t *data = (hoomd_data_t *)mydata;
01445
01446 meta->count = -1;
01447 meta->has_velocities = data->ts_meta.has_velocities;
01448 if (meta->has_velocities) {
01449 vmdcon_printf(VMDCON_INFO, "hoomdplugin) Importing velocities.\n");
01450 }
01451 return MOLFILE_SUCCESS;
01452 }
01453 #endif
01454
01455
01456
01457 #if vmdplugin_ABIVERSION > 14
01458 static int read_hoomd_bonds(void *v, int *nbonds, int **fromptr, int **toptr,
01459 float **bondorder, int **bondtype,
01460 int *nbondtypes, char ***bondtypename) {
01461 hoomd_data_t *data = (hoomd_data_t *)v;
01462
01463 *nbonds = data->numbonds;
01464 *bondorder = NULL;
01465
01466 if (data->numbonds > 0) {
01467 *fromptr = data->from;
01468 *toptr = data->to;
01469 *bondtype = data->bondtype;
01470 *nbondtypes = data->numbondtypes;
01471 *bondtypename = data->bondtypename;
01472 } else {
01473 *fromptr = NULL;
01474 *toptr = NULL;
01475 *bondtype = NULL;
01476 *nbondtypes = 0;
01477 *bondtypename = NULL;
01478 vmdcon_printf(VMDCON_WARN,
01479 "hoomdplugin) no bonds defined in data file.\n");
01480 }
01481 return MOLFILE_SUCCESS;
01482 }
01483 #else
01484 static int read_hoomd_bonds(void *v, int *nbonds, int **fromptr, int **toptr,
01485 float **bondorder) {
01486 hoomd_data_t *data = (hoomd_data_t *)v;
01487
01488 *nbonds = data->numbonds;
01489 *bondorder = NULL;
01490
01491 if (data->numbonds > 0) {
01492 *fromptr = data->from;
01493 *toptr = data->to;
01494 } else {
01495 *fromptr = NULL;
01496 *toptr = NULL;
01497 vmdcon_printf(VMDCON_WARN,
01498 "hoomdplugin) no bonds defined in data file.\n");
01499 }
01500 return MOLFILE_SUCCESS;
01501 }
01502 #endif
01503
01504 #if vmdplugin_ABIVERSION > 15
01505 static int read_hoomd_angles(void *v, int *numangles, int **angles,
01506 int **angletypes, int *numangletypes,
01507 char ***angletypenames, int *numdihedrals,
01508 int **dihedrals, int **dihedraltypes,
01509 int *numdihedraltypes, char ***dihedraltypenames,
01510 int *numimpropers, int **impropers,
01511 int **impropertypes, int *numimpropertypes,
01512 char ***impropertypenames, int *numcterms,
01513 int **cterms, int *ctermcols, int *ctermrows) {
01514
01515 hoomd_data_t *data = (hoomd_data_t *)v;
01516
01517 *numangles = 0;
01518 *angles = NULL;
01519 *angletypes = NULL;
01520 *numangletypes = 0;
01521 *angletypenames = NULL;
01522 *numdihedrals = 0;
01523 *dihedrals = NULL;
01524 *dihedraltypes = NULL;
01525 *numdihedraltypes = 0;
01526 *dihedraltypenames = NULL;
01527 *numimpropers = 0;
01528 *impropers = NULL;
01529 *impropertypes = NULL;
01530 *numimpropertypes = 0;
01531 *impropertypenames = NULL;
01532 *numcterms = 0;
01533 *cterms = NULL;
01534 *ctermrows = 0;
01535 *ctermcols = 0;
01536
01537 if (data->numangles > 0) {
01538 *numangles = data->numangles;
01539 *angles = data->angle;
01540 *angletypes = data->angletype;
01541 *numangletypes = data->numangletypes;
01542 *angletypenames = data->angletypename;
01543 } else {
01544 vmdcon_printf(VMDCON_INFO,
01545 "hoomdplugin) no angles defined in data file.\n");
01546 }
01547 if (data->numdihedrals > 0) {
01548 *numdihedrals = data->numdihedrals;
01549 *dihedrals = data->dihedral;
01550 *dihedraltypes = data->dihedraltype;
01551 *numdihedraltypes = data->numdihedraltypes;
01552 *dihedraltypenames = data->dihedraltypename;
01553 } else {
01554 vmdcon_printf(VMDCON_INFO,
01555 "hoomdplugin) no dihedrals defined in data file.\n");
01556 }
01557 if (data->numimpropers > 0) {
01558 *numimpropers = data->numimpropers;
01559 *impropers = data->improper;
01560 *impropertypes = data->impropertype;
01561 *numimpropertypes = data->numimpropertypes;
01562 *impropertypenames = data->impropertypename;
01563 } else {
01564 vmdcon_printf(VMDCON_INFO,
01565 "hoomdplugin) no impropers defined in data file.\n");
01566 }
01567 return MOLFILE_SUCCESS;
01568 }
01569 #else
01570 static int read_hoomd_angles(void *v, int *numangles, int **angles,
01571 double **angleforces, int *numdihedrals,
01572 int **dihedrals, double **dihedralforces,
01573 int *numimpropers, int **impropers,
01574 double **improperforces, int *numcterms,
01575 int **cterms, int *ctermcols, int *ctermrows,
01576 double **ctermforces) {
01577 hoomd_data_t *data = (hoomd_data_t *)v;
01578
01579 *numangles = 0;
01580 *angles = NULL;
01581 *angleforces = NULL;
01582 *numdihedrals = 0;
01583 *dihedrals = NULL;
01584 *dihedralforces = NULL;
01585 *numimpropers = 0;
01586 *impropers = NULL;
01587 *improperforces = NULL;
01588 *numcterms = 0;
01589 *cterms = NULL;
01590 *ctermrows = 0;
01591 *ctermcols = 0;
01592 *ctermforces = NULL;
01593
01594 if (data->numangles > 0) {
01595 *numangles = data->numangles;
01596 *angles = data->angle;
01597 } else {
01598 vmdcon_printf(VMDCON_INFO,
01599 "hoomdplugin) no angles defined in data file.\n");
01600 }
01601 if (data->numdihedrals > 0) {
01602 *numdihedrals = data->numdihedrals;
01603 *dihedrals = data->dihedral;
01604 } else {
01605 vmdcon_printf(VMDCON_INFO,
01606 "hoomdplugin) no dihedrals defined in data file.\n");
01607 }
01608 if (data->numimpropers > 0) {
01609 *numimpropers = data->numimpropers;
01610 *impropers = data->improper;
01611 } else {
01612 vmdcon_printf(VMDCON_INFO,
01613 "hoomdplugin) no impropers defined in data file.\n");
01614 }
01615 return MOLFILE_SUCCESS;
01616 }
01617 #endif
01618
01619 static int read_hoomd_timestep(void *mydata, int natoms,
01620 molfile_timestep_t *ts) {
01621 int i;
01622 hoomd_data_t *data;
01623
01624 data=(hoomd_data_t *)mydata;
01625
01626 if (data->parse_error) return MOLFILE_ERROR;
01627 if (data->p == NULL) return MOLFILE_ERROR;
01628
01629 if (data->numframe > 1) {
01630
01631 data->doneframe = 0;
01632 do {
01633 if (MOLFILE_ERROR == hoomd_parse_line(data)) {
01634 XML_ParserFree(data->p);
01635 data->p = NULL;
01636 return MOLFILE_ERROR;
01637 }
01638 if (data->parse_error) return MOLFILE_ERROR;
01639 } while (!feof(data->fp) && !data->doneframe);
01640 if (feof(data->fp)) return MOLFILE_ERROR;
01641 }
01642 data->numframe++;
01643
01644 if (ts != NULL) {
01645
01646
01647 for (i=0; i<natoms; ++i) {
01648 ts->coords[3*i+0] = data->ts.coords[3*i+0]
01649 + data->ts.A * data->imagecounts[3*i+0];
01650 ts->coords[3*i+1] = data->ts.coords[3*i+1]
01651 + data->ts.B * data->imagecounts[3*i+1];
01652 ts->coords[3*i+2] = data->ts.coords[3*i+2]
01653 + data->ts.C * data->imagecounts[3*i+2];
01654 }
01655 #if vmdplugin_ABIVERSION > 10
01656
01657 if (ts->velocities != NULL) {
01658 for (i=0; i<natoms; ++i) {
01659 ts->velocities[3*i+0] = data->ts.velocities[3*i+0];
01660 ts->velocities[3*i+1] = data->ts.velocities[3*i+1];
01661 ts->velocities[3*i+2] = data->ts.velocities[3*i+2];
01662 }
01663 }
01664 #endif
01665 ts->A = data->ts.A;
01666 ts->B = data->ts.B;
01667 ts->C = data->ts.C;
01668 ts->alpha = data->ts.alpha;
01669 ts->beta = data->ts.beta;
01670 ts->gamma = data->ts.gamma;
01671 }
01672 return MOLFILE_SUCCESS;
01673 }
01674
01675 static void close_hoomd_read(void *mydata) {
01676 hoomd_data_t *data = (hoomd_data_t *)mydata;
01677 int i;
01678
01679 if (data->p != NULL) XML_ParserFree(data->p);
01680
01681 fclose(data->fp);
01682
01683 SAFE_FREE(data->imagecounts);
01684 SAFE_FREE(data->from);
01685 SAFE_FREE(data->to);
01686 SAFE_FREE(data->bondtype);
01687 for (i=0; i < data->numbondtypes; ++i) {
01688 SAFE_FREE(data->bondtypename[i]);
01689 }
01690 SAFE_FREE(data->bondtypename);
01691 SAFE_FREE(data->angle);
01692 SAFE_FREE(data->angletype);
01693 for (i=0; i < data->numangletypes; ++i) {
01694 SAFE_FREE(data->angletypename[i]);
01695 }
01696 SAFE_FREE(data->angletypename);
01697 SAFE_FREE(data->dihedral);
01698 SAFE_FREE(data->dihedraltype);
01699 for (i=0; i < data->numdihedraltypes; ++i) {
01700 SAFE_FREE(data->dihedraltypename[i]);
01701 }
01702 SAFE_FREE(data->dihedraltypename);
01703 SAFE_FREE(data->improper);
01704 SAFE_FREE(data->impropertype);
01705 for (i=0; i < data->numimpropertypes; ++i) {
01706 SAFE_FREE(data->impropertypename[i]);
01707 }
01708 SAFE_FREE(data->impropertypename);
01709 SAFE_FREE(data->ts.coords);
01710 SAFE_FREE(data->ts.velocities);
01711 SAFE_FREE(data->filename);
01712 SAFE_FREE(data);
01713 }
01714
01718 static void *open_hoomd_write(const char *filename, const char *filetype,
01719 int natoms) {
01720 FILE *fd;
01721 hoomd_data_t *data;
01722 time_t mytime;
01723
01724 mytime = time(NULL);
01725
01726
01727 fd = fopen(filename, "w");
01728 if (!fd) {
01729 vmdcon_printf(VMDCON_ERROR, "hoomdplugin) Unable to open HOOMD-blue data file %s "
01730 "for writing\n", filename);
01731 return NULL;
01732 }
01733
01734 data = (hoomd_data_t *)malloc(sizeof(hoomd_data_t));
01735 data->numatoms = natoms;
01736 data->filename = strdup(filename);
01737 data->fp = fd;
01738 data->atomlist = NULL;
01739 data->numbonds = 0;
01740 data->to = NULL;
01741 data->from = NULL;
01742 data->bondtype = NULL;
01743 data->numbondtypes = 0;
01744 data->bondtypename = NULL;
01745 data->numangles = 0;
01746 data->angle = NULL;
01747 data->angletype = NULL;
01748 data->numangletypes = 0;
01749 data->angletypename = NULL;
01750 data->numdihedrals = 0;
01751 data->dihedral = NULL;
01752 data->dihedraltype = NULL;
01753 data->numdihedraltypes = 0;
01754 data->dihedraltypename = NULL;
01755 data->numimpropers = 0;
01756 data->improper = NULL;
01757 data->impropertype = NULL;
01758 data->numimpropertypes = 0;
01759 data->impropertypename = NULL;
01760 data->numframe = 0;
01761 data->imagecounts = (int *)malloc(sizeof(int)*3*(data->numatoms));
01762
01763 fputs("<?xml version =\"1.0\" encoding =\"UTF-8\" ?>\n",fd);
01764 fprintf(fd,"<%s version=\"%d.%d\">\n", hoomd_tag_name[HOOMD_XML],
01765 HOOMD_FORMAT_MAJV, HOOMD_FORMAT_MINV);
01766 fprintf(fd, "<!-- generated by VMD on: %s", ctime(&mytime));
01767 fprintf(fd, " %s plugin v%d.%d by %s -->\n", plugin.prettyname,
01768 plugin.majorv, plugin.minorv, plugin.author);
01769
01770 fprintf(data->fp,"<%s time_step=\"%d\" natoms=\"%d\">\n",
01771 hoomd_tag_name[HOOMD_CONFIG], data->numframe, data->numatoms);
01772
01773 return data;
01774 }
01775
01776
01777
01778 #define SECTION_OPEN(tag,num) \
01779 fprintf(data->fp, "<%s num=\"%d\">\n", hoomd_tag_name[tag], num)
01780
01781 #define SECTION_CLOSE(tag) \
01782 fprintf(data->fp, "</%s>\n", hoomd_tag_name[tag])
01783
01784 #define SECTION_WRITE_ATOMIC(tag,fmt,val) \
01785 SECTION_OPEN(tag,numatoms); \
01786 for (i=0; i < numatoms; ++i) \
01787 fprintf(data->fp,fmt, val); \
01788 SECTION_CLOSE(tag)
01789
01792 static int write_hoomd_structure(void *mydata, int optflags,
01793 const molfile_atom_t *atoms) {
01794 int i, numatoms;
01795 hoomd_data_t *data = (hoomd_data_t *)mydata;
01796 numatoms = data->numatoms;
01797
01798
01799
01800
01801 SECTION_WRITE_ATOMIC(HOOMD_TYPE,"%s\n",atoms[i].type);
01802 SECTION_WRITE_ATOMIC(HOOMD_BODY,"%d\n",atoms[i].resid);
01803
01804
01805 if (optflags & MOLFILE_RADIUS) {
01806 SECTION_WRITE_ATOMIC(HOOMD_DIAMETER,"%f\n",2.0f*atoms[i].radius);
01807 }
01808
01809 if (optflags & MOLFILE_MASS) {
01810 SECTION_WRITE_ATOMIC(HOOMD_MASS,"%f\n",atoms[i].mass);
01811 }
01812
01813 if (optflags & MOLFILE_CHARGE) {
01814 SECTION_WRITE_ATOMIC(HOOMD_CHARGE,"%f\n",atoms[i].charge);
01815 }
01816
01817 if ((data->numbonds > 0) && (data->from != NULL) && (data->to != NULL)) {
01818 SECTION_OPEN(HOOMD_BOND, data->numbonds);
01819 if (data->bondtype != NULL) {
01820 if (data->bondtypename != NULL) {
01821
01822 for (i=0; i < data->numbonds; ++i) {
01823 if (data->bondtype[i] < 0) {
01824
01825 fprintf(data->fp,"unkown %d %d\n", data->from[i]-1, data->to[i]-1);
01826 } else {
01827 fprintf(data->fp,"%s %d %d\n", data->bondtypename[data->bondtype[i]],
01828 data->from[i]-1, data->to[i]-1);
01829 }
01830 }
01831 } else {
01832
01833 for (i=0; i < data->numbonds; ++i) {
01834 fprintf(data->fp,"bondtype%d %d %d\n", data->bondtype[i],
01835 data->from[i]-1, data->to[i]-1);
01836 }
01837 }
01838 } else {
01839
01840 for (i=0; i < data->numbonds; ++i) {
01841 fprintf(data->fp,"bond %d %d\n", data->from[i]-1, data->to[i]-1);
01842 }
01843 }
01844 SECTION_CLOSE(HOOMD_BOND);
01845 }
01846
01847 if ( (data->numangles > 0) && (data->angle != NULL) ) {
01848 SECTION_OPEN(HOOMD_ANGLE,data->numangles);
01849 if (data->angletype != NULL) {
01850 if (data->angletypename != NULL) {
01851
01852 for (i=0; i < data->numangles; ++i) {
01853 if (data->angletype[i] < 0) {
01854
01855 fprintf(data->fp,"unkown %d %d %d\n", data->angle[3*i]-1,
01856 data->angle[3*i+1]-1, data->angle[3*i+2]-1);
01857 } else {
01858 fprintf(data->fp,"%s %d %d %d\n",
01859 data->angletypename[data->angletype[i]],
01860 data->angle[3*i]-1, data->angle[3*i+1]-1,
01861 data->angle[3*i+2]-1);
01862 }
01863 }
01864 } else {
01865
01866 for (i=0; i < data->numangles; ++i) {
01867 fprintf(data->fp,"angletype%d %d %d %d\n",
01868 data->angletype[i], data->angle[3*i]-1,
01869 data->angle[3*i+1]-1, data->angle[3*i+2]-1);
01870 }
01871 }
01872 } else {
01873
01874 for (i=0; i < data->numangles; ++i) {
01875 fprintf(data->fp,"angle %d %d %d\n", data->angle[3*i]-1,
01876 data->angle[3*i+1]-1, data->angle[3*i+2]-1);
01877 }
01878 }
01879 SECTION_CLOSE(HOOMD_ANGLE);
01880 }
01881
01882 if ( (data->numdihedrals > 0) && (data->dihedral != NULL) ) {
01883 SECTION_OPEN(HOOMD_DIHEDRAL, data->numdihedrals);
01884 if (data->dihedraltype != NULL) {
01885 if (data->dihedraltypename != NULL) {
01886
01887 for (i=0; i < data->numdihedrals; ++i) {
01888 if (data->dihedraltype[i] < 0) {
01889
01890 fprintf(data->fp,"unkown %d %d %d %d\n", data->dihedral[4*i]-1,
01891 data->dihedral[4*i+1]-1, data->dihedral[4*i+2]-1,
01892 data->dihedral[4*i+3]-1);
01893 } else {
01894 fprintf(data->fp,"%s %d %d %d %d\n",
01895 data->dihedraltypename[data->dihedraltype[i]],
01896 data->dihedral[4*i]-1, data->dihedral[4*i+1]-1,
01897 data->dihedral[4*i+2]-1, data->dihedral[4*i+3]-1);
01898 }
01899 }
01900 } else {
01901
01902 for (i=0; i < data->numdihedrals; ++i) {
01903 fprintf(data->fp,"dihedraltype%d %d %d %d %d\n",
01904 data->dihedraltype[i], data->dihedral[4*i]-1,
01905 data->dihedral[4*i+1]-1, data->dihedral[4*i+2]-1,
01906 data->dihedral[4*i+3]-1);
01907 }
01908 }
01909 } else {
01910
01911 for (i=0; i < data->numdihedrals; ++i) {
01912 fprintf(data->fp,"dihedral %d %d %d %d\n",
01913 data->dihedral[4*i]-1, data->dihedral[4*i+1]-1,
01914 data->dihedral[4*i+2]-1, data->dihedral[4*i+3]-1);
01915 }
01916 }
01917 SECTION_CLOSE(HOOMD_DIHEDRAL);
01918 }
01919
01920 if ( (data->numimpropers > 0) && (data->improper != NULL) ) {
01921 SECTION_OPEN(HOOMD_IMPROPER, data->numimpropers);
01922 if (data->impropertype != NULL) {
01923 if (data->impropertypename != NULL) {
01924
01925 for (i=0; i < data->numimpropers; ++i) {
01926 if (data->impropertype[i] < 0) {
01927
01928 fprintf(data->fp,"unkown %d %d %d %d\n", data->improper[4*i]-1,
01929 data->improper[4*i+1]-1, data->improper[4*i+2]-1,
01930 data->improper[4*i+3]-1);
01931 } else {
01932 fprintf(data->fp,"%s %d %d %d %d\n",
01933 data->impropertypename[data->impropertype[i]],
01934 data->improper[4*i]-1, data->improper[4*i+1]-1,
01935 data->improper[4*i+2]-1, data->improper[4*i+3]-1);
01936 }
01937 }
01938 } else {
01939
01940 for (i=0; i < data->numimpropers; ++i) {
01941 fprintf(data->fp,"impropertype%d %d %d %d %d\n",
01942 data->impropertype[i], data->improper[4*i]-1,
01943 data->improper[4*i+1]-1, data->improper[4*i+2]-1,
01944 data->improper[4*i+3]-1);
01945 }
01946 }
01947 } else {
01948
01949 for (i=0; i < data->numimpropers; ++i) {
01950 fprintf(data->fp,"improper %d %d %d %d\n",
01951 data->improper[4*i]-1, data->improper[4*i+1]-1,
01952 data->improper[4*i+2]-1, data->improper[4*i+3]-1);
01953 }
01954 }
01955 SECTION_CLOSE(HOOMD_IMPROPER);
01956 }
01957 return MOLFILE_SUCCESS;
01958 }
01959
01960 #undef SECTION_WRITE_ATOMIC
01961
01962
01963 #if vmdplugin_ABIVERSION > 14
01964 static int write_hoomd_bonds(void *v, int nbonds, int *fromptr, int *toptr,
01965 float *bondorder, int *bondtype,
01966 int nbondtypes, char **bondtypename) {
01967 hoomd_data_t *data = (hoomd_data_t *)v;
01968 int i;
01969
01970 data->numbonds=0;
01971 data->numbondtypes=0;
01972 data->bondtype=NULL;
01973
01974
01975 if ( (nbonds > 0) && (fromptr != NULL) && (toptr != NULL) ) {
01976 data->numbonds = nbonds;
01977 SAFE_CALLOC(data->from, int, nbonds);
01978 memcpy(data->from, fromptr, nbonds * sizeof(int));
01979 SAFE_CALLOC(data->to, int, nbonds);
01980 memcpy(data->to, toptr, nbonds * sizeof(int));
01981 if (bondtype != NULL) {
01982 SAFE_CALLOC(data->bondtype, int, nbonds);
01983 memcpy(data->bondtype, bondtype, nbonds * sizeof(int));
01984 }
01985 }
01986 if (nbondtypes > 0) {
01987 data->numbondtypes=nbondtypes;
01988 if (bondtypename != NULL) {
01989 SAFE_CALLOC(data->bondtypename, char *, nbondtypes+1);
01990 for (i=0; i < nbondtypes; ++i) {
01991 data->bondtypename[i] = strdup(bondtypename[i]);
01992 }
01993 }
01994 }
01995
01996 return MOLFILE_SUCCESS;
01997 }
01998 #else
01999 static int write_hoomd_bonds(void *v, int nbonds, int *fromptr, int *toptr, float *bondorder) {
02000 hoomd_data_t *data = (hoomd_data_t *)v;
02001
02002 data->numbonds=0;
02003 data->numbondtypes=0;
02004 data->bondtype=NULL;
02005
02006
02007 if ( (nbonds > 0) && (fromptr != NULL) && (toptr != NULL) ) {
02008 data->numbonds = nbonds;
02009 SAFE_CALLOC(data->from, int, nbonds);
02010 memcpy(data->from, fromptr, nbonds * sizeof(int));
02011 SAFE_CALLOC(data->to, int, nbonds);
02012 memcpy(data->to, toptr, nbonds * sizeof(int));
02013 }
02014
02015 return MOLFILE_SUCCESS;
02016 }
02017 #endif
02018
02019 #if vmdplugin_ABIVERSION > 15
02020 static int write_hoomd_angles(void * v, int numangles, const int *angles,
02021 const int *angletypes, int numangletypes,
02022 const char **angletypenames, int numdihedrals,
02023 const int *dihedrals, const int *dihedraltypes,
02024 int numdihedraltypes, const char **dihedraltypenames,
02025 int numimpropers, const int *impropers,
02026 const int *impropertypes, int numimpropertypes,
02027 const char **impropertypenames, int numcterms,
02028 const int *cterm, int ctermcols, int ctermrows) {
02029 hoomd_data_t *data = (hoomd_data_t *)v;
02030 int i;
02031
02032
02033 data->numangles = numangles;
02034 data->numdihedrals = numdihedrals;
02035 data->numimpropers = numimpropers;
02036
02037 if (data->numangles > 0) {
02038 SAFE_CALLOC(data->angle, int, 3*data->numangles);
02039 memcpy(data->angle, angles, 3*(data->numangles)*sizeof(int));
02040 }
02041 if (data->numdihedrals > 0) {
02042 SAFE_CALLOC(data->dihedral, int, 4*data->numdihedrals);
02043 memcpy(data->dihedral, dihedrals, 4*(data->numdihedrals)*sizeof(int));
02044 }
02045 if (data->numimpropers > 0) {
02046 SAFE_CALLOC(data->improper, int, 4*data->numimpropers);
02047 memcpy(data->improper, impropers, 4*(data->numimpropers)*sizeof(int));
02048 }
02049
02050 if (angletypes != NULL) {
02051 SAFE_CALLOC(data->angletype, int, data->numangles);
02052 memcpy(data->angletype, angletypes, (data->numangles)*sizeof(int));
02053 }
02054 if (dihedraltypes != NULL) {
02055 SAFE_CALLOC(data->dihedraltype, int, data->numdihedrals);
02056 memcpy(data->dihedraltype, dihedraltypes, (data->numdihedrals)*sizeof(int));
02057 }
02058 if (impropertypes != NULL) {
02059 SAFE_CALLOC(data->impropertype, int, data->numimpropers);
02060 memcpy(data->impropertype, impropertypes, (data->numimpropers)*sizeof(int));
02061 }
02062
02063 data->numangletypes = numangletypes;
02064 if (data->numangletypes > 0) {
02065 if (angletypenames != NULL) {
02066 SAFE_CALLOC(data->angletypename, char *, (data->numangletypes)+1);
02067 for (i=0; i < data->numangletypes; ++i) {
02068 data->angletypename[i] = strdup(angletypenames[i]);
02069 }
02070 }
02071 }
02072 data->numdihedraltypes = numdihedraltypes;
02073 if (data->numdihedraltypes > 0) {
02074 if (dihedraltypenames != NULL) {
02075 SAFE_CALLOC(data->dihedraltypename, char *, (data->numdihedraltypes)+1);
02076 for (i=0; i < data->numdihedraltypes; ++i) {
02077 data->dihedraltypename[i] = strdup(dihedraltypenames[i]);
02078 }
02079 }
02080 }
02081
02082 data->numimpropertypes = numimpropertypes;
02083 if (data->numimpropertypes > 0) {
02084 if (impropertypenames != NULL) {
02085 SAFE_CALLOC(data->impropertypename, char *, (data->numimpropertypes)+1);
02086 for (i=0; i < data->numimpropertypes; ++i) {
02087 data->impropertypename[i] = strdup(impropertypenames[i]);
02088 }
02089 }
02090 }
02091
02092 #if 0
02093 data->numcterms = numcterms;
02094 if (data->numcterms > 0) {
02095
02096 data->cterms = (int *) malloc(8*data->numcterms*sizeof(int));
02097 memcpy(data->cterms, cterms, 8*data->numcterms*sizeof(int));
02098 }
02099 #endif
02100 return MOLFILE_SUCCESS;
02101 }
02102
02103 #else
02104 static int write_hoomd_angles(void * v, int numangles, const int *angles,
02105 const double *angleforces, int numdihedrals,
02106 const int *dihedrals, const double *dihedralforces,
02107 int numimpropers, const int *impropers,
02108 const double *improperforces, int numcterms,
02109 const int *cterms, int ctermcols, int ctermrows,
02110 const double *ctermforces) {
02111 hoomd_data_t *data = (hoomd_data_t *)v;
02112
02113 data->numangles = numangles;
02114 if (data->numangles > 0) {
02115 data->angle = (int *) malloc(3*data->numangles*sizeof(int));
02116 memcpy(data->angle, angles, 3*data->numangles*sizeof(int));
02117 }
02118
02119 data->numdihedrals = numdihedrals;
02120 data->numimpropers = numimpropers;
02121 data->numcterms = numcterms;
02122
02123 if (data->numdihedrals > 0) {
02124 data->dihedrals = (int *) malloc(4*data->numdihedrals*sizeof(int));
02125 memcpy(data->dihedrals, dihedrals, 4*data->numdihedrals*sizeof(int));
02126 }
02127
02128 data->impropers = (int *) malloc(4*data->numimpropers*sizeof(int));
02129 memcpy(data->impropers, impropers, 4*data->numimpropers*sizeof(int));
02130 }
02131 #if 0
02132 if (data->numcterms > 0) {
02133
02134 data->cterms = (int *) malloc(8*data->numcterms*sizeof(int));
02135 memcpy(data->cterms, cterms, 8*data->numcterms*sizeof(int));
02136 }
02137 #endif
02138 return MOLFILE_SUCCESS;
02139 }
02140 #endif
02141
02142
02143 static int write_hoomd_timestep(void *mydata, const molfile_timestep_t *ts) {
02144 hoomd_data_t *data = (hoomd_data_t *)mydata;
02145 const float *pos;
02146 float px, py, pz;
02147 int i, *img, numatoms;
02148 numatoms = data->numatoms;
02149
02150 if (data->numframe != 0) {
02151 fprintf(data->fp,"<%s time_step=\"%d\" natoms=\"%d\" dimensions=\"3\">\n",
02152 hoomd_tag_name[HOOMD_CONFIG],
02153 #if vmdplugin_ABIVERSION > 10
02154 (int)(ts->physical_time+0.5f),
02155 #else
02156 data->numframe,
02157 #endif
02158 numatoms);
02159 }
02160
02161 fprintf(data->fp, "<%s lx=\"%f\" ly=\"%f\" lz=\"%f\" />\n",
02162 hoomd_tag_name[HOOMD_BOX], ts->A, ts->B, ts->C);
02163
02164 #define HOOMD_CONV_POS(crd, box, p, i) \
02165 if (fabsf(box) > SMALL) { \
02166 float tmp; \
02167 tmp = floorf(crd / box + 0.5f); \
02168 p = crd - tmp*box; \
02169 i = (int) tmp; \
02170 } else { \
02171 p = crd; \
02172 i = 0.0f; \
02173 }
02174
02175 SECTION_OPEN(HOOMD_POS, numatoms);
02176 pos = ts->coords;
02177 img = data->imagecounts;
02178 for (i = 0; i < numatoms; ++i) {
02179 HOOMD_CONV_POS(*pos, ts->A, px, *img);
02180 ++pos; ++img;
02181 HOOMD_CONV_POS(*pos, ts->B, py, *img);
02182 ++pos; ++img;
02183 HOOMD_CONV_POS(*pos, ts->C, pz, *img);
02184 ++pos; ++img;
02185
02186 fprintf(data->fp, "%.6f %.6f %.6f\n", px, py, pz);
02187 }
02188 SECTION_CLOSE(HOOMD_POS);
02189
02190
02191 if ( (fabsf(ts->A)+fabsf(ts->B)+fabsf(ts->C)) > SMALL) {
02192 SECTION_OPEN(HOOMD_IMAGE,numatoms);
02193 img = data->imagecounts;
02194 for (i = 0; i < numatoms; ++i) {
02195 fprintf(data->fp, "%d %d %d\n", img[0], img[1], img[2]);
02196 img += 3;
02197 }
02198 SECTION_CLOSE(HOOMD_IMAGE);
02199 }
02200
02201 #if vmdplugin_ABIVERSION > 10
02202 if (ts->velocities != NULL) {
02203 const float *vel;
02204
02205 SECTION_OPEN(HOOMD_VEL, numatoms);
02206 vel = ts->velocities;
02207 for (i = 0; i < numatoms; ++i) {
02208 fprintf(data->fp, "%.6f %.6f %.6f\n", vel[0], vel[1], vel[2]);
02209 vel += 3;
02210 }
02211 SECTION_CLOSE(HOOMD_VEL);
02212 }
02213 #endif
02214
02215 SECTION_CLOSE(HOOMD_CONFIG);
02216 data->numframe ++;
02217
02218 return MOLFILE_SUCCESS;
02219 }
02220
02221
02222 static void close_hoomd_write(void *mydata) {
02223 hoomd_data_t *data = (hoomd_data_t *)mydata;
02224 int i;
02225
02226
02227 if (data->numframe == 0) {
02228 SECTION_CLOSE(HOOMD_CONFIG);
02229 }
02230
02231
02232 SECTION_CLOSE(HOOMD_XML);
02233 fclose(data->fp);
02234
02235 if (data->numbonds > 0) {
02236 free(data->from);
02237 free(data->to);
02238 if (data->bondtype != NULL) {
02239 free(data->bondtype);
02240 }
02241 }
02242 if (data->numbondtypes > 0) {
02243 if (data->bondtypename != NULL) {
02244 for (i=0; i < data->numbondtypes; ++i) {
02245 free(data->bondtypename[i]);
02246 }
02247 free(data->bondtypename);
02248 }
02249 }
02250
02251 if (data->numangles > 0) {
02252 free(data->angle);
02253 if (data->angletype != NULL) {
02254 free(data->angletype);
02255 }
02256 }
02257 if (data->numangletypes > 0) {
02258 if (data->angletypename != NULL) {
02259 for (i=0; i < data->numangletypes; ++i) {
02260 free(data->angletypename[i]);
02261 }
02262 free(data->angletypename);
02263 }
02264 }
02265
02266 if (data->numdihedrals > 0) {
02267 free(data->dihedral);
02268 if (data->dihedraltype != NULL) {
02269 free(data->dihedraltype);
02270 }
02271 }
02272 if (data->numdihedraltypes > 0) {
02273 if (data->dihedraltypename != NULL) {
02274 for (i=0; i < data->numdihedraltypes; ++i) {
02275 free(data->dihedraltypename[i]);
02276 }
02277 free(data->dihedraltypename);
02278 }
02279 }
02280
02281 if (data->numimpropers > 0) {
02282 free(data->improper);
02283 if (data->impropertype != NULL) {
02284 free(data->impropertype);
02285 }
02286 }
02287 if (data->numimpropertypes > 0) {
02288 if (data->impropertypename != NULL) {
02289 for (i=0; i < data->numimpropertypes; ++i) {
02290 free(data->impropertypename[i]);
02291 }
02292 free(data->impropertypename);
02293 }
02294 }
02295
02296
02297 if (data->atomlist) {
02298 free(data->atomlist);
02299 }
02300 free(data->imagecounts);
02301 free(data->filename);
02302 free(data);
02303 }
02304
02305 #undef SECTION_OPEN
02306 #undef SECTION_CLOSE
02307
02308
02309
02310 VMDPLUGIN_API int VMDPLUGIN_init() {
02311 memset(&plugin, 0, sizeof(molfile_plugin_t));
02312 plugin.abiversion = vmdplugin_ABIVERSION;
02313 plugin.type = MOLFILE_PLUGIN_TYPE;
02314 plugin.name = "hoomd";
02315 plugin.prettyname = "HOOMD-blue XML File";
02316 plugin.author = "Axel Kohlmeyer";
02317 plugin.majorv = 0;
02318 plugin.minorv = 10;
02319 plugin.is_reentrant = VMDPLUGIN_THREADUNSAFE;
02320 plugin.filename_extension = "xml";
02321
02322 plugin.open_file_read = open_hoomd_read;
02323 plugin.read_structure = read_hoomd_structure;
02324 plugin.read_bonds = read_hoomd_bonds;
02325 plugin.read_next_timestep = read_hoomd_timestep;
02326 plugin.close_file_read = close_hoomd_read;
02327 #if vmdplugin_ABIVERSION > 10
02328 plugin.read_timestep_metadata = read_timestep_metadata;
02329 #endif
02330
02331 plugin.open_file_write = open_hoomd_write;
02332 plugin.write_bonds = write_hoomd_bonds;
02333 #if vmdplugin_ABIVERSION > 9
02334 plugin.read_angles = read_hoomd_angles;
02335 plugin.write_angles = write_hoomd_angles;
02336 #endif
02337 plugin.write_structure = write_hoomd_structure;
02338 plugin.write_timestep = write_hoomd_timestep;
02339 plugin.close_file_write = close_hoomd_write;
02340 return VMDPLUGIN_SUCCESS;
02341 }
02342
02343 VMDPLUGIN_API int VMDPLUGIN_register(void *v, vmdplugin_register_cb cb) {
02344 (*cb)(v, (vmdplugin_t *)&plugin);
02345 return VMDPLUGIN_SUCCESS;
02346 }
02347
02348 VMDPLUGIN_API int VMDPLUGIN_fini() {
02349 return VMDPLUGIN_SUCCESS;
02350 }
02351
02352
02353 #ifdef TEST_PLUGIN
02354
02355 int main(int argc, char *argv[]) {
02356 molfile_timestep_t timestep;
02357 molfile_atom_t *atoms;
02358 void *v, *w;
02359 int i, natoms, nbonds, optflags;
02360 #if vmdplugin_ABIVERSION > 10
02361 molfile_timestep_metadata_t ts_meta;
02362 #endif
02363 int nbondtypes, nangletypes, ndihtypes, nimptypes;
02364 int nangles, ndihedrals, nimpropers, ncterms;
02365 int ctermrows, ctermcols;
02366 int *from, *to, *bondtype, *angles, *angletype, *dihedrals, *dihedraltype;
02367 int *impropers, *impropertype, *cterms;
02368 #if vmdplugin_ABIVERSION < 16
02369 double *angleforces, *dihedralforces, *improperforces, *ctermforces;
02370 #endif
02371 float *order;
02372 char **btnames, **atnames, **dtnames, **itnames;
02373
02374 bondtype=NULL;
02375 btnames=NULL;
02376 nbondtypes=0;
02377 nangletypes=0;
02378
02379 VMDPLUGIN_init();
02380
02381 while (--argc) {
02382 ++argv;
02383 v = open_hoomd_read(*argv, "hoomd", &natoms);
02384 if (!v) {
02385 fprintf(stderr, "open_hoomd_read failed for file %s\n", *argv);
02386 return 1;
02387 }
02388 fprintf(stderr, "open_hoomd_read succeeded for file %s\n", *argv);
02389 fprintf(stderr, "number of atoms: %d\n", natoms);
02390
02391 atoms=(molfile_atom_t *)malloc(natoms*sizeof(molfile_atom_t));
02392 if (read_hoomd_structure(v, &optflags, atoms) != MOLFILE_SUCCESS) return 1;
02393 #if vmdplugin_ABIVERSION > 14
02394 read_hoomd_bonds(v, &nbonds, &from, &to, &order, &bondtype, &nbondtypes, &btnames);
02395
02396 read_hoomd_angles(v, &nangles, &angles, &angletype, &nangletypes, &atnames,
02397 &ndihedrals, &dihedrals, &dihedraltype, &ndihtypes, &dtnames,
02398 &nimpropers, &impropers, &impropertype, &nimptypes, &itnames,
02399 &ncterms, &cterms, &ctermcols, &ctermrows);
02400 #else
02401 read_hoomd_bonds(v, &nbonds, &from, &to, &order);
02402 #if vmdplugin_ABIVERSION > 9
02403 read_hoomd_angles(v, &nangles, &angles, &angleforces, &ndihedrals,
02404 &dihedrals, &dihedralforces, &nimpropers, &impropers,
02405 &improperforces, &ncterms, &cterms, &ctermcols,
02406 &ctermrows, &ctermforces);
02407 #endif
02408 #endif
02409 fprintf(stderr, "found: %d atoms, %d bonds, %d bondtypes, %d angles, "
02410 "%d angletypes\nfound: %d dihedrals, %d impropers, %d cterms\n",
02411 natoms, nbonds, nbondtypes, nangles, nangletypes, ndihedrals,
02412 nimpropers, ncterms);
02413
02414 fputs("ATOMS:\n", stderr);
02415 for(i=0; (i<20) && (i<natoms); ++i) {
02416 fprintf(stderr,"%05d: %s/%s %d\n",i+1,atoms[i].name,
02417 atoms[i].type, atoms[i].atomicnumber);
02418 }
02419
02420 fputs("BONDS:\n", stderr);
02421 if (nbonds > 0) {
02422 if (bondtype && nbondtypes > 0) {
02423 for(i=0; (i<20) && (i<nbonds);++i) {
02424 fprintf(stderr,"%05d: %s/%d %d %d\n", i+1, btnames[bondtype[i]],
02425 bondtype[i], from[i], to[i]);
02426 }
02427 } else {
02428 for(i=0;(i<20) && (i<nbonds);++i) {
02429 fprintf(stderr,"%05d: %d %d\n",i+1,from[i], to[i]);
02430 }
02431 }
02432 }
02433
02434 #if vmdplugin_ABIVERSION > 9
02435 fputs("ANGLES:\n", stderr);
02436 if (nangles > 0) {
02437 if (angletype && nangletypes > 0) {
02438 for(i=0; (i<20) && (i<nangles);++i) {
02439 fprintf(stderr,"%05d: %s/%d %d %d %d\n", i+1, atnames[angletype[i]],
02440 angletype[i], angles[3*i], angles[3*i+1], angles[3*i+2]);
02441 }
02442 } else {
02443 for(i=0;(i<20) && (i<nangles);++i) {
02444 fprintf(stderr,"%05d: %d %d %d\n",i+1, angles[3*i],
02445 angles[3*i+1], angles[3*i+2]);
02446 }
02447 }
02448 }
02449 #endif
02450
02451 i = 0;
02452 timestep.coords = (float *)malloc(3*sizeof(float)*natoms);
02453 #if vmdplugin_ABIVERSION > 10
02454
02455 timestep.velocities = NULL;
02456 read_timestep_metadata(v,&ts_meta);
02457 if (ts_meta.has_velocities) {
02458 fprintf(stderr, "found timestep velocities metadata.\n");
02459 }
02460 timestep.velocities = (float *) malloc(3*natoms*sizeof(float));
02461 #endif
02462
02463 while (!read_hoomd_timestep(v, natoms, ×tep)) {
02464 i++;
02465 }
02466 fprintf(stderr, "ended read_next_timestep on frame %d\n", i);
02467 if (argc < 2) {
02468 w = open_hoomd_write("test.xml","hoomd",natoms);
02469 #if vmdplugin_ABIVERSION > 14
02470 write_hoomd_bonds(w, nbonds, from, to, NULL, bondtype, nbondtypes, btnames);
02471 #else
02472 write_hoomd_bonds(w, nbonds, from, to, NULL);
02473 #endif
02474 #if vmdplugin_ABIVERSION > 9
02475 #if vmdplugin_ABIVERSION > 15
02476 write_hoomd_angles(w, nangles, angles, angletype, nangletypes, atnames,
02477 ndihedrals, dihedrals, dihedraltype, ndihtypes, dtnames,
02478 nimpropers, impropers, impropertype, nimptypes, itnames,
02479 ncterms, cterms, ctermcols, ctermrows);
02480 #else
02481 write_hoomd_angles(w, nangles, angles, NULL, ndihedrals, dihedrals, NULL,
02482 nimpropers, impropers, NULL, ncterms, cterms,
02483 ctermcols, ctermrows, NULL);
02484 #endif
02485 #endif
02486 write_hoomd_structure(w, optflags, atoms);
02487 write_hoomd_timestep(w, ×tep);
02488 close_hoomd_write(w);
02489 fprintf(stderr, "done with writing hoomd file test.xml.\n");
02490 }
02491 close_hoomd_read(v);
02492
02493 free(atoms);
02494 free(timestep.coords);
02495 #if vmd_ABIVERSION > 10
02496 free(timestep.velocities);
02497 #endif
02498 }
02499 VMDPLUGIN_fini();
02500 return 0;
02501 }
02502
02503 #endif
02504
