vaspxmlplugin.c Source File

00001 /***************************************************************************
00002  * RCS INFORMATION:
00003  *
00004  *      $RCSfile: vaspxmlplugin.c,v $
00005  *      $Author: johns $       $Locker:  $             $State: Exp $
00006  *      $Revision: 1.9 $       $Date: 2016/11/06 17:43:13 $
00007  *
00008  ***************************************************************************/
00009 
00010 /*
00011  *  VASP plugins for VMD
00012  *  Sung Sakong, Dept. of Phys., Univsity Duisburg-Essen
00013  *  
00014  *  VASP manual   
00015  *  http://cms.mpi.univie.ac.at/vasp/
00016  * 
00017  *  LINUX
00018  *  g++ -O2 -Wall -I. -I$VMDBASEDIR/plugins/include -c vaspxmlplugin.c
00019  *  ld -shared -o vaspxmlplugin.so vaspxmlplugin.o
00020  *
00021  *  MACOSX
00022  *  c++ -O2 -Wall -I. -I$VMDBASEDIR/plugins/include -c vaspxmlplugin.c
00023  *  c++ -bundle -o vaspxmlplugin.so vaspxmlplugin.o
00024  *
00025  *  Install
00026  *  copy vaspxmlplugin.so $VMDBASEDIR/plugins/$ARCH/molfile
00027  */
00028 
00029  /* Be aware that the XML file has 2 extra frames in addition to the
00030   * ionic time steps. There is a 'initial position' frame, then the
00031   * frames of each time step, and in the end once again a 'final position'
00032   * frame. Hence, the xml file will load two frames more than the ones
00033   * in XDATCAR or OUTCAR; the first two and last two are thus identical.
00034   */
00035 
00036 #include <stdio.h>
00037 #include <stdlib.h>
00038 #include <string.h>
00039 
00040 #include "molfile_plugin.h"
00041 #include "vaspplugin.h"
00042 #include "periodic_table.h"
00043 
00044 
00045 static void *open_vaspxml_read(const char *filename, const char *filetype, int *natoms)
00046 {
00047   vasp_plugindata_t *data;
00048   char lineptr[LINESIZE];
00049   int cellcoords, finished;
00050 
00051   /* Verify that input is OK */
00052   if (!filename || !natoms) return NULL;
00053 
00054   /* Start with undefined value; set it after successful read */
00055   *natoms = MOLFILE_NUMATOMS_UNKNOWN;
00056 
00057   data = vasp_plugindata_malloc();
00058   if (!data) return NULL;
00059 
00060   data->file = fopen(filename, "rb");
00061   if (!data->file) {
00062     vasp_plugindata_free(data);
00063     return NULL;
00064   }
00065 
00066   data->filename = strdup(filename);
00067 
00068   /* Scan xml file */
00069   data->numatoms = cellcoords = finished = 0;
00070   while (fgets(lineptr, LINESIZE, data->file) && !finished) {
00071 
00072     if (strstr(lineptr, "SYSTEM") != NULL && data->titleline == NULL) {
00073        /* Extract title line */
00074        char *begin = strstr(lineptr, ">") + 1;
00075        char *end = strstr(lineptr, "</i>");
00076        if (end) *end = '\0';
00077        if (begin) data->titleline = strdup(begin);
00078 
00079     } else if (strstr(lineptr, "atominfo") != NULL && data->numatoms == 0) {
00080        /* Extract number of atoms */
00081        fgets(lineptr, LINESIZE, data->file);
00082        sscanf(lineptr, " <atoms> %d </atoms>", &data->numatoms);
00083 
00084     } else if (strstr(lineptr, "crystal") != NULL && cellcoords == 0) {
00085        /* Extract lattice vectors */
00086        int i;
00087        fgets(lineptr, LINESIZE, data->file);
00088        for (i = 0; i < 3 && fgets(lineptr, LINESIZE, data->file); ++i) cellcoords += sscanf(lineptr, " <v> %f %f %f </v>", &data->cell[i][0], &data->cell[i][1], &data->cell[i][2]);
00089     }
00090 
00091     finished = data->titleline != NULL && data->numatoms != 0 && cellcoords != 0;
00092   }
00093 
00094   if (data->numatoms <= 0) {
00095      vasp_plugindata_free(data);
00096      fprintf(stderr, "\n\nVASP xml read) ERROR: file '%s' does not contain the number of atoms.\n", filename);
00097      return NULL;
00098   }
00099 
00100   if (cellcoords != 9) {
00101      vasp_plugindata_free(data);
00102      fprintf(stderr, "\n\nVASP xml read) ERROR: file '%s' does not contain lattice vectors.\n", filename);
00103      return NULL;
00104   }
00105 
00106   vasp_buildrotmat(data);
00107 
00108   *natoms = data->numatoms;
00109   rewind(data->file);
00110 
00111   return data;
00112 }
00113 
00114 
00115 static int read_vaspxml_structure(void *mydata, int *optflags, molfile_atom_t *atoms)
00116 {
00117   vasp_plugindata_t *data = (vasp_plugindata_t *)mydata;
00118   int atomcount, coordscount, finished;
00119   char lineptr[LINESIZE];
00120 
00121   /* Verify that input is OK */
00122   if (!data || !optflags || !atoms) return MOLFILE_ERROR;
00123 
00124   *optflags = MOLFILE_MASS; /* we set atom mass from the PTE. */
00125   *optflags |= MOLFILE_ATOMICNUMBER | MOLFILE_RADIUS; 
00126 
00127   /* Scan xml file */
00128   atomcount = coordscount = finished = 0;
00129   while (fgets(lineptr, LINESIZE, data->file) && !finished) {
00130 
00131     /* Extract atom types */
00132     if (strstr(lineptr, "atomtype") != NULL && atomcount == 0) {
00133       int i;
00134       fgets(lineptr, LINESIZE, data->file);
00135       for (i = 0; i < data->numatoms; ++i, ++atomcount) {
00136         molfile_atom_t *atom = &(atoms[i]);
00137         char atomtype[5];
00138         int idx;
00139         fgets(lineptr, LINESIZE, data->file);
00140         if (1 != sscanf(lineptr, " <rc><c> %4s </c>", atomtype)) break;
00141 
00142         idx = get_pte_idx(atomtype);
00143 
00144         /* Required settings */
00145         strncpy(atom->name, get_pte_label(idx), sizeof(atom->name));
00146         strncpy(atom->type, atom->name, sizeof(atom->type));
00147         atom->resname[0] = '\0';
00148         atom->resid = 1;
00149         atom->segid[0]  ='\0';
00150         atom->chain[0] = '\0';
00151 
00152         /* Optional flags (as defined in *optflags) */
00153         atom->mass = get_pte_mass(idx);
00154         atom->radius = get_pte_vdw_radius(idx);
00155         atom->atomicnumber = idx;
00156       }
00157 
00158     /* Verify presence of coordinates for all atoms */
00159     } else if (strstr(lineptr, "positions") != NULL && coordscount == 0) {
00160         int i;
00161         for (i = 0; i < data->numatoms && fgets(lineptr, LINESIZE, data->file); ++i) {
00162            float x, y, z;
00163            if (3 != sscanf(lineptr, " <v> %f %f %f <\v>", &x, &y, &z)) break;
00164         }
00165         coordscount = 3 * i;
00166     }
00167 
00168     finished = atomcount != 0 && coordscount != 0;
00169   }
00170 
00171   if (atomcount != data->numatoms) {
00172     fprintf(stderr, "\n\nVASP xml read) ERROR: file '%s' does not have list of atom names.\n", data->filename);
00173     return MOLFILE_ERROR;
00174   }
00175 
00176   if (coordscount != 3 * data->numatoms) {
00177     fprintf(stderr, "\n\nVASP xml read)  file '%s' does not contain coordinates of all atoms.\n", data->filename);
00178     return MOLFILE_ERROR;
00179   }
00180 
00181   rewind(data->file);
00182 
00183   return MOLFILE_SUCCESS;
00184 }
00185 
00186 
00187 static int read_vaspxml_timestep(void *mydata, int natoms, molfile_timestep_t *ts)
00188 {
00189   vasp_plugindata_t *data = (vasp_plugindata_t *)mydata;
00190   char lineptr[LINESIZE];
00191 
00192   /* only save coords if we're given a timestep pointer, */
00193   /* otherwise assume that VMD wants us to skip past it. */
00194   if (!data || !ts) return MOLFILE_EOF;
00195 
00196   /* Scan xml file */
00197   while (fgets(lineptr, LINESIZE, data->file)) {
00198 
00199     /* Extract coordinates of all atoms */
00200     if (strstr(lineptr, "positions") != NULL) {
00201       int i;
00202       for (i = 0; i < data->numatoms && fgets(lineptr, LINESIZE, data->file); ++i) {
00203         float x, y, z, rotx, roty, rotz;
00204         if (3 != sscanf(lineptr, " <v> %f %f %f </v>", &x, &y, &z)) return MOLFILE_EOF;
00205 
00206         rotx = x*data->cell[0][0] + y*data->cell[1][0] + z*data->cell[2][0];
00207         roty = x*data->cell[0][1] + y*data->cell[1][1] + z*data->cell[2][1];
00208         rotz = x*data->cell[0][2] + y*data->cell[1][2] + z*data->cell[2][2];
00209 
00210         ts->coords[3*i  ] = data->rotmat[0][0]*rotx + data->rotmat[0][1]*roty + data->rotmat[0][2]*rotz;
00211         ts->coords[3*i+1] = data->rotmat[1][0]*rotx + data->rotmat[1][1]*roty + data->rotmat[1][2]*rotz;
00212         ts->coords[3*i+2] = data->rotmat[2][0]*rotx + data->rotmat[2][1]*roty + data->rotmat[2][2]*rotz;
00213       }
00214       vasp_timestep_unitcell(ts, data);
00215       return MOLFILE_SUCCESS;
00216     }
00217   }
00218 
00219   return MOLFILE_EOF;
00220 }
00221 
00222 
00223 static void close_vaspxml_read(void *mydata)
00224 {
00225   vasp_plugindata_t *data = (vasp_plugindata_t *)mydata;
00226   vasp_plugindata_free(data);
00227 }
00228 
00229 
00230 /* registration stuff */
00231 static molfile_plugin_t plugin;
00232 
00233 int VMDPLUGIN_init() {
00234   memset(&plugin, 0, sizeof(molfile_plugin_t));
00235   plugin.abiversion = vmdplugin_ABIVERSION;
00236   plugin.type = MOLFILE_PLUGIN_TYPE;
00237   plugin.name = "xml";
00238   plugin.prettyname = "VASP_xml";
00239   plugin.author = "Sung Sakong";
00240   plugin.majorv = 0;
00241   plugin.minorv = 7;
00242   plugin.is_reentrant = VMDPLUGIN_THREADUNSAFE;
00243   plugin.filename_extension = "xml";
00244   plugin.open_file_read = open_vaspxml_read;
00245   plugin.read_structure = read_vaspxml_structure;
00246   plugin.read_next_timestep = read_vaspxml_timestep;
00247   plugin.close_file_read = close_vaspxml_read;
00248   return VMDPLUGIN_SUCCESS;
00249 }
00250 
00251 int VMDPLUGIN_register(void *v, vmdplugin_register_cb cb) {
00252   (*cb)(v, (vmdplugin_t *)&plugin);
00253   return VMDPLUGIN_SUCCESS;
00254 }
00255 
00256 int VMDPLUGIN_fini() {
00257   return VMDPLUGIN_SUCCESS;
00258 }