Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members   Related Pages  

hoomdplugin.c

Go to the documentation of this file.
00001 /***************************************************************************
00002  *
00003  * HOOMD and HOOMD-blue style XML format data/topology file reader and writer
00004  *
00005  * Copyright (c) 2009 Axel Kohlmeyer <akohlmey@gmail.com>
00006  *
00007  ***************************************************************************/
00008 
00009 /***************************************************************************
00010  * RCS INFORMATION:
00011  *
00012  *      $RCSfile: hoomdplugin.c,v $
00013  *      $Author: akohlmey $       $Locker:  $             $State: Exp $
00014  *      $Revision: 1.16 $       $Date: 2011/12/06 04:57:12 $
00015  *
00016  ***************************************************************************/
00017 
00018 /*
00019  * HOOMD-blue xml topology/trajectory file format:
00020  *
00021  * http://codeblue.umich.edu/hoomd-blue/doc/page_xml_file_format.html
00022  *
00023  * NOTES: 
00024  * - all XML elements are converted to lowercase before parsing, so case 
00025  *   will be ignored on reading. on writing only lowercase will be written.
00026  * - attributes are in general optional, but it is recommended to provide 
00027  *   them all. the plugin will print an informational message if it inserts 
00028  *   a missing attribute with a default value and a warning message if it 
00029  *   encounters a yet unknown attribute.
00030  * - units attributes are deprecated and always ignored.
00031  * - the angle, dihedral, improper nodes are only supported for 
00032  *   version 1.1 files. in files with no version attribute to hoomd_xml 
00033  *   or version 1.0 those nodes will be ignored.
00034  * - accelerations nodes appear in version 1.2. The molfile API
00035  *   has no method to pass them on, so they are currently ignored.
00036  * - charge nodes are only supported for version 1.3 and later files. 
00037  *   in files with no version attribute to hoomd_xml, version 1.0,
00038  *   or version 1.2, those nodes will be ignored.
00039  * - body nodes are supported starting with version 1.4. The body 
00040  *   value is mapped to the "resid" field in molfile.
00041  * - as of version 1.4 also the "vizsigma" attribute to the 
00042  *   configuration node is honored. It provides a scaling factor
00043  *   for the radius of particles and can be overridded by the 
00044  *   VMDHOOMDSIGMA environment variable. Both are ignored, however,
00045  *   if a diameter node is present.
00046  *
00047  ***************************************************
00048 
00049 <?xml version ="1.0" encoding ="UTF-8" ?>
00050 <hoomd_xml version="1.4">
00051   <!-- comments -->
00052   <configuration time_step="0" dimensions="3" vizsigma="0.5" ntypes="2">
00053     <box  lx="49.914" ly= "49.914" lz="49.914" />
00054     <position>
00055       3.943529 4.884501 12.317140
00056       3.985539 5.450660 13.107670
00057       3.521055 6.212838 13.516340
00058     </position>
00059     <type num="6">
00060       CT
00061       CM
00062       BB
00063       TRP
00064       LEU
00065       WS
00066     </type>
00067     <bond num="10">
00068       CT-CM 0 1
00069       CM-CT 1 2
00070       BB-TRP 3 4
00071       BB-BB  3 5
00072       BB-LEU 5 6
00073       BB-BB  5 7
00074       BB-TRP 7 8
00075       BB-BB  7 9
00076       BB-LEU 9 10
00077       BB-BB  9 3
00078     </bond>
00079     <angle num="9">
00080       CT-CM-CT   0  1  2
00081       BB-BB-TRP  9  3  4
00082       BB-BB-TRP  5  3  4
00083       BB-BB-LEU  3  5  6
00084       BB-BB-BB   3  5  7
00085       BB-BB-TRP  5  7  8
00086       BB-BB-BB   5  7  9
00087       BB-BB-LEU  7  9 10
00088       BB-BB-BB   7  9  3
00089     </angle>
00090     <diameter num="3">
00091       1.0 1.0 2.0
00092     </diameter>
00093     <mass num="3">
00094       1.0 1.0 2.0
00095     </mass>
00096     <charge num="3">
00097       -1.0 -1.0 2.0
00098     </charge>
00099     <body num="3">
00100       -1 0 0
00101     </charge>
00102     <wall>
00103       <coord ox="1.0" oy="2.0" oz="3.0" nx="4.0" ny="5.0" nz="6.0"/>
00104       <coord ox="7.0" oy="8.0" oz="9.0" nx="10.0" ny="11.0" nz="-12.0"/>
00105     </wall>
00106   </configuration>
00107 </hoomd>
00108 
00109 ****************************************************/
00110 
00111 #include "largefiles.h"   /* platform dependent 64-bit file I/O defines */
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 /* XXX: it is a bit annoying to have to resolve to an external XML 
00127  * library for parsing. XML is messy and expat is plain c and makes
00128  * life much easier. perhaps a portable internal xml parser would 
00129  * do as well, but the one used in hoomd is in c++ and requires to 
00130  * the the whole file into memory as one gian string. :-( */
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 /* numeric and symbolic representations of the supported xml nodes */
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 /* This list has to stay consistend with the defines above */
00168 static const char *hoomd_tag_name[HOOMD_NUMTAGS] = {
00169   "(none)", "hoomd_xml", "configuration", "box", /*  0 -  3 */
00170   "position", "image", "velocity",               /*  4 -  6 */
00171   "type", "body", "mass", "charge", "diameter",  /*  7 - 11 */
00172   "bond", "angle", "dihedral", "improper",       /* 12 - 15 */
00173   "wall", "coord", "acceleration", "orientation",/* 16 - 19 */
00174   "(error)"                                      /* 20 -    */
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 /* forward declaration */
00227 static molfile_plugin_t plugin;
00228 
00229 /*
00230  * Error reporting macro for use in DEBUG mode
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 /* make sure pointers are NULLed after free(3)ing them. */
00244 #define SAFE_FREE(ptr) free(ptr); ptr=NULL
00245 /* calloc with test of success */
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   /* convert text mode tags to numbers and parse attributes. */
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:             /* root node */
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:          /* new frame */
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:             /* new box dimensions */
00404       /* hoomd only supports orthogonal boxes at the moment */
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:             /* coordinates */
00423       for (i=0; attr[i]; i+=2) {
00424         XML_ATTR_COPY(attr, i, key, value);
00425         if (strcmp(key,"units") == 0) {
00426           ;                     /* ignore */
00427         } else if (strcmp(key,"num") == 0) {
00428           if (d->numatoms < 1) {
00429             d->numatoms = atoi(value); /* number of positions. */
00430           }
00431         } else {
00432           XML_UNK_ATTR(element, key, value);
00433         }
00434       }
00435       break;
00436 
00437     case HOOMD_VEL:             /* velocities */
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           ;                     /* ignore */
00444         } else if (strcmp(key,"num") == 0) {
00445           ;                     /* XXX: number of velocities. use for check. */
00446         } else {
00447           XML_UNK_ATTR(element, key, value);
00448         }
00449       }
00450 #endif
00451       break;
00452 
00453     case HOOMD_BODY:            /* body/resid tag */
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             ;                     /* XXX: number of body tags. use for check. */
00463           } else {
00464             XML_UNK_ATTR(element, key, value);
00465           }
00466         }
00467       }
00468       break;
00469 
00470     case HOOMD_MASS:            /* atom 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           ;                     /* ignore */
00476         } else if (strcmp(key,"num") == 0) {
00477           ;                     /* XXX: number of masses. use for check. */
00478         } else {
00479           XML_UNK_ATTR(element, key, value);
00480         }
00481       }
00482       break;
00483 
00484     case HOOMD_CHARGE:          /* atom 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             ;                     /* ignore */
00495           } else if (strcmp(key,"num") == 0) {
00496             ;                     /* XXX: number of masses. use for check. */
00497           } else {
00498             XML_UNK_ATTR(element, key, value);
00499           }
00500         }
00501       }
00502       break;
00503 
00504     case HOOMD_DIAMETER:          /* particle 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           ;                     /* ignore */
00510         } else if (strcmp(key,"num") == 0) {
00511           ;                     /* XXX: number of diameters. use for check. */
00512         } else {
00513           XML_UNK_ATTR(element, key, value);
00514         }
00515       }
00516       break;
00517 
00518     case HOOMD_IMAGE:           /* fallthrough */
00519     case HOOMD_TYPE:            /* fallthrough */
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           ;                 /* XXX: number of atom types. use for check. */
00525         } else {
00526           XML_UNK_ATTR(element, key, value);
00527         }
00528       }
00529       break;
00530 
00531     /* angle type definitions */
00532     case HOOMD_ANGLE:           /* fallthrough */
00533     case HOOMD_DIHEDRAL:        /* fallthrough */
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               ;                 /* XXX: number of angles. use for check. */
00547             } else {
00548               XML_UNK_ATTR(element, key, value);
00549             }
00550           }
00551         }
00552       }
00553       break;
00554 
00555     case HOOMD_WALL:            /* wall section */
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:           /* wall coordinate */
00563       for (i=0; attr[i]; i+=2) {
00564         XML_ATTR_COPY(attr, i, key, value);
00565         if (strcmp(key,"units") == 0) {
00566           ;                       /* ignore */
00567         } else if(strcmp(key,"ox") == 0) {
00568           ; /* do nothing. we don't do anything with walls, yet.*/
00569         } else if(strcmp(key,"oy") == 0) {
00570           ; /* do nothing. we don't do anything with walls, yet.*/
00571         } else if(strcmp(key,"oz") == 0) {
00572           ; /* do nothing. we don't do anything with walls, yet.*/
00573         } else if(strcmp(key,"nx") == 0) {
00574           ; /* do nothing. we don't do anything with walls, yet.*/
00575         } else if(strcmp(key,"ny") == 0) {
00576           ; /* do nothing. we don't do anything with walls, yet.*/
00577         } else if(strcmp(key,"nz") == 0) {
00578           ; /* do nothing. we don't do anything with walls, yet.*/
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               ;                 /* XXX: number of accelerations. use for check. */
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:            /* orientation tag */
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             ;                     /* XXX: number of orientation tags. use for check. */
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       /* store number of atoms the first time we read this section and
00645        * if it has not yet been passed on at the <configuration> tag.
00646        * VMD (currently) assumes that this doesn't change, so we
00647        * don't have to parse this more than once. */
00648       if (d->numatoms < 0) {    
00649         d->numatoms = d->counter;
00650       }
00651       break;
00652 
00653     case HOOMD_XML:             /* fallthrough. */
00654     case HOOMD_BOX:             /* fallthrough. */
00655     case HOOMD_IMAGE:           /* fallthrough. */
00656     case HOOMD_WALL:            /* fallthrough. */
00657     case HOOMD_COORD:           /* fallthrough. */
00658     case HOOMD_ACCEL:           /* fallthrough. */
00659     case HOOMD_ORIENT:          /* fallthrough. */
00660       ;                         /* nothing to do */
00661       break;
00662 
00663     case HOOMD_POS:
00664       /* assign the number of atoms if it has not been done before */
00665       if ( (d->counter > 0) && (d->numatoms < 0) ) {
00666         d->numatoms = d->counter/3;
00667       }                         /* fallthrough */
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:            /* fallthrough. */
00678     case HOOMD_CHARGE:          /* fallthrough. */
00679     case HOOMD_BODY:            /* fallthrough. */
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       /* store number of bonds the first time we read this section.
00691        * VMD (currently) assumes that this doesn't change, so we
00692        * don't have to do this more than once. */
00693       if (d->numbonds < 0) {
00694         d->numbonds = d->counter/3; /* three 'words' per bond */
00695       }
00696       break;
00697 
00698     case HOOMD_ANGLE:
00699       /* store number of angles the first time we read this section.
00700        * VMD (currently) assumes that this doesn't change, so we
00701        * don't have to do this more than once. */
00702       if (d->numangles < 0) {
00703         if (d->majv == 1 && d->minv > 0) { /* angle only valid from version 1.1 on */
00704           d->numangles = d->counter/4;    /* four 'words' per angle */
00705         } else {
00706           d->numangles = 0;
00707         }
00708       }
00709       break;
00710 
00711     case HOOMD_DIHEDRAL:
00712       /* store number of dihedrals the first time we read this section.
00713        * VMD (currently) assumes that this doesn't change, so we
00714        * don't have to do this more than once. */
00715       if (d->numdihedrals < 0) {
00716         if (d->majv == 1 && d->minv > 0) { /* dihedral only valid from version 1.1 on */
00717           d->numdihedrals = d->counter/5;    /* four 'words' per dihedral */
00718         } else {
00719           d->numdihedrals = 0;
00720         }
00721       }
00722       break;
00723 
00724     case HOOMD_IMPROPER:
00725       /* store number of impropers the first time we read this section.
00726        * VMD (currently) assumes that this doesn't change, so we
00727        * don't have to do this more than once. */
00728       if (d->numimpropers < 0) {
00729         if (d->majv == 1 && d->minv > 0) { /* improper only valid from version 1.1 on */
00730           d->numimpropers = d->counter/5;    /* four 'words' per improper */
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;          /* XXX: bug in program */
00764   if (d->parsedepth < 1) return;  /* not yet within XML block */
00765   if (d->parsedepth > HOOMD_MAXDEPTH) return;  /* too much nesting */
00766   if (len < 1) return;
00767 
00768   mytag = d->currtag[d->parsedepth];
00769   switch (mytag) {
00770     
00771     case HOOMD_TYPE:
00772 
00773       /* count the number of atoms by counting the types */
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         /* assign atom types the second time we parse the first frame */
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       /* count the number of bonds by counting the types */
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               /* CAVEAT: atom indices start at 1 here. psf style. XXX */
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) {  /* angle is only valid from version 1.1 on */
00858         /* count the number of angles by counting the types */
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                 /* CAVEAT: atom indices have to start at 1 here. psf style. XXX */
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) {  /* dihedral is only valid from version 1.1 on */
00901         /* count the number of dihedrals by counting the types */
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                 /* CAVEAT: atom indices start at 1 here. psf style. XXX */
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) {  /* improper is only valid from version 1.1 on */
00946         /* count the number of impropers by counting the types */
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                 /* CAVEAT: atom indices start at 1 here. psf style. XXX */
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       /* set radius, from diameter block. radius = 0.5*diameter. */
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       /* set mass. */
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       /* set resid from <body> tag. */
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       /* set charge. */
01069       if (d->majv == 1 && d->minv > 2) {  /* charge is only valid from version 1.2 on */
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       /* only try to read coordinates. if there is storage for them. */
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       /* scaling factor for guessed diameters */
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       /* loop through file until we have parsed the first configuration */
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       /* reset parsing */
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);  /* is grown on demand */
01354   SAFE_CALLOC(data->angletypename, char*, 1);  /* is grown on demand */
01355   SAFE_CALLOC(data->dihedraltypename, char*, 1);  /* is grown on demand */
01356   SAFE_CALLOC(data->impropertypename, char*, 1);  /* is grown on demand */
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   /* initialize atomdata with typical defaults */
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   /* read the first configuration again, but keep the data this time. */
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   /* allow overriding hardcoded sigma values from environment */
01416   envvar = getenv("VMDHOOMDSIGMA");
01417   if (envvar) data->mysigma = atof(envvar);
01418 
01419   /* fix up settings */
01420   for (i=0, a=atoms; i < data->numatoms; ++i, ++a) {
01421 
01422     if (!(data->optflags & MOLFILE_RADIUS)) {
01423       /* use vizsigma parameter to adjust guessed particle radii,
01424          but only for xml files that don't have a diameter section. */
01425       a->radius = data->mysigma * get_pte_vdw_radius(a->atomicnumber);
01426     }
01427 
01428     if (!(data->optflags & MOLFILE_MASS)) {
01429       /* guess atom mass, if not provided */
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;           /* hoomd does not support CMAP */
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;          /* we currently only support angles. */
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     /* read the next configuration. */
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     /* only save coords if we're given a timestep pointer, */
01646     /* otherwise assume that VMD wants us to skip past it. */
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     /* copy velocities */
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 /* some macros for consistency and convenience */
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   /* write fields we know about */
01800   /* required by molfile: */
01801   SECTION_WRITE_ATOMIC(HOOMD_TYPE,"%s\n",atoms[i].type);
01802   SECTION_WRITE_ATOMIC(HOOMD_BODY,"%d\n",atoms[i].resid);
01803 
01804   /* optional: */
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         /* case 1: symbolic bond types available */
01822         for (i=0; i < data->numbonds; ++i) {
01823           if (data->bondtype[i] < 0) {
01824             /* case 1a: symbolic bond types not assigned */
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         /* case 2: only numerical bond types available */
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       /* case 3: no bond type info available. */
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         /* case 1: symbolic angle types available */
01852         for (i=0; i < data->numangles; ++i) {
01853           if (data->angletype[i] < 0) {
01854             /* case 1a: symbolic angle types not assigned */
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         /* case 2: only numerical angle types available */
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       /* case 3: no angle type info available. */
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         /* case 1: symbolic dihedral types available */
01887         for (i=0; i < data->numdihedrals; ++i) {
01888           if (data->dihedraltype[i] < 0) {
01889             /* case 1a: symbolic dihedral types not assigned */
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         /* case 2: only numerical dihedral types available */
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       /* case 3: no dihedral type info available. */
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         /* case 1: symbolic improper types available */
01925         for (i=0; i < data->numimpropers; ++i) {
01926           if (data->impropertype[i] < 0) {
01927             /* case 1a: symbolic improper types not assigned */
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         /* case 2: only numerical improper types available */
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       /* case 3: no improper type info available. */
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 /* clean up, part 1 */
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   /* save info until we actually write out the structure data */
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   /* save info until we actually write out the structure data */
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   /* save info until we actually write out the structure file */
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                           /* the rest is not yet supported. */
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   /* save info until we actually write out the structure file */
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                           /* XXX not yet supported. */
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   /* only print image section if we have (some) box info. */
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   /* we have not yet written any positions, so the configuration tag is still open. */
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 /* cleanup, part 2 */
02305 #undef SECTION_OPEN
02306 #undef SECTION_CLOSE
02307 
02308 /* registration stuff */
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     /* XXX: this should be tested at each time step. */
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, &timestep)) {
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, &timestep);
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 

Generated on Fri Aug 14 03:05:21 2020 for VMD Plugins (current) by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002