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

jsplugin.c

Go to the documentation of this file.
00001 /***************************************************************************
00002  *cr                                                                       
00003  *cr            (C) Copyright 1995-2016 The Board of Trustees of the           
00004  *cr                        University of Illinois                       
00005  *cr                         All Rights Reserved                        
00006  *cr                                                                   
00007  ***************************************************************************/
00008 
00009 /***************************************************************************
00010  * RCS INFORMATION:
00011  *
00012  *      $RCSfile: jsplugin.c,v $
00013  *      $Author: johns $       $Locker:  $             $State: Exp $
00014  *      $Revision: 1.91 $       $Date: 2020/10/21 16:01:57 $
00015  *
00016  ***************************************************************************
00017  * DESCRIPTION:
00018  *   This plugin implements a high-performance binary molecular structure
00019  *   and trajectory storage format.  This file format currently uses a simple 
00020  *   non-redundant hash table approach for compression of per-atom
00021  *   character strings, properties, and tags.  Trajectory data is stored 
00022  *   with a file structure that avoids the need to transpose or convert
00023  *   dense blocks of cartesian coordinates into the most commonly used
00024  *   interleaved  x/y/z ordering.  The file structure also enables zero-copy 
00025  *   vectorized I/O methods to be used high-performance reads for 
00026  *   visualization and analysis with reduced operating system overhead.
00027  *   The plugin optionally supports the use of a block-based file structure
00028  *   and block-aligned memory buffers for direct I/O that bypasses the 
00029  *   OS filesystem buffer caches for multi-gigabyte-per-second read rates
00030  *   from SSD RAID arrays.
00031  *
00032  *   At present, only VMD, NAMD, and psfgen make use of this format.
00033  *   It started out as a test/example code and is slowly becoming
00034  *   more robust and featureful.
00035  *
00036  *   We should be able to implement a selective read approach that gathers
00037  *   discontiguous blocks of the trajectory using the POSIX lio_listio()
00038  *   APIs.  On Unix we currently use I/O wrappers that are based on the 
00039  *   lseek() and readv() APIs.  By using lio_listio() we could eliminate
00040  *   the separate lseek calls and service multiple timestep reads in a 
00041  *   single request, even included cases with discontiguous requests.
00042  *
00043  *   VMD test results for Linux host with an 8-way RAID0 of commodity 
00044  *   Intel 510 SSDs with SATA III 6Gbit/sec interfaces:
00045  *     Non-direct I/O using standard readv(): 1203 MB/sec
00046  *     Direct I/O, readv(), 4KB blocked file: 2130 MB/sec
00047  ***************************************************************************/
00048 /***
00049  ** Standalone test binary compilation flags for 64-bit Linux:
00050     cc -O3 -m64 -I../../include -DTEST_JSPLUGIN jsplugin.c \
00051       -o ~/bin/readjs -lm
00052  **
00053  ** Standalone test binary compilation flags for 64-bit Linux w/ CUDA:
00054     cc -O3 -m64 -I../../include -I/usr/local/cuda/include \
00055       -DTEST_JSPLUGIN -DENABLECUDATESTS jsplugin.c \
00056       -o ~/bin/readjs -L/usr/local/cuda/lib64 -lcudart -lm
00057  **
00058  ** Standalone test binary compilation flags for 64-bit Linux w/ CUDA+GDS:
00059     cc -O3 -m64 -I../../include -I/usr/local/cuda/include \
00060       -I/Projects/vmd/cuda/johns/gds/gds-alpha/lib \
00061       -DTEST_JSPLUGIN -DENABLECUDATESTS -DENABLECUDAGDS jsplugin.c \
00062       -o ~/bin/readjs -L/Projects/vmd/cuda/johns/gds/gds-alpha/lib -lcufile \
00063       -L/usr/local/cuda/lib64 -lcudart -lm
00064  **
00065  ** Standalone test binary compilation flags for Solaris:
00066     cc -fast -xarch=v9a -I../../include -DTEST_JSPLUGIN jsplugin.c \
00067       -o ~/bin/readjs -lm
00068  **
00069  ** Profiling flags for Solaris:
00070     cc -xpg -fast -xarch=v9a -g -I../../include -DTEST_JSPLUGIN jsplugin.c \
00071       -o ~/bin/readjs -lm
00072  **
00073  ** Test run for DGX-2:
00074     ~/bin/readjs /raid/scratch/solvatebar1204kb.js
00075  ***************************************************************************/
00076 
00077 #define INFOMSGS  1
00078 
00079 #if 1
00080 #define ENABLEJSSHORTREADS 1
00081 #endif
00082 
00083 #define VMDPLUGIN_STATIC
00084 #include "largefiles.h"   /* platform dependent 64-bit file I/O defines */
00085 #include "fastio.h"       /* must come before others, for O_DIRECT...   */
00086 
00087 #include <sys/stat.h>
00088 #include <sys/types.h>
00089 #include <stdio.h>
00090 #include <stdlib.h>
00091 #include <stddef.h>
00092 #include <string.h>
00093 #include <math.h>
00094 
00095 #include "hash.h"
00096 #include "endianswap.h"
00097 #include "molfile_plugin.h"
00098 
00099 
00100 /* allocate memory and return a pointer that is aligned on a given   */
00101 /* byte boundary, to be used for page- or sector-aligned I/O buffers */
00102 /* We use this if posix_memalign() is not available...               */
00103 #if defined(_WIN64)  /* sizeof(size_t) == sizeof(void*) */
00104 #define myintptrtype size_t
00105 #elif 1 /* sizeof(unsigned long) == sizeof(void*) */
00106 #define myintptrtype unsigned long
00107 #else
00108 #define myintptrtype uintptr_t  /* C99 */
00109 #endif
00110 /*
00111  * XXX On MSVC we get warnings about type conversions for 
00112  *     size_t vs. fio_size_t
00113  */
00114 static void *alloc_aligned_ptr(size_t sz, size_t blocksz, void **unalignedptr) {
00115   /* pad the allocation to an even multiple of the block size */
00116   size_t padsz = (sz + (blocksz - 1)) & (~(blocksz - 1));
00117   void * ptr = malloc(padsz + blocksz);
00118   *unalignedptr = ptr;
00119   return (void *) ((((myintptrtype) ptr) + (blocksz-1)) & (~(blocksz-1)));
00120 }
00121 
00122 
00123 #ifndef M_PI_2
00124 #define M_PI_2 1.57079632679489661922
00125 #endif
00126 
00127 #define JSHEADERSTRING   "JS Binary Structure and Trajectory File Format"                
00128 #define JSMAGICNUMBER    0x31337
00129 #define JSENDIANISM      0x12345678
00130 
00131 #define JSMAJORVERSION   2
00132 #define JSMINORVERSION   19
00133 
00134 #define JSNFRAMESOFFSET  (strlen(JSHEADERSTRING) + 20)
00135 
00136 #define JSNOERR             0
00137 #define JSBADFILE           1
00138 #define JSBADFORMAT         2
00139 
00140 
00141 /* Threshold atom count beyond which block-based I/O is used by default */
00142 /* The overhead from block-alignment padding bytes becomes essentially  */
00143 /* inconsequential (< 1%) for structures with more than 50,000 atoms.   */
00144 #define JSBLOCKIO_THRESH    50000
00145 
00146 
00147 /* Option flag macros and their meanings */
00148 #define JSOPT_NOOPTIONS     0x00000000  /* no structure, only coords    */
00149 
00150 /* Timesteps are block-size padded and page- or sector-aligned for      */
00151 /* direct I/O, using  OS-specific APIs that completely bypass the OS    */
00152 /* kernel filesystem buffer cache.                                      */
00153 /* The use of direct I/O APIs can raise performance tremendously on     */
00154 /* high-end RAIDs.  Tests on an 8-way RAID0 of Intel 510 SSDs raise the */
00155 /* peak I/O rate from 1100 MB/sec up to 2020 MB/sec with direct I/O.    */
00156 #define JSOPT_TS_BLOCKIO    0x10000000
00157 
00158 /* large data blocks */
00159 #define JSOPT_STRUCTURE     0x00000001  /* file contains structure data */
00160 #define JSOPT_BONDS         0x00000002  /* file contains bond info      */
00161 #define JSOPT_BONDORDERS    0x00000004  /* file contains bond orders    */
00162 #define JSOPT_ANGLES        0x00000008  /* file contains angle info     */
00163 #define JSOPT_CTERMS        0x00000010  /* file contains cross-terms    */
00164 
00165 /* optional per-atom fields */
00166 #define JSOPT_OCCUPANCY     0x00000100  /* file contains occupancy      */
00167 #define JSOPT_BFACTOR       0x00000200  /* file contains b-factor       */
00168 #define JSOPT_MASS          0x00000400  /* file contains masses         */
00169 #define JSOPT_CHARGE        0x00000800  /* file contains charges        */
00170 #define JSOPT_RADIUS        0x00001000  /* file contains radii          */
00171 #define JSOPT_ATOMICNUMBER  0x00002000  /* file contains atomic numbers */
00172 
00173 typedef struct {
00174   int verbose;                 /* flag to enable console info output    */
00175   fio_fd fd;                   /* main file descriptor                  */
00176   ptrdiff_t natoms;            /* handle uses a long type for natoms to */
00177                                /* help force promotion of file offset   */
00178                                /* arithmetic to long types              */
00179 
00180 #if JSMAJORVERSION > 1
00181   int parsed_structure;        /* flag indicating structure is parsed   */
00182   char *path;                  /* path to file                          */
00183 
00184   /* info for block-based direct I/O */ 
00185   int directio_pgsize_queried; /* caller has queried page/blocksz       */
00186   int directio_enabled;        /* block-based direct I/O is available   */
00187   fio_fd directio_fd;          /* block-based direct I/O using O_DIRECT */
00188   int directio_block_size;     /* block size to use for direct ts I/O   */
00189   void *directio_ucell_ptr;    /* unaligned unit cell buffer ptr        */
00190   void *directio_ucell_blkbuf; /* block-aligned unit cell buffer pt r   */
00191 
00192   /* timestep file offset, block padding, and stride information */
00193   fio_size_t ts_file_offset;   /* file offset to first timestep         */
00194   fio_size_t ts_crd_sz;        /* size of TS coordinates                */
00195   fio_size_t ts_crd_padsz;     /* size of TS block-padded coordinates   */
00196   fio_size_t ts_ucell_sz;      /* size of TS unit cell                  */
00197   fio_size_t ts_ucell_padsz;   /* size of TS block-padded unit cell     */
00198   
00199   /* structure info */
00200   int optflags;
00201   molfile_atom_t *atomlist;
00202   molfile_metadata_t *meta;
00203 
00204   /* bond info */
00205   int nbonds;
00206   int *bondfrom;
00207   int *bondto;
00208   float *bondorders;
00209 
00210   /* angle/dihedral/improper/cross-term info */
00211   int numangles, *angles;
00212   int numdihedrals, *dihedrals;
00213   int numimpropers, *impropers;
00214   int numcterms, *cterms;
00215 #endif
00216 
00217   /* trajectory info */
00218   int nframes;
00219   double tsdelta;
00220   int reverseendian;
00221   int with_unitcell;
00222 
00223   /* convenient buffer full of zeros for block-multiple padding */
00224   unsigned char blockpad[MOLFILE_DIRECTIO_MAX_BLOCK_SIZE];
00225 } jshandle;
00226 
00227 
00228 /* report the block size required to read this JS file */
00229 static int read_js_timestep_pagealign_size(void *v, int *pagealignsz) {
00230   jshandle *js = (jshandle *)v;
00231 
00232   // mark that the caller has queried the page alignment size
00233   js->directio_pgsize_queried = 1;
00234 
00235   /* assigne page alignment size based on file contents */
00236   if (js->optflags & JSOPT_TS_BLOCKIO) 
00237     *pagealignsz = js->directio_block_size;
00238   else 
00239     *pagealignsz = 1;
00240 
00241   return 0;
00242 }
00243 
00244 
00245 /* Use block-based I/O by default when writing structures larger */
00246 /* than JSBLOCKIO_THRESH atoms, or when directed by the user and */
00247 /* not otherwise prohibited...                                   */
00248 static void js_blockio_check_and_set(jshandle *js) {
00249   if ((getenv("VMDJSNOBLOCKIO") == NULL) && 
00250       ((js->natoms > JSBLOCKIO_THRESH) || getenv("VMDJSBLOCKIO"))) {
00251     js->optflags |= JSOPT_TS_BLOCKIO;
00252     js->directio_block_size = MOLFILE_DIRECTIO_MIN_BLOCK_SIZE; 
00253   }
00254 }
00255 
00256 
00257 static void *open_js_read(const char *path, const char *filetype, int *natoms) {
00258   jshandle *js;
00259   int jsmagicnumber, jsendianism, jsmajorversion, jsminorversion;
00260   struct stat stbuf;
00261   char strbuf[1024];
00262   int tmpnatoms=0;
00263 
00264   if (!path) return NULL;
00265 
00266   /* See if the file exists, and get its size */
00267   memset(&stbuf, 0, sizeof(struct stat));
00268   if (stat(path, &stbuf)) {
00269     printf("jsplugin) Could not access file '%s'.\n", path);
00270     perror("jsplugin) stat: ");
00271 /*    return NULL; */
00272   }
00273 
00274   js = (jshandle *)malloc(sizeof(jshandle));
00275   memset(js, 0, sizeof(jshandle));
00276   js->verbose = (getenv("VMDJSVERBOSE") != NULL);
00277 #if defined(_WIN64)
00278   js->verbose = 1;
00279 #endif
00280 
00281 #if JSMAJORVERSION > 1
00282   js->parsed_structure=0;
00283   js->directio_block_size=1;
00284   js->directio_ucell_ptr = NULL;
00285   js->directio_ucell_blkbuf = NULL;
00286 
00287   js->directio_pgsize_queried=0;
00288   js->directio_enabled=0;
00289   js->ts_file_offset=0;
00290   js->ts_crd_sz=0;
00291   js->ts_ucell_sz=0;
00292   js->ts_crd_padsz=0;
00293   js->ts_ucell_padsz=0;
00294 #endif
00295 
00296   if (fio_open(path, FIO_READ, &js->fd) < 0) {
00297     printf("jsplugin) Could not open file '%s' for reading.\n", path);
00298     free(js);
00299     return NULL;
00300   }
00301 
00302   /* emit header information */
00303   fio_fread(strbuf, strlen(JSHEADERSTRING), 1, js->fd);
00304   strbuf[strlen(JSHEADERSTRING)] = '\0';
00305   if (strcmp(strbuf, JSHEADERSTRING)) {
00306     printf("jsplugin) Bad trajectory header!\n");
00307     printf("jsplugin) Read string: %s\n", strbuf);
00308     fio_fclose(js->fd);
00309     free(js);
00310     return NULL;
00311   }
00312 
00313   fio_read_int32(js->fd, &jsmagicnumber);
00314   fio_read_int32(js->fd, &jsendianism);
00315   fio_read_int32(js->fd, &jsmajorversion);
00316   fio_read_int32(js->fd, &jsminorversion);
00317   fio_read_int32(js->fd, &tmpnatoms); /* handle-internal natoms is a long */
00318   fio_read_int32(js->fd, &js->nframes);
00319   if ((jsmagicnumber != JSMAGICNUMBER) || (jsendianism != JSENDIANISM)) {
00320 #if defined(INFOMSGS)
00321     if (js->verbose)
00322       printf("jsplugin) opposite endianism file, enabling byte swapping\n");
00323 #endif
00324     js->reverseendian = 1;
00325     swap4_aligned(&jsmagicnumber, 1);
00326     swap4_aligned(&jsendianism, 1);
00327     swap4_aligned(&jsmajorversion, 1);
00328     swap4_aligned(&jsminorversion, 1);
00329     swap4_aligned(&tmpnatoms, 1);
00330     swap4_aligned(&js->nframes, 1);
00331   } else {
00332 #if defined(INFOMSGS)
00333     if (js->verbose)
00334       printf("jsplugin) native endianism file\n");
00335 #endif
00336   }
00337 
00338   if ((jsmagicnumber != JSMAGICNUMBER) || (jsendianism != JSENDIANISM)) {
00339     fio_fclose(js->fd);
00340     free(js);
00341     return NULL;
00342   }
00343  
00344   if (jsmajorversion != JSMAJORVERSION) {
00345     printf("jsplugin) major version mismatch\n");
00346     printf("jsplugin)   file version: %d\n", jsmajorversion);
00347     printf("jsplugin)   plugin version: %d\n", JSMAJORVERSION);
00348     fio_fclose(js->fd);
00349     free(js);
00350     return NULL;
00351   }
00352 
00353   /* Copy integer natoms to handle natoms, could be a long. */
00354   /* The handle natoms uses long to help force promotion of */
00355   /* integer file offset calculations to long types...      */
00356   js->natoms = tmpnatoms;
00357   *natoms = tmpnatoms;
00358 
00359   /* copy path if we succeeded in opening the file */
00360   js->path = (char *) calloc(strlen(path)+1, 1);
00361   strcpy(js->path, path);
00362 
00363 #if 1
00364   /* read flags data from the file */
00365   fio_read_int32(js->fd, &js->optflags); 
00366   if (js->reverseendian)
00367     swap4_aligned(&js->optflags, 1);
00368 
00369 #if defined(INFOMSGS)
00370   if (js->verbose)
00371     printf("jsplugin) read option flags: %0x08x\n", js->optflags);
00372 #endif
00373 
00374   /* Check to see if block-based trajectory I/O is used  */
00375   /* and read in the block size for this file.           */
00376   if (js->optflags & JSOPT_TS_BLOCKIO) {
00377     fio_fread(&js->directio_block_size, sizeof(int), 1, js->fd);
00378     if (js->reverseendian)
00379       swap4_aligned(&js->directio_block_size, 1);
00380 
00381 #if defined(INFOMSGS)
00382     if (js->verbose) {
00383       printf("jsplugin) File uses direct I/O block size: %d bytes\n", 
00384              js->directio_block_size);
00385     }
00386 #endif
00387 
00388     /* Check to ensure that we can handle the block size used by the */
00389     /* file we are reading.  We may use variable block sizes in      */
00390     /* the future as more high-end filesystems begin to support      */
00391     /* 8KB, 16KB, or larger block sizes for enhanced sequential I/O  */
00392     /* performance on very large files.                              */
00393     if (js->directio_block_size > MOLFILE_DIRECTIO_MAX_BLOCK_SIZE) {
00394       printf("jsplugin) File block size exceeds jsplugin block size limit.\n");
00395       printf("jsplugin) Direct I/O unavailable for file '%s'\n", js->path);
00396     } else {
00397       if (fio_open(js->path, FIO_READ | FIO_DIRECT, &js->directio_fd) < 0) {
00398         printf("jsplugin) Direct I/O unavailable for file '%s'\n", js->path);
00399       } else {
00400         js->directio_enabled = 1;
00401       } 
00402     }
00403   }
00404 
00405 #if defined(ENABLEJSSHORTREADS)
00406   /* test code for an implementation that does short reads that */
00407   /* skip bulk solvent, useful for faster loading of very large */
00408   /* structures                                                 */
00409   if (getenv("VMDJSMAXATOMIDX") != NULL) {
00410     ptrdiff_t maxatomidx = atoi(getenv("VMDJSMAXATOMIDX"));
00411     if (maxatomidx < 0)
00412       maxatomidx = 0;
00413     if (maxatomidx >= js->natoms)
00414       maxatomidx = js->natoms - 1;
00415 
00416     printf("jsplugin) Short-reads of timesteps enabled: %td / %td atoms (%.2f%%)\n",
00417            maxatomidx, js->natoms, 100.0*(maxatomidx+1) / ((double) js->natoms));
00418   }
00419 #endif
00420 #endif
00421 
00422   return js;
00423 }
00424 
00425 
00426 #if JSMAJORVERSION > 1
00427 
00428 /* Compute the file offset for the first timestep and move */
00429 /* the file pointer to the correct position to read/write  */
00430 /* the first timestep.  Takes care of block alignment when */
00431 /* needed.                                                 */ 
00432 static int js_calc_timestep_blocking_info(void *mydata) {
00433   fio_size_t ts_block_offset, bszmask;
00434   jshandle *js = (jshandle *) mydata;
00435   int iorc=0;
00436 
00437   /* Record the current file offset so we can use it to */
00438   /* compute the absolute offset to the first timestep. */
00439   js->ts_file_offset = fio_ftell(js->fd);
00440 
00441   /* pad current offset to the start of the next block  */ 
00442   bszmask = js->directio_block_size - 1;
00443   ts_block_offset = (js->ts_file_offset + bszmask) & (~bszmask);
00444 
00445 #if defined(INFOMSGS)
00446   if (js->verbose) {
00447     printf("jsplugin) TS block size %td  curpos: %td  blockpos: %td\n", 
00448            (ptrdiff_t) js->directio_block_size, 
00449            (ptrdiff_t) js->ts_file_offset, 
00450            (ptrdiff_t) ts_block_offset);
00451   }
00452 #endif
00453 
00454   /* seek to the first block of the first timestep */
00455   js->ts_file_offset = ts_block_offset;
00456   if (js->directio_enabled)
00457     iorc = fio_fseek(js->directio_fd, js->ts_file_offset, FIO_SEEK_SET);
00458   else
00459     iorc = fio_fseek(js->fd, js->ts_file_offset, FIO_SEEK_SET);
00460   if (iorc < 0) {
00461     perror("jsplugin) fseek(): ");
00462   }
00463 
00464   /* compute timestep block padding/skipping for both */
00465   /* coordinate blocks and unit cell blocks           */
00466   js->ts_crd_sz = js->natoms * 3L * sizeof(float);
00467   js->ts_crd_padsz = (js->ts_crd_sz + bszmask) & (~bszmask);
00468 
00469   js->ts_ucell_sz = 6L * sizeof(double);
00470   js->ts_ucell_padsz = (js->ts_ucell_sz + bszmask) & (~bszmask);
00471 
00472   /* allocate TS unit cell buffer in an aligned, block-size-multiple buffer */
00473   /* unaligned unit cell buffer ptr */
00474 #if defined(USE_POSIX_MEMALIGN)
00475   if (posix_memalign((void**) &js->directio_ucell_ptr, 
00476       js->directio_block_size, js->ts_ucell_padsz)) {
00477     printf("jsplugin) Couldn't allocate aligned unit cell block buffer!\n");
00478   }
00479   /* the returned pointer is already block-aligned, and can free() */
00480   js->directio_ucell_blkbuf = js->directio_ucell_ptr;
00481 #else
00482   js->directio_ucell_blkbuf = (float *) 
00483     alloc_aligned_ptr(js->ts_ucell_padsz, js->directio_block_size, 
00484                       (void**) &js->directio_ucell_ptr);
00485 #endif
00486 
00487 #if defined(INFOMSGS)
00488   if (js->verbose) {
00489     printf("jsplugin) TS crds sz: %td psz: %td  ucell sz: %td psz: %td\n",
00490            (ptrdiff_t) js->ts_crd_sz,
00491            (ptrdiff_t) js->ts_crd_padsz, 
00492            (ptrdiff_t) js->ts_ucell_sz, 
00493            (ptrdiff_t) js->ts_ucell_padsz);
00494   }
00495 #endif
00496 
00497   return MOLFILE_SUCCESS;
00498 }
00499 
00500 
00501 static int read_js_structure(void *mydata, int *optflags,
00502                              molfile_atom_t *atoms) {
00503   jshandle *js = (jshandle *) mydata;
00504   ptrdiff_t i;
00505 
00506   if (optflags != NULL)
00507     *optflags = MOLFILE_NOOPTIONS; /* set to no options until we read them */
00508 
00509 #if 0
00510   /* read flags data from the file */
00511   fio_read_int32(js->fd, &js->optflags); 
00512   if (js->reverseendian)
00513     swap4_aligned(&js->optflags, 1);
00514 
00515 #if defined(INFOMSGS)
00516   if (js->verbose)
00517     printf("jsplugin) read option flags: %0x08x\n", js->optflags);
00518 #endif
00519 
00520   /* Check to see if block-based trajectory I/O is used  */
00521   /* and read in the block size for this file.           */
00522   if (js->optflags & JSOPT_TS_BLOCKIO) {
00523     fio_fread(&js->directio_block_size, sizeof(int), 1, js->fd);
00524     if (js->reverseendian)
00525       swap4_aligned(&js->directio_block_size, 1);
00526 
00527 #if defined(INFOMSGS)
00528     if (js->verbose) {
00529       printf("jsplugin) File uses direct I/O block size: %d bytes\n", 
00530              js->directio_block_size);
00531     }
00532 #endif
00533 
00534     /* Check to ensure that we can handle the block size used by the */
00535     /* file we are reading.  We may use variable block sizes in      */
00536     /* the future as more high-end filesystems begin to support      */
00537     /* 8KB, 16KB, or larger block sizes for enhanced sequential I/O  */
00538     /* performance on very large files.                              */
00539     if (js->directio_block_size > MOLFILE_DIRECTIO_MAX_BLOCK_SIZE) {
00540       printf("jsplugin) File block size exceeds jsplugin block size limit.\n");
00541       printf("jsplugin) Direct I/O unavailable for file '%s'\n", js->path);
00542     } else {
00543       if (fio_open(js->path, FIO_READ | FIO_DIRECT, &js->directio_fd) < 0) {
00544         printf("jsplugin) Direct I/O unavailable for file '%s'\n", js->path);
00545       } else {
00546         js->directio_enabled = 1;
00547       } 
00548     }
00549   }
00550 #endif
00551 
00552 
00553   /* emit warning message if the caller didn't check the required */
00554   /* alignment size, but the file supports block based direct I/O */
00555   if (js->directio_enabled && !js->directio_pgsize_queried) {
00556     printf("jsplugin) Warning: File supports block-based direct I/O, but\n");
00557     printf("jsplugin)          caller failed to query required alignment.\n");
00558     printf("jsplugin)          Block-based direct I/O is now disabled.\n");
00559      
00560     js->directio_enabled=0; // ensure we disable direct I/O early on
00561   }
00562 
00563 #if defined(INFOMSGS)
00564   if (js->verbose) {
00565     printf("jsplugin) Direct I/O %sabled for file '%s'\n", 
00566            (js->directio_enabled) ? "en" : "dis", js->path);
00567   }
00568 #endif
00569 
00570 
00571 #if 0
00572 #if defined(ENABLEJSSHORTREADS)
00573   /* test code for an implementation that does short reads that */
00574   /* skip bulk solvent, useful for faster loading of very large */
00575   /* structures                                                 */
00576   if (getenv("VMDJSMAXATOMIDX") != NULL) {
00577     ptrdiff_t maxatomidx = atoi(getenv("VMDJSMAXATOMIDX"));
00578     if (maxatomidx < 0)
00579       maxatomidx = 0;
00580     if (maxatomidx >= js->natoms)
00581       maxatomidx = js->natoms - 1;
00582 
00583     printf("jsplugin) Short-reads of timesteps enabled: %ld / %ld atoms (%.2f%%)\n",
00584            maxatomidx, js->natoms, 100.0*(maxatomidx+1) / ((double) js->natoms));
00585   }
00586 #endif
00587 #endif
00588 
00589   /* Mark the handle to indicate we've parsed the structure.             */
00590   /* If any errors occur after this point, they are likely fatal anyway. */
00591   js->parsed_structure = 1;
00592 
00593   /* determine whether or not this file contains structure info or not */
00594   if (js->optflags & JSOPT_STRUCTURE) {
00595     int numatomnames, numatomtypes, numresnames, numsegids, numchains;
00596     char **atomnames = NULL;
00597     char **atomtypes = NULL;
00598     char **resnames = NULL;
00599     char **segids = NULL;
00600     char **chains = NULL;
00601     short *shortbuf = NULL; /* temp buf for decoding atom records */
00602     int *intbuf = NULL;     /* temp buf for decoding atom records */
00603     float *fltbuf = NULL;   /* temp buf for decoding atom records */
00604  
00605     /* read in block of name string table sizes */
00606     fio_read_int32(js->fd, &numatomnames); 
00607     fio_read_int32(js->fd, &numatomtypes); 
00608     fio_read_int32(js->fd, &numresnames);
00609     fio_read_int32(js->fd, &numsegids);
00610     fio_read_int32(js->fd, &numchains); 
00611     if (js->reverseendian) {
00612       swap4_aligned(&numatomnames, 1);
00613       swap4_aligned(&numatomtypes, 1);
00614       swap4_aligned(&numresnames, 1);
00615       swap4_aligned(&numsegids, 1);
00616       swap4_aligned(&numchains, 1);
00617     }
00618 
00619 #if defined(INFOMSGS)
00620     if (js->verbose) {
00621       printf("jsplugin) reading string tables...\n");
00622       printf("jsplugin) %d %d %d %d %d\n",
00623              numatomnames, numatomtypes, numresnames, numsegids, numchains);
00624     }
00625 #endif
00626 
00627     /* skip forward to first TS if the caller gives us NULL ptrs */
00628     if (optflags == NULL && atoms == NULL) {
00629       size_t offset=0;
00630       offset += numatomnames * (16L * sizeof(char));
00631       offset += numatomtypes * (16L * sizeof(char));
00632       offset += numresnames  * (8L * sizeof(char));
00633       offset += numsegids    * (8L * sizeof(char));
00634       offset += numchains    * (2L * sizeof(char));
00635       offset += js->natoms * sizeof(short); /* atom name indices    */
00636       offset += js->natoms * sizeof(short); /* atom type indices    */
00637       offset += js->natoms * sizeof(short); /* residue name indices */
00638       offset += js->natoms * sizeof(short); /* segment name indices */
00639       offset += js->natoms * sizeof(short); /* chain name indices   */
00640       offset += js->natoms * sizeof(int);   /* residue indices      */
00641       
00642       /* optional per-atom fields */
00643       if (js->optflags & JSOPT_OCCUPANCY)
00644         offset += js->natoms * sizeof(float); 
00645       if (js->optflags & JSOPT_BFACTOR)
00646         offset += js->natoms * sizeof(float); 
00647       if (js->optflags & JSOPT_MASS)
00648         offset += js->natoms * sizeof(float); 
00649       if (js->optflags & JSOPT_CHARGE)
00650         offset += js->natoms * sizeof(float); 
00651       if (js->optflags & JSOPT_RADIUS)
00652         offset += js->natoms * sizeof(float); 
00653       if (js->optflags & JSOPT_ATOMICNUMBER)
00654         offset += js->natoms * sizeof(int);
00655 
00656       fio_fseek(js->fd, offset, FIO_SEEK_CUR);
00657       offset=0;
00658 
00659       /* these require actually seeking as we process... */
00660       if (js->optflags & JSOPT_BONDS) {
00661         fio_fread(&js->nbonds, sizeof(int), 1, js->fd);
00662         if (js->reverseendian)
00663           swap4_aligned(&js->nbonds, 1);
00664 #if defined(INFOMSGS)
00665         if (js->verbose) {
00666           printf("jsplugin)   %d bonds...\n", js->nbonds);
00667         }
00668 #endif
00669 
00670         offset += 2L * js->nbonds * sizeof(int);
00671         if (js->optflags & JSOPT_BONDORDERS)
00672           offset += js->nbonds * sizeof(float);
00673 
00674         fio_fseek(js->fd, offset, FIO_SEEK_CUR);
00675         offset=0;
00676       }
00677 
00678       if (js->optflags & JSOPT_ANGLES) {
00679         fio_fread(&js->numangles, sizeof(int), 1, js->fd);
00680         if (js->reverseendian)
00681           swap4_aligned(&js->numangles, 1);
00682 #if defined(INFOMSGS)
00683         if (js->verbose) {
00684           printf("jsplugin)   %d angles...\n", js->numangles);
00685         }
00686 #endif
00687         fio_fseek(js->fd, sizeof(int)*3L*js->numangles, FIO_SEEK_CUR);
00688 
00689         fio_fread(&js->numdihedrals, sizeof(int), 1, js->fd);
00690         if (js->reverseendian)
00691           swap4_aligned(&js->numdihedrals, 1);
00692 #if defined(INFOMSGS)
00693         if (js->verbose) {
00694           printf("jsplugin)   %d dihedrals...\n", js->numdihedrals);
00695         }
00696 #endif
00697         fio_fseek(js->fd, sizeof(int)*4L*js->numdihedrals, FIO_SEEK_CUR);
00698 
00699         fio_fread(&js->numimpropers, sizeof(int), 1, js->fd);
00700         if (js->reverseendian)
00701           swap4_aligned(&js->numimpropers, 1);
00702 #if defined(INFOMSGS)
00703         if (js->verbose) {
00704           printf("jsplugin)   %d impropers...\n", js->numimpropers);
00705         }
00706 #endif
00707         fio_fseek(js->fd, sizeof(int)*4L*js->numimpropers, FIO_SEEK_CUR);
00708       }
00709 
00710       if (js->optflags & JSOPT_CTERMS) {
00711         fio_fread(&js->numcterms, sizeof(int), 1, js->fd);
00712         if (js->reverseendian)
00713           swap4_aligned(&js->numcterms, 1);
00714 #if defined(INFOMSGS)
00715         if (js->verbose) {
00716           printf("jsplugin)   %d cterms...\n", js->numcterms);
00717         }
00718 #endif
00719         fio_fseek(js->fd, sizeof(int)*8L*js->numcterms, FIO_SEEK_CUR);
00720       }
00721   
00722       /* record the file offset for the first timestep */
00723       js_calc_timestep_blocking_info(js);
00724 
00725       return MOLFILE_SUCCESS;
00726     }
00727 
00728 
00729     /* allocate string tables */
00730     atomnames = (char **) malloc(numatomnames * sizeof(char *));
00731     atomtypes = (char **) malloc(numatomtypes * sizeof(char *));
00732     resnames  = (char **) malloc(numresnames  * sizeof(char *));
00733     segids    = (char **) malloc(numsegids    * sizeof(char *));
00734     chains    = (char **) malloc(numchains    * sizeof(char *));
00735 
00736 #if defined(INFOMSGS)
00737     if (js->verbose)
00738       printf("jsplugin)   atom names...\n");
00739 #endif
00740 
00741     /* read in the string tables */
00742     for (i=0; i<numatomnames; i++) {
00743       atomnames[i] = (char *) malloc(16L * sizeof(char));
00744       fio_fread(atomnames[i], 16L * sizeof(char), 1, js->fd);
00745     }
00746 
00747 #if defined(INFOMSGS)
00748     if (js->verbose)
00749       printf("jsplugin)   atom types...\n");
00750 #endif
00751     for (i=0; i<numatomtypes; i++) {
00752       atomtypes[i] = (char *) malloc(16L * sizeof(char));
00753       fio_fread(atomtypes[i], 16L * sizeof(char), 1, js->fd);
00754     }
00755 
00756 #if defined(INFOMSGS)
00757     if (js->verbose)
00758       printf("jsplugin)   residue names...\n");
00759 #endif
00760     for (i=0; i<numresnames; i++) {
00761       resnames[i] = (char *) malloc(8L * sizeof(char));
00762       fio_fread(resnames[i], 8L * sizeof(char), 1, js->fd);
00763     }
00764 
00765 #if defined(INFOMSGS)
00766     if (js->verbose)
00767       printf("jsplugin)   segment names...\n");
00768 #endif
00769     for (i=0; i<numsegids; i++) {
00770       segids[i] = (char *) malloc(8L * sizeof(char));
00771       fio_fread(segids[i], 8L * sizeof(char), 1, js->fd);
00772     }
00773 
00774 #if defined(INFOMSGS)
00775     if (js->verbose)
00776       printf("jsplugin)   chain names...\n");
00777 #endif
00778     for (i=0; i<numchains; i++) {
00779       chains[i] = (char *) malloc(2L * sizeof(char));
00780       fio_fread(chains[i], 2L * sizeof(char), 1, js->fd);
00781     }
00782 
00783 #if defined(INFOMSGS)
00784     if (js->verbose)
00785       printf("jsplugin) reading numeric field tables...\n");
00786 #endif
00787     /* read in all of the atom fields */
00788     shortbuf = (short *) malloc(js->natoms * sizeof(short));
00789 
00790 #if defined(INFOMSGS)
00791     if (js->verbose)
00792       printf("jsplugin)   atom name indices...\n");
00793 #endif
00794     /* read in atom names */
00795     fio_fread(shortbuf, js->natoms * sizeof(short), 1, js->fd);
00796     if (js->reverseendian)
00797       swap2_aligned(shortbuf, js->natoms);
00798     for (i=0; i<js->natoms; i++) {
00799       strcpy(atoms[i].name, atomnames[shortbuf[i]]);
00800     }
00801     for (i=0; i<numatomnames; i++)
00802       free(atomnames[i]);
00803     free(atomnames);
00804 
00805 #if defined(INFOMSGS)
00806     if (js->verbose)
00807       printf("jsplugin)   atom type indices...\n");
00808 #endif
00809     /* read in atom types */
00810     fio_fread(shortbuf, js->natoms * sizeof(short), 1, js->fd);
00811     if (js->reverseendian)
00812       swap2_aligned(shortbuf, js->natoms);
00813     for (i=0; i<js->natoms; i++) {
00814       strcpy(atoms[i].type, atomtypes[shortbuf[i]]);
00815     }
00816     for (i=0; i<numatomtypes; i++)
00817       free(atomtypes[i]);
00818     free(atomtypes);
00819 
00820 #if defined(INFOMSGS)
00821     if (js->verbose)
00822       printf("jsplugin)   residue name indices...\n");
00823 #endif
00824     /* read in resnames */
00825     fio_fread(shortbuf, js->natoms * sizeof(short), 1, js->fd);
00826     if (js->reverseendian)
00827       swap2_aligned(shortbuf, js->natoms);
00828     for (i=0; i<js->natoms; i++) {
00829       strcpy(atoms[i].resname, resnames[shortbuf[i]]);
00830     }
00831     for (i=0; i<numresnames; i++)
00832       free(resnames[i]);
00833     free(resnames);
00834     
00835 #if defined(INFOMSGS)
00836     if (js->verbose)
00837       printf("jsplugin)   segment name indices...\n");
00838 #endif
00839     /* read in segids */
00840     fio_fread(shortbuf, js->natoms * sizeof(short), 1, js->fd);
00841     if (js->reverseendian)
00842       swap2_aligned(shortbuf, js->natoms);
00843     for (i=0; i<js->natoms; i++) {
00844       strcpy(atoms[i].segid, segids[shortbuf[i]]);
00845     }
00846     for (i=0; i<numsegids; i++)
00847       free(segids[i]);
00848     free(segids);
00849 
00850 #if defined(INFOMSGS)
00851     if (js->verbose)
00852       printf("jsplugin)   chain name indices...\n");
00853 #endif
00854     /* read in chains */
00855     fio_fread(shortbuf, js->natoms * sizeof(short), 1, js->fd);
00856     if (js->reverseendian)
00857       swap2_aligned(shortbuf, js->natoms);
00858     for (i=0; i<js->natoms; i++) {
00859       strcpy(atoms[i].chain, chains[shortbuf[i]]);
00860     }
00861     for (i=0; i<numchains; i++)
00862       free(chains[i]);
00863     free(chains);
00864 
00865     if (shortbuf != NULL) {
00866       free(shortbuf);
00867       shortbuf=NULL;
00868     }
00869 
00870     /* 
00871      * read in integer data blocks 
00872      */
00873     intbuf = (int *) malloc(js->natoms * sizeof(int));
00874 
00875 #if defined(INFOMSGS)
00876     if (js->verbose)
00877       printf("jsplugin)   residue indices...\n");
00878 #endif
00879     /* read in resid */
00880     fio_fread(intbuf, js->natoms * sizeof(int), 1, js->fd);
00881     if (js->reverseendian)
00882       swap4_aligned(intbuf, js->natoms);
00883     for (i=0; i<js->natoms; i++) {
00884       atoms[i].resid = intbuf[i];
00885     }    
00886      
00887     if (intbuf != NULL) {
00888       free(intbuf);
00889       intbuf = NULL;
00890     }
00891 
00892 
00893 #if defined(INFOMSGS)
00894     if (js->verbose)
00895       printf("jsplugin) reading optional per-atom tables...\n");
00896 #endif
00897     /*
00898      * read in optional single-precision float data blocks
00899      */ 
00900     if (js->optflags & (JSOPT_OCCUPANCY | JSOPT_BFACTOR | 
00901         JSOPT_MASS | JSOPT_RADIUS | JSOPT_CHARGE)) 
00902       fltbuf = (float *) malloc(js->natoms * sizeof(float));
00903 
00904     /* read in optional data if it exists */
00905     if (js->optflags & JSOPT_OCCUPANCY) {
00906 #if defined(INFOMSGS)
00907       if (js->verbose)
00908         printf("jsplugin)   occupancy...\n");
00909 #endif
00910       *optflags |= MOLFILE_OCCUPANCY;
00911       fio_fread(fltbuf, js->natoms * sizeof(float), 1, js->fd);
00912       if (js->reverseendian)
00913         swap4_aligned(fltbuf, js->natoms);
00914       for (i=0; i<js->natoms; i++) {
00915         atoms[i].occupancy = fltbuf[i];
00916       }    
00917     }
00918 
00919     if (js->optflags & JSOPT_BFACTOR) {
00920 #if defined(INFOMSGS)
00921       if (js->verbose)
00922         printf("jsplugin)   bfactor...\n");
00923 #endif
00924       *optflags |= MOLFILE_BFACTOR;
00925       fio_fread(fltbuf, js->natoms * sizeof(float), 1, js->fd);
00926       if (js->reverseendian)
00927         swap4_aligned(fltbuf, js->natoms);
00928       for (i=0; i<js->natoms; i++) {
00929         atoms[i].bfactor = fltbuf[i];
00930       }    
00931     }
00932 
00933     if (js->optflags & JSOPT_MASS) { 
00934 #if defined(INFOMSGS)
00935       if (js->verbose)
00936         printf("jsplugin)   mass...\n");
00937 #endif
00938       *optflags |= MOLFILE_MASS;
00939       fio_fread(fltbuf, js->natoms * sizeof(float), 1, js->fd);
00940       if (js->reverseendian)
00941         swap4_aligned(fltbuf, js->natoms);
00942       for (i=0; i<js->natoms; i++) {
00943         atoms[i].mass = fltbuf[i];
00944       }    
00945     }
00946 
00947     if (js->optflags & JSOPT_CHARGE) { 
00948 #if defined(INFOMSGS)
00949       if (js->verbose)
00950         printf("jsplugin)   charge...\n");
00951 #endif
00952       *optflags |= MOLFILE_CHARGE;
00953       fio_fread(fltbuf, js->natoms * sizeof(float), 1, js->fd);
00954       if (js->reverseendian)
00955         swap4_aligned(fltbuf, js->natoms);
00956       for (i=0; i<js->natoms; i++) {
00957         atoms[i].charge = fltbuf[i];
00958       }    
00959     }
00960 
00961     if (js->optflags & JSOPT_RADIUS) { 
00962 #if defined(INFOMSGS)
00963       if (js->verbose)
00964         printf("jsplugin)   radius...\n");
00965 #endif
00966       *optflags |= MOLFILE_RADIUS;
00967       fio_fread(fltbuf, js->natoms * sizeof(float), 1, js->fd);
00968       if (js->reverseendian)
00969         swap4_aligned(fltbuf, js->natoms);
00970       for (i=0; i<js->natoms; i++) {
00971         atoms[i].radius = fltbuf[i];
00972       }    
00973     }
00974 
00975     if (fltbuf != NULL) {
00976       free(fltbuf);
00977       fltbuf=NULL;
00978     }
00979 
00980     /*
00981      * read in optional integer data blocks
00982      */ 
00983     if (js->optflags & JSOPT_ATOMICNUMBER)
00984       intbuf = (int *) malloc(js->natoms * sizeof(int));
00985 
00986     if (js->optflags & JSOPT_ATOMICNUMBER) { 
00987 #if defined(INFOMSGS)
00988       if (js->verbose)
00989         printf("jsplugin)   atomic number...\n");
00990 #endif
00991       *optflags |= MOLFILE_ATOMICNUMBER;
00992       fio_fread(intbuf, js->natoms * sizeof(int), 1, js->fd);
00993       if (js->reverseendian)
00994         swap4_aligned(intbuf, js->natoms);
00995       for (i=0; i<js->natoms; i++) {
00996         atoms[i].atomicnumber = intbuf[i];
00997       }    
00998     }
00999 
01000     if (intbuf != NULL) {
01001       free(intbuf);
01002       intbuf = NULL;
01003     }
01004 
01005 
01006     /*
01007      * read in bonds and fractional bond orders
01008      */ 
01009     if (js->optflags & JSOPT_BONDS) {
01010       fio_fread(&js->nbonds, sizeof(int), 1, js->fd);
01011       if (js->reverseendian)
01012         swap4_aligned(&js->nbonds, 1);
01013 #if defined(INFOMSGS)
01014       if (js->verbose)
01015         printf("jsplugin)   %d bonds...\n", js->nbonds);
01016 #endif
01017 
01018       js->bondfrom = (int *) malloc(js->nbonds * sizeof(int));
01019       js->bondto = (int *) malloc(js->nbonds * sizeof(int));
01020       fio_fread(js->bondfrom, js->nbonds * sizeof(int), 1, js->fd);
01021       fio_fread(js->bondto, js->nbonds * sizeof(int), 1, js->fd);
01022       if (js->reverseendian) {
01023         swap4_aligned(js->bondfrom, js->nbonds);
01024         swap4_aligned(js->bondto, js->nbonds);
01025       }
01026 
01027       if (js->optflags & JSOPT_BONDORDERS) {
01028 #if defined(INFOMSGS)
01029         if (js->verbose)
01030           printf("jsplugin)   bond orders...\n");
01031 #endif
01032         js->bondorders = (float *) malloc(js->nbonds * sizeof(float));
01033         fio_fread(js->bondorders, js->nbonds * sizeof(float), 1, js->fd);
01034         if (js->reverseendian)
01035           swap4_aligned(js->bondorders, js->nbonds);
01036       }
01037     }
01038 
01039     if (js->optflags & JSOPT_ANGLES) {
01040       fio_fread(&js->numangles, sizeof(int), 1, js->fd);
01041       if (js->reverseendian)
01042         swap4_aligned(&js->numangles, 1);
01043 #if defined(INFOMSGS)
01044       if (js->verbose)
01045         printf("jsplugin)   %d angles...\n", js->numangles);
01046 #endif
01047       js->angles = (int *) malloc(3L * js->numangles * sizeof(int));
01048       fio_fread(js->angles, sizeof(int)*3L*js->numangles, 1, js->fd);
01049       if (js->reverseendian)
01050         swap4_aligned(js->angles, 3L*js->numangles);
01051 
01052       fio_fread(&js->numdihedrals, sizeof(int), 1, js->fd);
01053       if (js->reverseendian)
01054         swap4_aligned(&js->numdihedrals, 1);
01055 #if defined(INFOMSGS)
01056       if (js->verbose)
01057         printf("jsplugin)   %d dihedrals...\n", js->numdihedrals);
01058 #endif
01059       js->dihedrals = (int *) malloc(4L * js->numdihedrals * sizeof(int));
01060       fio_fread(js->dihedrals, sizeof(int)*4L*js->numdihedrals, 1, js->fd);
01061       if (js->reverseendian)
01062         swap4_aligned(js->dihedrals, 4L*js->numdihedrals);
01063 
01064       fio_fread(&js->numimpropers, sizeof(int), 1, js->fd);
01065       if (js->reverseendian)
01066         swap4_aligned(&js->numimpropers, 1);
01067       js->impropers = (int *) malloc(4L * js->numimpropers * sizeof(int));
01068 #if defined(INFOMSGS)
01069       if (js->verbose)
01070         printf("jsplugin)   %d impropers...\n", js->numimpropers);
01071 #endif
01072       fio_fread(js->impropers, sizeof(int)*4L*js->numimpropers, 1, js->fd);
01073       if (js->reverseendian)
01074         swap4_aligned(js->impropers, 4L*js->numimpropers);
01075     }
01076 
01077     if (js->optflags & JSOPT_CTERMS) {
01078       fio_fread(&js->numcterms, sizeof(int), 1, js->fd);
01079       if (js->reverseendian)
01080         swap4_aligned(&js->numcterms, 1);
01081       js->cterms = (int *) malloc(8L * js->numcterms * sizeof(int));
01082 #if defined(INFOMSGS)
01083       if (js->verbose)
01084         printf("jsplugin)   %d cterms...\n", js->numcterms);
01085 #endif
01086       fio_fread(js->cterms, sizeof(int)*8L*js->numcterms, 1, js->fd);
01087       if (js->reverseendian)
01088         swap4_aligned(js->cterms, 8L*js->numcterms);
01089     }
01090 
01091 #if defined(INFOMSGS)
01092     if (js->verbose) {
01093       printf("jsplugin) final optflags: %08x\n", *optflags);
01094       printf("jsplugin) structure information complete\n");
01095     }
01096 #endif
01097 
01098     /* record the file offset for the first timestep */
01099     js_calc_timestep_blocking_info(js);
01100 
01101     return MOLFILE_SUCCESS;
01102   }
01103 
01104 #if defined(INFOMSGS)
01105   if (js->verbose)
01106     printf("jsplugin) no structure information available\n");
01107 #endif
01108 
01109   /* record the file offset for the first timestep */
01110   js_calc_timestep_blocking_info(js);
01111 
01112   /* else, we have no structure information */
01113   return MOLFILE_NOSTRUCTUREDATA;
01114 }
01115 
01116 
01117 static int read_js_bonds(void *v, int *nbonds, int **fromptr, int **toptr, 
01118                          float **bondorder, int **bondtype, 
01119                          int *nbondtypes, char ***bondtypename) {
01120   jshandle *js = (jshandle *)v;
01121 
01122   *nbonds = 0;
01123   *fromptr = NULL;
01124   *toptr = NULL;
01125   *bondorder = NULL;
01126   *bondtype = NULL;
01127   *nbondtypes = 0;
01128   *bondtypename = NULL;
01129 
01130   if (js->optflags & JSOPT_BONDS) {
01131     *nbonds = js->nbonds;
01132     *fromptr = js->bondfrom;
01133     *toptr = js->bondto;
01134 
01135     if (js->optflags & JSOPT_BONDORDERS) {
01136       *bondorder = js->bondorders;
01137     }
01138   }
01139 
01140   return MOLFILE_SUCCESS;
01141 }
01142 
01143 #if vmdplugin_ABIVERSION > 14
01144 static int read_js_angles(void *v, int *numangles, int **angles, 
01145                           int **angletypes, int *numangletypes, 
01146                           char ***angletypenames, int *numdihedrals,
01147                           int **dihedrals, int **dihedraltypes, 
01148                           int *numdihedraltypes, char ***dihedraltypenames,
01149                           int *numimpropers, int **impropers, 
01150                           int **impropertypes, int *numimpropertypes, 
01151                           char ***impropertypenames, int *numcterms, 
01152                           int **cterms, int *ctermcols, int *ctermrows) {
01153   jshandle *js = (jshandle *)v;
01154 
01155   /* initialize data to zero */
01156   *numangles         = 0;
01157   *angles            = NULL;
01158   *angletypes        = NULL;
01159   *numangletypes     = 0;
01160   *angletypenames    = NULL;
01161   *numdihedrals      = 0;
01162   *dihedrals         = NULL;
01163   *dihedraltypes     = NULL;
01164   *numdihedraltypes  = 0;
01165   *dihedraltypenames = NULL;
01166   *numimpropers      = 0;
01167   *impropers         = NULL;
01168   *impropertypes     = NULL;
01169   *numimpropertypes  = 0;
01170   *impropertypenames = NULL;
01171   *numcterms         = 0;
01172   *cterms            = NULL;
01173   *ctermrows         = 0;
01174   *ctermcols         = 0;
01175 
01176   *numangles = js->numangles;
01177   *angles = js->angles;
01178 
01179   *numdihedrals = js->numdihedrals;
01180   *dihedrals = js->dihedrals;
01181 
01182   *numimpropers = js->numimpropers;
01183   *impropers = js->impropers;
01184 
01185   *numcterms = js->numcterms;
01186   *cterms = js->cterms;
01187   *ctermcols = 0;
01188   *ctermrows = 0;
01189 
01190   return MOLFILE_SUCCESS;
01191 }
01192 #else
01193 static int read_js_angles(void *v,
01194                int *numangles,    int **angles,    double **angleforces,
01195                int *numdihedrals, int **dihedrals, double **dihedralforces,
01196                int *numimpropers, int **impropers, double **improperforces,
01197                int *numcterms,    int **cterms,
01198                int *ctermcols,    int *ctermrows,  double **ctermforces) {
01199   jshandle *js = (jshandle *)v;
01200 
01201   *numangles = js->numangles;
01202   *angles = js->angles;
01203   *angleforces = NULL;
01204 
01205   *numdihedrals = js->numdihedrals;
01206   *dihedrals = js->dihedrals;
01207   *dihedralforces = NULL;
01208 
01209   *numimpropers = js->numimpropers;
01210   *impropers = js->impropers;
01211   *improperforces = NULL;
01212 
01213   *numcterms = js->numcterms;
01214   *cterms = js->cterms;
01215   *ctermcols = 0;
01216   *ctermrows = 0;
01217   *ctermforces = NULL;
01218 
01219   return MOLFILE_SUCCESS;
01220 }
01221 #endif
01222 
01223 #endif
01224 
01225 
01226 #if 1 
01227 // XXX prototypical out-of-core trajectory analysis API
01228 static int read_js_timestep_index_offsets(void *v, int natoms, 
01229                                           ptrdiff_t frameindex,
01230                                           int firstatom, int numatoms,
01231                                           fio_fd *directio_fd,
01232                                           ptrdiff_t *startoffset,
01233                                           ptrdiff_t *fileoffset,
01234                                           ptrdiff_t *readlen) {
01235   jshandle *js = (jshandle *)v;
01236   fio_size_t framelen;
01237 
01238 #if JSMAJORVERSION > 1
01239   /* If we haven't yet read (or skipped) the structure data, then we    */
01240   /* need to begin by skipping past it before we try to read the        */
01241   /* first timestep.  In the case of files with block-aligned timesteps,*/
01242   /* this will also get our file pointer to the right block-aligned     */
01243   /* location.                                                          */
01244   if (!js->parsed_structure)
01245     read_js_structure(v, NULL, NULL);
01246 #endif
01247 
01248   /* compute total read/seek size of timestep */
01249   framelen = js->ts_crd_padsz + js->ts_ucell_padsz;
01250 
01251   if (directio_fd != NULL)
01252     *directio_fd = js->directio_fd;
01253 
01254   /* compute file offset for requested timestep */
01255   if (fileoffset != NULL)
01256     *fileoffset = (frameindex * framelen) + js->ts_file_offset;
01257 
01258   /* compute startoffset for first requested atom */
01259   if (startoffset != NULL)
01260     *startoffset = firstatom * 3L * sizeof(float);
01261  
01262   /* compute required read size */ 
01263   if (readlen != NULL)
01264     *readlen = framelen;
01265 
01266   return MOLFILE_SUCCESS;
01267 }
01268 
01269 
01270 #if 0
01271 static int read_js_timestep_index(void *v, int natoms, 
01272                                   ptrdiff_t frameindex,
01273                                   molfile_timestep_t *ts) {
01274 }
01275 #endif
01276 
01277 #endif
01278 
01279 
01280 
01281 static int read_js_timestep(void *v, int natoms, molfile_timestep_t *ts) {
01282   jshandle *js = (jshandle *)v;
01283   fio_size_t framelen;
01284 
01285 #if JSMAJORVERSION > 1
01286   /* If we haven't yet read (or skipped) the structure data, then we    */
01287   /* need to begin by skipping past it before we try to read the        */
01288   /* first timestep.  In the case of files with block-aligned timesteps,*/
01289   /* this will also get our file pointer to the right block-aligned     */
01290   /* location.                                                          */
01291   if (!js->parsed_structure)
01292     read_js_structure(v, NULL, NULL);
01293 #endif
01294 
01295   /* compute total read/seek size of timestep */
01296   framelen = js->ts_crd_padsz + js->ts_ucell_padsz;
01297 
01298   /* if we have a valid ts pointer, read the timestep, otherwise skip it */ 
01299   if (ts != NULL) {
01300     fio_size_t readlen=0;
01301     fio_iovec iov[2];
01302 
01303     /* set unit cell pointer to the TS block-aligned buffer area */
01304     double *unitcell = (double *) js->directio_ucell_blkbuf;
01305 
01306     unitcell[0] = unitcell[2] = unitcell[5] = 1.0f;
01307     unitcell[1] = unitcell[3] = unitcell[4] = 90.0f;
01308 
01309 #if defined(ENABLEJSSHORTREADS)
01310     /* test code for an implementation that does short reads that */
01311     /* skip bulk solvent, useful for faster loading of very large */
01312     /* structures                                                 */
01313     if (getenv("VMDJSMAXATOMIDX") != NULL) {
01314       fio_size_t bszmask;
01315       ptrdiff_t maxatompadsz, skipatompadsz;
01316 
01317       ptrdiff_t maxatomidx = atoi(getenv("VMDJSMAXATOMIDX"));
01318       if (maxatomidx < 0)
01319         maxatomidx = 0;
01320       if (maxatomidx >= js->natoms)
01321         maxatomidx = js->natoms - 1;
01322 
01323       /* pad max read to the start of the next block  */
01324       bszmask = js->directio_block_size - 1;
01325       maxatompadsz = ((maxatomidx*3L*sizeof(float)) + bszmask) & (~bszmask);
01326       skipatompadsz = js->ts_crd_padsz - maxatompadsz;
01327 
01328       readlen=0;
01329       if (js->directio_enabled) {
01330         if (fio_fread(ts->coords, maxatompadsz, 1, js->directio_fd) == 1)
01331           readlen = maxatompadsz;
01332         if (fio_fseek(js->directio_fd, skipatompadsz, FIO_SEEK_CUR) == 0)
01333           readlen += skipatompadsz;
01334         if (fio_fread(unitcell, js->ts_ucell_padsz, 1, js->directio_fd) == 1)
01335           readlen += js->ts_ucell_padsz;
01336       } else {
01337         if (fio_fread(ts->coords, maxatompadsz, 1, js->fd) == 1)
01338           readlen = maxatompadsz;
01339         if (fio_fseek(js->fd, skipatompadsz, FIO_SEEK_CUR) == 0)
01340           readlen += skipatompadsz;
01341         if (fio_fread(unitcell, js->ts_ucell_padsz, 1, js->fd) == 1)
01342           readlen += js->ts_ucell_padsz;
01343       }
01344 
01345 #if 0
01346       /* clear all non-read atom coords to zeros */
01347       memset(ts->coords+3L*maxatomidx,0,3L*sizeof(float)*(js->natoms-maxatomidx));
01348 #endif
01349 
01350     }  else {
01351 #endif
01352  
01353     /* setup the I/O vector */
01354     iov[0].iov_base = (fio_caddr_t) ts->coords;   /* read coordinates    */
01355     iov[1].iov_base = (fio_caddr_t) unitcell;     /* read PBC unit cell  */
01356 
01357     if (js->directio_enabled) {
01358       iov[0].iov_len  = js->ts_crd_padsz;
01359       iov[1].iov_len  = js->ts_ucell_padsz;
01360     } else {
01361       iov[0].iov_len  = js->ts_crd_sz;
01362       iov[1].iov_len  = js->ts_ucell_sz;
01363     }
01364    
01365 #if 1
01366     /* Use fall-back code instead of readv():                            */
01367     /*  Some platforms implement readv() as user level code in libc,     */
01368     /*  and due to POSIX atomicity requirements for readv()/writev(),    */
01369     /*  they may copy data to internal temp buffers, which can kill      */
01370     /*  performance, and in cases when doing single I/O ops on large,    */
01371     /*  buffers, e.g. > 2GB, can fail with shorts reads or writes...     */
01372     /*  On such platforms it is best to avoid using readv()/writev()...  */
01373     {
01374       int readcnt = 0;
01375       readlen = 0;
01376       if (js->directio_enabled) {
01377         readcnt =  fio_fread(iov[0].iov_base, iov[0].iov_len, 1, js->directio_fd);
01378         readcnt += fio_fread(iov[1].iov_base, iov[1].iov_len, 1, js->directio_fd);
01379       } else {
01380         fio_size_t seeklen=0;
01381       
01382         readcnt =  fio_fread(iov[0].iov_base, iov[0].iov_len, 1, js->fd);
01383         seeklen = js->ts_crd_padsz - js->ts_crd_sz;
01384         if (seeklen > 0)
01385           fio_fseek(js->fd, seeklen, FIO_SEEK_CUR);
01386         readcnt += fio_fread(iov[1].iov_base, iov[1].iov_len, 1, js->fd);
01387         seeklen = js->ts_ucell_padsz - js->ts_ucell_sz;
01388         if (seeklen > 0)
01389           fio_fseek(js->fd, seeklen, FIO_SEEK_CUR);
01390       }
01391 
01392       /* if both records read correctly, then the reads are okay */
01393       if (readcnt == 2)
01394         readlen = framelen;
01395     }
01396 #else
01397     /* Do all of the reads with a single syscall, for peak efficiency.   */
01398     /* On smart kernels, readv() causes only one context switch, and     */
01399     /* can effeciently scatter the reads to the various buffers.         */
01400     if (js->directio_enabled) {
01401       readlen = fio_readv(js->directio_fd, &iov[0], 2); 
01402     } else {
01403       // XXX we can't use readv() when not using direct I/O since we 
01404       // can't make intervening seek calls required if the caller
01405       // doesn't provide appropriate buffers.
01406       // readlen = fio_readv(js->fd, &iov[0], 2); 
01407 
01408       fio_size_t seeklen=0;
01409       readcnt =  fio_fread(iov[0].iov_base, iov[0].iov_len, 1, js->fd);
01410       seeklen = js->ts_crd_padsz - js->ts_crd_sz;
01411       if (seeklen > 0)
01412         fio_fseek(js->fd, seeklen, FIO_SEEK_CUR);
01413       readcnt += fio_fread(iov[1].iov_base, iov[1].iov_len, 1, js->fd);
01414       seeklen = js->ts_ucell_padsz - js->ts_ucell_sz;
01415       if (seeklen > 0)
01416         fio_fseek(js->fd, seeklen, FIO_SEEK_CUR);
01417     }
01418 #endif
01419 
01420 #if defined(ENABLEJSSHORTREADS)
01421    }
01422 #endif 
01423  
01424     /* check the number of read bytes versus what we expected */
01425     if (readlen != framelen) {
01426       if (readlen < 0) {
01427         perror("jsplugin) fio_readv(): ");
01428       } else if (readlen != 0) {
01429         printf("jsplugin) mismatched read: %td, expected %td\n", 
01430                (ptrdiff_t) readlen, (ptrdiff_t) framelen);
01431       }
01432 
01433       return MOLFILE_EOF;
01434     }
01435 
01436     /* perform byte swapping if necessary */
01437     if (js->reverseendian) {
01438       swap4_aligned(ts->coords, js->natoms * 3L);
01439       swap8_aligned(unitcell, 6);
01440     }
01441 
01442     /* copy unit cell values into VMD */
01443     ts->A = unitcell[0];
01444     ts->B = unitcell[1];
01445     ts->C = unitcell[2];
01446     ts->alpha = 90.0 - asin(unitcell[3]) * 90.0 / M_PI_2;
01447     ts->beta  = 90.0 - asin(unitcell[4]) * 90.0 / M_PI_2;
01448     ts->gamma = 90.0 - asin(unitcell[5]) * 90.0 / M_PI_2;
01449   } else {
01450     /* skip this frame, seek to the next frame */
01451     if (js->directio_enabled) {
01452       if (fio_fseek(js->directio_fd, framelen, FIO_SEEK_CUR)) 
01453         return MOLFILE_EOF;
01454     } else {
01455       if (fio_fseek(js->fd, framelen, FIO_SEEK_CUR)) 
01456         return MOLFILE_EOF;
01457     }
01458   }
01459  
01460   return MOLFILE_SUCCESS;
01461 }
01462 
01463 
01464 static void close_js_read(void *v) {
01465   jshandle *js = (jshandle *)v;
01466   fio_fclose(js->fd);
01467 
01468 #if JSMAJORVERSION > 1
01469   if (js->path)
01470     free(js->path);
01471 
01472   if (js->directio_enabled)
01473     fio_fclose(js->directio_fd);
01474 
01475   if (js->directio_ucell_ptr)
01476     free(js->directio_ucell_ptr);
01477 
01478   if (js->bondfrom)
01479     free(js->bondfrom);
01480   if (js->bondto)
01481     free(js->bondto);
01482   if (js->bondorders)
01483     free(js->bondorders);
01484 
01485   /* free angle data */
01486   if (js->angles != NULL)
01487     free(js->angles);
01488   if (js->dihedrals != NULL)
01489     free(js->dihedrals);
01490   if (js->impropers != NULL)
01491     free(js->impropers);
01492   if (js->cterms)
01493     free(js->cterms);
01494 #endif
01495 
01496   free(js);
01497 }
01498 
01499 
01500 static void *open_js_write(const char *path, const char *filetype, int natoms) {
01501   jshandle *js;
01502 
01503   js = (jshandle *) malloc(sizeof(jshandle));
01504   memset(js, 0, sizeof(jshandle));
01505 #if JSMAJORVERSION > 1
01506   js->parsed_structure=0;
01507   js->directio_block_size=1;
01508   js->directio_ucell_ptr = NULL;
01509   js->directio_ucell_blkbuf = NULL;
01510 
01511   js->directio_enabled=0;
01512   js->ts_file_offset=0;
01513   js->ts_crd_sz=0;
01514   js->ts_ucell_sz=0;
01515   js->ts_crd_padsz=0;
01516   js->ts_ucell_padsz=0;
01517 #endif
01518 
01519   if (fio_open(path, FIO_WRITE, &js->fd) < 0) {
01520     printf("jsplugin) Could not open file %s for writing\n", path);
01521     free(js);
01522     return NULL;
01523   }
01524 
01525   js->natoms = natoms;
01526   js->with_unitcell = 1;
01527 
01528   /* emit header information */
01529   fio_write_str(js->fd, JSHEADERSTRING);
01530   fio_write_int32(js->fd, JSMAGICNUMBER);
01531   fio_write_int32(js->fd, JSENDIANISM);
01532   fio_write_int32(js->fd, JSMAJORVERSION);
01533   fio_write_int32(js->fd, JSMINORVERSION);
01534 
01535   /* write number of atoms */
01536   fio_write_int32(js->fd, natoms);
01537 
01538   /* write number of frames, to be updated later */
01539   js->nframes = 0;
01540   fio_write_int32(js->fd, js->nframes);
01541 
01542   return js;
01543 }
01544 
01545 
01546 #if JSMAJORVERSION > 1
01547 
01548 static int write_js_structure(void *mydata, int optflags,
01549                               const molfile_atom_t *atoms) {
01550   jshandle *js = (jshandle *) mydata;
01551   ptrdiff_t i;
01552 
01553   /* use block-based I/O by default when writing structures larger */
01554   /* than JSBLOCKIO_THRESH atoms, or when directed by the user     */
01555   js_blockio_check_and_set(js);
01556 
01557   js->optflags |= JSOPT_STRUCTURE;
01558 
01559   if (optflags & MOLFILE_OCCUPANCY)
01560     js->optflags |= JSOPT_OCCUPANCY;
01561 
01562   if (optflags & MOLFILE_BFACTOR)
01563     js->optflags |= JSOPT_BFACTOR;
01564 
01565   if (optflags & MOLFILE_BFACTOR)
01566     js->optflags |= JSOPT_BFACTOR;
01567 
01568   if (optflags & MOLFILE_MASS)
01569     js->optflags |= JSOPT_MASS;
01570 
01571   if (optflags & MOLFILE_CHARGE)
01572     js->optflags |= JSOPT_CHARGE;
01573  
01574   if (optflags & MOLFILE_RADIUS)
01575     js->optflags |= JSOPT_RADIUS;
01576 
01577   if (optflags & MOLFILE_ATOMICNUMBER)
01578     js->optflags |= JSOPT_ATOMICNUMBER;
01579 
01580   /* write flags data to the file */
01581   fio_write_int32(js->fd, js->optflags); 
01582 printf("jsplugin) writing option flags: %0x08x\n", js->optflags);
01583 
01584   /* Check to see if block-based trajectory I/O is used  */
01585   /* and write out the block size for this file.         */
01586   if (js->optflags & JSOPT_TS_BLOCKIO) {
01587     fio_fwrite(&js->directio_block_size, sizeof(int), 1, js->fd);
01588     printf("jsplugin) Block-based I/O enabled: block size %d bytes\n", 
01589            js->directio_block_size);
01590   }
01591 
01592 printf("jsplugin) writing structure...\n");
01593   /* determine whether or not this file contains structure info or not */
01594   if (js->optflags & JSOPT_STRUCTURE) {
01595     int numatomnames, numatomtypes, numresnames, numsegids, numchains;
01596     char **atomnames = NULL;
01597     char **atomtypes = NULL;
01598     char **resnames = NULL;
01599     char **segids = NULL;
01600     char **chains = NULL;
01601     short *shortbuf = NULL; /* temp buf for encoding atom records */
01602     int *intbuf = NULL;     /* temp buf for encoding atom records */
01603     float *fltbuf = NULL;   /* temp buf for encoding atom records */
01604 
01605     hash_t tmphash;         /* temporary hash table */
01606     hash_t atomnamehash;
01607     hash_t atomtypehash;
01608     hash_t resnamehash;
01609     hash_t segidhash;
01610     hash_t chainhash;
01611     int hashcnt;
01612 
01613 
01614 printf("jsplugin) counting atom names, types, etc...\n");
01615     /* generate hash tables to count the number of unique strings */
01616     hash_init(&tmphash, 127);
01617     for (i=0; i<js->natoms; i++)
01618       hash_insert(&tmphash, atoms[i].name, 0);
01619     numatomnames = hash_entries(&tmphash);
01620     hash_destroy(&tmphash);
01621 
01622     hash_init(&tmphash, 127);
01623     for (i=0; i<js->natoms; i++)
01624       hash_insert(&tmphash, atoms[i].type, 0);
01625     numatomtypes = hash_entries(&tmphash);
01626     hash_destroy(&tmphash);
01627 
01628     hash_init(&tmphash, 127);
01629     for (i=0; i<js->natoms; i++)
01630       hash_insert(&tmphash, atoms[i].resname, 0);
01631     numresnames = hash_entries(&tmphash);
01632     hash_destroy(&tmphash);
01633 
01634     hash_init(&tmphash, 127);
01635     for (i=0; i<js->natoms; i++)
01636       hash_insert(&tmphash, atoms[i].segid, 0);
01637     numsegids = hash_entries(&tmphash);
01638     hash_destroy(&tmphash);
01639 
01640     hash_init(&tmphash, 127);
01641     for (i=0; i<js->natoms; i++)
01642       hash_insert(&tmphash, atoms[i].chain, 0);
01643     numchains = hash_entries(&tmphash);
01644     hash_destroy(&tmphash);
01645  
01646 printf("jsplugin) writing unique string counts...\n");
01647 printf("jsplugin) %d %d %d %d %d\n",
01648        numatomnames, numatomtypes, numresnames, numsegids, numchains);
01649 
01650     /* write block of name string table sizes */
01651     fio_write_int32(js->fd, numatomnames); 
01652     fio_write_int32(js->fd, numatomtypes); 
01653     fio_write_int32(js->fd, numresnames);
01654     fio_write_int32(js->fd, numsegids);
01655     fio_write_int32(js->fd, numchains); 
01656 
01657 printf("jsplugin) writing string tables...\n");
01658 
01659     atomnames = (char **) malloc(numatomnames * sizeof(char *));
01660     atomtypes = (char **) malloc(numatomtypes * sizeof(char *));
01661     resnames = (char **) malloc(numresnames * sizeof(char *));
01662     segids = (char **) malloc(numsegids * sizeof(char *));
01663     chains = (char **) malloc(numchains * sizeof(char *));
01664 
01665 printf("jsplugin)   atom names...\n");
01666     /* generate and write out the string tables */
01667     hash_init(&atomnamehash, 127);
01668     for (hashcnt=0,i=0; i<js->natoms; i++) {
01669       /* add a new string table entry for hash inserts that don't yet exist */
01670       if (hash_insert(&atomnamehash, atoms[i].name, hashcnt) == HASH_FAIL) {
01671         atomnames[hashcnt] = (char *) calloc(1, 16L * sizeof(char));
01672         strcpy(atomnames[hashcnt], atoms[i].name);
01673         hashcnt++;
01674       }
01675     }
01676     for (i=0; i<numatomnames; i++) {
01677       fio_fwrite(atomnames[i], 16L * sizeof(char), 1, js->fd);
01678     }
01679 
01680 
01681 printf("jsplugin)   atom types...\n");
01682     hash_init(&atomtypehash, 127);
01683     for (hashcnt=0,i=0; i<js->natoms; i++) {
01684       /* add a new string table entry for hash inserts that don't yet exist */
01685       if (hash_insert(&atomtypehash, atoms[i].type, hashcnt) == HASH_FAIL) {
01686         atomtypes[hashcnt] = (char *) calloc(1, 16L * sizeof(char));
01687         strcpy(atomtypes[hashcnt], atoms[i].type);
01688         hashcnt++;
01689       }
01690     }
01691     for (i=0; i<numatomtypes; i++) {
01692       fio_fwrite(atomtypes[i], 16L * sizeof(char), 1, js->fd);
01693     }
01694 
01695 
01696 printf("jsplugin)   residue names...\n");
01697     hash_init(&resnamehash, 127);
01698     for (hashcnt=0,i=0; i<js->natoms; i++) {
01699       /* add a new string table entry for hash inserts that don't yet exist */
01700       if (hash_insert(&resnamehash, atoms[i].resname, hashcnt) == HASH_FAIL) {
01701         resnames[hashcnt] = (char *) calloc(1, 8L * sizeof(char));
01702         strcpy(resnames[hashcnt], atoms[i].resname);
01703         hashcnt++;
01704       }
01705     }
01706     for (i=0; i<numresnames; i++) {
01707       fio_fwrite(resnames[i], 8L * sizeof(char), 1, js->fd);
01708     }
01709 
01710 
01711 printf("jsplugin)   segment names...\n");
01712     hash_init(&segidhash, 127);
01713     for (hashcnt=0,i=0; i<js->natoms; i++) {
01714       /* add a new string table entry for hash inserts that don't yet exist */
01715       if (hash_insert(&segidhash, atoms[i].segid, hashcnt) == HASH_FAIL) {
01716         segids[hashcnt] = (char *) calloc(1, 8L * sizeof(char));
01717         strcpy(segids[hashcnt], atoms[i].segid);
01718         hashcnt++;
01719       }
01720     }
01721     for (i=0; i<numsegids; i++) {
01722       fio_fwrite(segids[i], 8L * sizeof(char), 1, js->fd);
01723     }
01724 
01725 
01726 printf("jsplugin)   chain names...\n");
01727     hash_init(&chainhash, 127);
01728     for (hashcnt=0,i=0; i<js->natoms; i++) {
01729       /* add a new string table entry for hash inserts that don't yet exist */
01730       if (hash_insert(&chainhash, atoms[i].chain, hashcnt) == HASH_FAIL) {
01731         chains[hashcnt] = (char *) calloc(1, 2L * sizeof(char));
01732         strcpy(chains[hashcnt], atoms[i].chain);
01733         hashcnt++;
01734       }
01735     }
01736     for (i=0; i<numchains; i++) {
01737       fio_fwrite(chains[i], 2L * sizeof(char), 1, js->fd);
01738     }
01739 
01740 
01741 printf("jsplugin) writing numeric field tables...\n");
01742     /* write out all of the atom fields */
01743     shortbuf = (short *) malloc(js->natoms * sizeof(short));
01744 
01745     /* write out atom names */
01746     for (i=0; i<js->natoms; i++) {
01747       shortbuf[i] = hash_lookup(&atomnamehash, atoms[i].name);
01748     }    
01749     fio_fwrite(shortbuf, js->natoms * sizeof(short), 1, js->fd);
01750 
01751     /* write out atom types */
01752     for (i=0; i<js->natoms; i++) {
01753       shortbuf[i] = hash_lookup(&atomtypehash, atoms[i].type);
01754     }    
01755     fio_fwrite(shortbuf, js->natoms * sizeof(short), 1, js->fd);
01756 
01757     /* write out resnames */
01758     for (i=0; i<js->natoms; i++) {
01759       shortbuf[i] = hash_lookup(&resnamehash, atoms[i].resname);
01760     }    
01761     fio_fwrite(shortbuf, js->natoms * sizeof(short), 1, js->fd);
01762     
01763     /* write out segids */
01764     for (i=0; i<js->natoms; i++) {
01765       shortbuf[i] = hash_lookup(&segidhash, atoms[i].segid);
01766     }    
01767     fio_fwrite(shortbuf, js->natoms * sizeof(short), 1, js->fd);
01768 
01769     /* write out chains */
01770     for (i=0; i<js->natoms; i++) {
01771       shortbuf[i] = hash_lookup(&chainhash, atoms[i].chain);
01772     }    
01773     fio_fwrite(shortbuf, js->natoms * sizeof(short), 1, js->fd);
01774 
01775     if (shortbuf != NULL) {
01776       free(shortbuf);
01777       shortbuf=NULL;
01778     }
01779 
01780     /* done with hash tables */
01781     hash_destroy(&atomnamehash);
01782     hash_destroy(&atomtypehash);
01783     hash_destroy(&resnamehash);
01784     hash_destroy(&segidhash);
01785     hash_destroy(&chainhash);
01786 
01787 
01788     /* 
01789      * write out integer data blocks 
01790      */
01791     intbuf = (int *) malloc(js->natoms * sizeof(int));
01792 
01793 printf("jsplugin)   residue indices...\n");
01794     /* write out resid */
01795     for (i=0; i<js->natoms; i++) {
01796       intbuf[i] = atoms[i].resid;
01797     }    
01798     fio_fwrite(intbuf, js->natoms * sizeof(int), 1, js->fd);
01799      
01800     if (intbuf != NULL) {
01801       free(intbuf);
01802       intbuf = NULL;
01803     }
01804 
01805 printf("jsplugin) writing optional per-atom tables...\n");
01806     /*
01807      * write out optional single-precision float data blocks
01808      */ 
01809     if (js->optflags & (JSOPT_OCCUPANCY | JSOPT_BFACTOR | 
01810         JSOPT_MASS | JSOPT_RADIUS | JSOPT_CHARGE)) 
01811       fltbuf = (float *) malloc(js->natoms * sizeof(float));
01812 
01813     /* write out optional data if it exists */
01814 
01815     if (js->optflags & JSOPT_OCCUPANCY) {
01816 printf("jsplugin)   writing occupancy...\n");
01817       for (i=0; i<js->natoms; i++) {
01818         fltbuf[i] = atoms[i].occupancy;
01819       }    
01820       fio_fwrite(fltbuf, js->natoms * sizeof(float), 1, js->fd);
01821     }
01822 
01823     if (js->optflags & JSOPT_BFACTOR) {
01824 printf("jsplugin)   writing bfactor...\n");
01825       for (i=0; i<js->natoms; i++) {
01826         fltbuf[i] = atoms[i].bfactor;
01827       }    
01828       fio_fwrite(fltbuf, js->natoms * sizeof(float), 1, js->fd);
01829     }
01830 
01831     if (js->optflags & JSOPT_MASS) { 
01832 printf("jsplugin)   writing mass...\n");
01833       for (i=0; i<js->natoms; i++) {
01834         fltbuf[i] = atoms[i].mass;
01835       }    
01836       fio_fwrite(fltbuf, js->natoms * sizeof(float), 1, js->fd);
01837     }
01838 
01839     if (js->optflags & JSOPT_CHARGE) { 
01840 printf("jsplugin)   writing charge...\n");
01841       for (i=0; i<js->natoms; i++) {
01842         fltbuf[i] = atoms[i].charge;
01843       }    
01844       fio_fwrite(fltbuf, js->natoms * sizeof(float), 1, js->fd);
01845     }
01846 
01847     if (js->optflags & JSOPT_RADIUS) { 
01848 printf("jsplugin)   writing radius...\n");
01849       for (i=0; i<js->natoms; i++) {
01850         fltbuf[i] = atoms[i].radius;
01851       }    
01852       fio_fwrite(fltbuf, js->natoms * sizeof(float), 1, js->fd);
01853     }
01854 
01855     if (fltbuf != NULL) {
01856       free(fltbuf);
01857       fltbuf=NULL;
01858     }
01859 
01860 
01861     /*
01862      * write out optional integer data blocks
01863      */ 
01864     if (js->optflags & JSOPT_ATOMICNUMBER)
01865       intbuf = (int *) malloc(js->natoms * sizeof(int));
01866 
01867     if (js->optflags & JSOPT_ATOMICNUMBER) { 
01868 printf("jsplugin)   writing atomic number...\n");
01869       for (i=0; i<js->natoms; i++) {
01870         intbuf[i] = atoms[i].atomicnumber;
01871       }    
01872       fio_fwrite(intbuf, js->natoms * sizeof(int), 1, js->fd);
01873     }
01874 
01875     if (intbuf != NULL) {
01876       free(intbuf);
01877       intbuf = NULL;
01878     }
01879 
01880 
01881     /*
01882      * write out bonds and fractional bond orders
01883      */ 
01884     if (js->optflags & JSOPT_BONDS) {
01885 printf("jsplugin) writing bonds...\n");
01886       fio_fwrite(&js->nbonds, sizeof(int), 1, js->fd);
01887       fio_fwrite(js->bondfrom, js->nbonds * sizeof(int), 1, js->fd);
01888       fio_fwrite(js->bondto, js->nbonds * sizeof(int), 1, js->fd);
01889 
01890       if (js->optflags & JSOPT_BONDORDERS) {
01891 printf("jsplugin) writing bond orders...\n");
01892         fio_fwrite(js->bondorders, js->nbonds * sizeof(float), 1, js->fd);
01893       }
01894     }
01895 
01896     /*
01897      * write out angles/dihedrals/impropers/cross-terms
01898      */
01899     if (js->optflags & JSOPT_ANGLES) {
01900 printf("jsplugin) writing angles/dihedrals/impropers...\n");
01901       fio_fwrite(&js->numangles, sizeof(int), 1, js->fd);
01902       fio_fwrite(js->angles, sizeof(int)*3L*js->numangles, 1, js->fd);
01903 
01904       fio_fwrite(&js->numdihedrals, sizeof(int), 1, js->fd);
01905       fio_fwrite(js->dihedrals, sizeof(int)*4L*js->numdihedrals, 1, js->fd);
01906 
01907       fio_fwrite(&js->numimpropers, sizeof(int), 1, js->fd);
01908       fio_fwrite(js->impropers, sizeof(int)*4L*js->numimpropers, 1, js->fd);
01909     }
01910     if (js->optflags & JSOPT_CTERMS) {
01911 printf("jsplugin) writing cross-terms\n");
01912       fio_fwrite(&js->numcterms, sizeof(int), 1, js->fd);
01913       fio_fwrite(js->cterms, sizeof(int)*8L*js->numcterms, 1, js->fd);
01914     }
01915 
01916     /* update the file offset for the first timestep */
01917     js_calc_timestep_blocking_info(js);
01918 
01919     return MOLFILE_SUCCESS;
01920   }
01921 
01922   /* update the file offset for the first timestep */
01923   js_calc_timestep_blocking_info(js);
01924 
01925   /* else, we have no structure information */
01926   return MOLFILE_NOSTRUCTUREDATA;
01927 }
01928 
01929 
01930 static int write_js_bonds(void *mydata, int nbonds, int *fromptr, int *toptr, 
01931                           float *bondorder,  int *bondtype, 
01932                           int nbondtypes, char **bondtypename) {
01933   jshandle *js = (jshandle *) mydata;
01934 
01935 #if defined(INFOMSGS)
01936     if (js->verbose) {
01937       printf("jsplugin) write_js_bonds():\n");
01938       printf("jsplugin) storing bond info for writing...\n");
01939       printf("jsplugin) %d %d\n", nbonds, nbondtypes);
01940     }
01941 #endif
01942 
01943   if (nbonds > 0 && fromptr != NULL && toptr != NULL) {
01944     js->optflags |= JSOPT_BONDS; 
01945 
01946     /* save bond info until we actually write out the structure file */
01947     js->nbonds = nbonds;
01948     js->bondfrom = (int *) malloc(nbonds * sizeof(int));
01949     memcpy(js->bondfrom, fromptr, nbonds * sizeof(int));
01950     js->bondto = (int *) malloc(nbonds * sizeof(int));
01951     memcpy(js->bondto, toptr, nbonds * sizeof(int));
01952 
01953     if (bondorder != NULL) {
01954       js->optflags |= JSOPT_BONDORDERS;
01955       js->bondorders = (float *) malloc(nbonds * sizeof(float));
01956       memcpy(js->bondorders, bondorder, nbonds * sizeof(float));
01957     }
01958   }
01959 
01960   return MOLFILE_SUCCESS;
01961 }
01962 
01963 #if vmdplugin_ABIVERSION > 14
01964 static int write_js_angles(void * v, int numangles, const int *angles,
01965                            const int *angletypes, int numangletypes,
01966                            const char **angletypenames, int numdihedrals, 
01967                            const int *dihedrals, const int *dihedraltype,
01968                            int numdihedraltypes, const char **dihedraltypenames,
01969                            int numimpropers, const int *impropers, 
01970                            const int *impropertypes, int numimpropertypes, 
01971                            const char **impropertypenames, int numcterms, 
01972                            const int *cterms, int ctermcols, int ctermrows) {
01973   jshandle *js = (jshandle *) v;
01974 
01975   /* save info until we actually write out the structure file */
01976   js->numangles = numangles;
01977   js->numdihedrals = numdihedrals;
01978   js->numimpropers = numimpropers;
01979   js->numcterms = numcterms;
01980 
01981 #if defined(INFOMSGS)
01982   if (js->verbose) {
01983     printf("jsplugin) write_js_angles():\n");
01984     printf("jsplugin) storing angles/dihedrals/impropers for writing...\n");
01985     printf("jsplugin) %d %d %d %d\n",
01986            numangles, numdihedrals, numimpropers, numcterms);
01987   }
01988 #endif
01989 
01990   if (js->numangles > 0 || js->numdihedrals > 0 || js->numimpropers > 0) {
01991     js->optflags |= JSOPT_ANGLES;
01992 
01993     js->angles = (int *) malloc(3L*js->numangles*sizeof(int));
01994     memcpy(js->angles, angles, 3L*js->numangles*sizeof(int));
01995     js->dihedrals = (int *) malloc(4L*js->numdihedrals*sizeof(int));
01996     memcpy(js->dihedrals, dihedrals, 4L*js->numdihedrals*sizeof(int));
01997     js->impropers = (int *) malloc(4L*js->numimpropers*sizeof(int));
01998     memcpy(js->impropers, impropers, 4L*js->numimpropers*sizeof(int));
01999   }
02000   if (js->numcterms > 0) {
02001     js->optflags |= JSOPT_CTERMS;
02002 
02003     js->cterms = (int *) malloc(8L*js->numcterms*sizeof(int));
02004     memcpy(js->cterms, cterms, 8L*js->numcterms*sizeof(int));
02005   }
02006 
02007   return MOLFILE_SUCCESS;
02008 }
02009 #else
02010 static int write_js_angles(void * v,
02011         int numangles,    const int *angles,    const double *angleforces,
02012         int numdihedrals, const int *dihedrals, const double *dihedralforces,
02013         int numimpropers, const int *impropers, const double *improperforces,
02014         int numcterms,   const int *cterms,
02015         int ctermcols, int ctermrows, const double *ctermforces) {
02016   jshandle *js = (jshandle *) v;
02017 
02018   /* save info until we actually write out the structure file */
02019   js->numangles = numangles;
02020   js->numdihedrals = numdihedrals;
02021   js->numimpropers = numimpropers;
02022   js->numcterms = numcterms;
02023 
02024   if (js->numangles > 0 || js->numdihedrals > 0 || js->numimpropers > 0) {
02025     js->optflags |= JSOPT_ANGLES;
02026 
02027     js->angles = (int *) malloc(3L*js->numangles*sizeof(int));
02028     memcpy(js->angles, angles, 3L*js->numangles*sizeof(int));
02029     js->dihedrals = (int *) malloc(4L*js->numdihedrals*sizeof(int));
02030     memcpy(js->dihedrals, dihedrals, 4L*js->numdihedrals*sizeof(int));
02031     js->impropers = (int *) malloc(4L*js->numimpropers*sizeof(int));
02032     memcpy(js->impropers, impropers, 4L*js->numimpropers*sizeof(int));
02033   }
02034   if (js->numcterms > 0) {
02035     js->optflags |= JSOPT_CTERMS;
02036 
02037     js->cterms = (int *) malloc(8L*js->numcterms*sizeof(int));
02038     memcpy(js->cterms, cterms, 8L*js->numcterms*sizeof(int));
02039   }
02040 
02041   return MOLFILE_SUCCESS;
02042 }
02043 #endif
02044 #endif
02045 
02046 
02047 static int write_js_timestep(void *v, const molfile_timestep_t *ts) { 
02048   jshandle *js = (jshandle *)v;
02049   double *unitcell=NULL;
02050   ptrdiff_t zeropadsz=0;
02051 
02052   /* If no structure data was written and this is the first timestep */
02053   /* we must complete writing the file header and performing the     */
02054   /* seek to the next filesystem block and VM-page boundary when     */
02055   /* using direct I/O APIs...                                        */
02056   if (js->directio_ucell_blkbuf == NULL) {
02057     printf("jsplugin) no structure data, writing timesteps only...\n");
02058 
02059     /* use block-based I/O by default when writing structures larger */
02060     /* than JSBLOCKIO_THRESH atoms, or when directed by the user     */
02061     js_blockio_check_and_set(js);
02062 
02063     /* write flags data to the file */
02064     fio_write_int32(js->fd, js->optflags); 
02065     printf("jsplugin) writing option flags: %0x08x\n", js->optflags);
02066 
02067     /* Check to see if block-based trajectory I/O is used  */
02068     /* and write out the block size for this file.         */
02069     if (js->optflags & JSOPT_TS_BLOCKIO) {
02070       fio_fwrite(&js->directio_block_size, sizeof(int), 1, js->fd);
02071       printf("jsplugin) Block-based I/O enabled: block size %d bytes\n", 
02072              js->directio_block_size);
02073     }
02074 
02075     /* update the file offset for the first timestep */
02076     js_calc_timestep_blocking_info(js);
02077   }
02078 
02079   /* set unit cell pointer to the TS block-aligned buffer area */
02080   unitcell = (double *) js->directio_ucell_blkbuf;
02081 
02082   js->nframes++; /* increment frame count written to the file so far */
02083 
02084   unitcell[0] = ts->A;
02085   unitcell[1] = ts->B;
02086   unitcell[2] = ts->C;
02087   unitcell[3] = sin((M_PI_2 / 90.0) * (90.0 - ts->alpha));
02088   unitcell[4] = sin((M_PI_2 / 90.0) * (90.0 - ts->beta));
02089   unitcell[5] = sin((M_PI_2 / 90.0) * (90.0 - ts->gamma));
02090 
02091   /* coordinates for all atoms */
02092   if (fio_fwrite(ts->coords, js->ts_crd_sz, 1, js->fd) != 1) {
02093     printf("jsplugin) Error writing timestep coords!\n");
02094     return MOLFILE_ERROR;
02095   }
02096 
02097   /* correctly handle block-based direct-I/O output format       */
02098   /* write out coord padding bytes using zero buffer in jshandle */
02099   zeropadsz = js->ts_crd_padsz - js->ts_crd_sz;
02100   if (zeropadsz > 0) {
02101     if ((zeropadsz > MOLFILE_DIRECTIO_MAX_BLOCK_SIZE) ||
02102         (fio_fwrite(js->blockpad, zeropadsz, 1, js->fd) != 1)) {
02103       printf("jsplugin) Error writing timestep coord padding!\n");
02104       return MOLFILE_ERROR;
02105     } 
02106   }
02107 
02108   /* PBC unit cell info */ 
02109   if (fio_fwrite(unitcell, js->ts_ucell_sz, 1, js->fd) != 1) {
02110     printf("jsplugin) Error writing timestep unit cell!\n");
02111     return MOLFILE_ERROR;
02112   }
02113 
02114   /* correctly handle block-based direct-I/O output format       */
02115   /* write out PBC padding bytes using zero buffer in jshandle */
02116   zeropadsz = js->ts_ucell_padsz - js->ts_ucell_sz;
02117   if (zeropadsz > 0) {
02118     if ((zeropadsz > MOLFILE_DIRECTIO_MAX_BLOCK_SIZE) ||
02119         (fio_fwrite(js->blockpad, zeropadsz, 1, js->fd) != 1)) {
02120       printf("jsplugin) Error writing timestep PBC padding!\n");
02121       return MOLFILE_ERROR;
02122     } 
02123   }
02124 
02125   return MOLFILE_SUCCESS;
02126 }
02127 
02128 
02129 static void close_js_write(void *v) {
02130   jshandle *js = (jshandle *)v;
02131 
02132   /* update the trajectory header information */
02133   fio_fseek(js->fd, JSNFRAMESOFFSET, FIO_SEEK_SET);
02134   fio_write_int32(js->fd, js->nframes);
02135   fio_fseek(js->fd, 0, FIO_SEEK_END);
02136 
02137   fio_fclose(js->fd);
02138 
02139 #if JSMAJORVERSION > 1
02140   if (js->directio_ucell_ptr)
02141     free(js->directio_ucell_ptr);
02142 
02143   if (js->bondfrom)
02144     free(js->bondfrom);
02145   if (js->bondto)
02146     free(js->bondto);
02147   if (js->bondorders)
02148     free(js->bondorders);
02149 
02150   if (js->angles)
02151     free(js->angles);
02152   if (js->dihedrals)
02153     free(js->dihedrals);
02154   if (js->impropers)
02155     free(js->impropers);
02156   if (js->cterms)
02157     free(js->cterms);
02158 #endif
02159 
02160   free(js);
02161 }
02162 
02163 
02164 /*
02165  * Initialization stuff here
02166  */
02167 static molfile_plugin_t plugin;
02168 
02169 #if !defined(VMDJSPLUGININCLUDESRC)
02170 
02171 VMDPLUGIN_API int VMDPLUGIN_init() {
02172   memset(&plugin, 0, sizeof(molfile_plugin_t));
02173   plugin.abiversion = vmdplugin_ABIVERSION;
02174   plugin.type = MOLFILE_PLUGIN_TYPE;
02175   plugin.name = "js";
02176   plugin.prettyname = "js";
02177   plugin.author = "John Stone";
02178   plugin.majorv = JSMAJORVERSION;
02179   plugin.minorv = JSMINORVERSION;
02180   plugin.is_reentrant = VMDPLUGIN_THREADSAFE;
02181   plugin.filename_extension = "js";
02182   plugin.open_file_read = open_js_read;
02183 #if JSMAJORVERSION > 1
02184   plugin.read_structure = read_js_structure;
02185   plugin.read_bonds = read_js_bonds;
02186   plugin.read_angles = read_js_angles;
02187 #endif
02188   plugin.read_next_timestep = read_js_timestep;
02189   plugin.close_file_read = close_js_read;
02190   plugin.open_file_write = open_js_write;
02191 #if JSMAJORVERSION > 1
02192   plugin.write_structure = write_js_structure;
02193   plugin.write_bonds = write_js_bonds;
02194   plugin.write_angles = write_js_angles;
02195 #endif
02196   plugin.write_timestep = write_js_timestep;
02197   plugin.close_file_write = close_js_write;
02198 #if vmdplugin_ABIVERSION > 17
02199   plugin.read_timestep_pagealign_size = read_js_timestep_pagealign_size;
02200 #endif
02201   return VMDPLUGIN_SUCCESS;
02202 }
02203 
02204 VMDPLUGIN_API int VMDPLUGIN_register(void *v, vmdplugin_register_cb cb) {
02205   (*cb)(v, (vmdplugin_t *)&plugin);
02206   return VMDPLUGIN_SUCCESS;
02207 }
02208 
02209 VMDPLUGIN_API int VMDPLUGIN_fini() {
02210   return VMDPLUGIN_SUCCESS;
02211 }
02212 
02213 #endif
02214   
02215 #ifdef TEST_JSPLUGIN
02216 
02217 #include <sys/time.h>
02218 
02219 #if defined(ENABLECUDATESTS)
02220 #include <cuda_runtime.h>
02221 
02222 #if defined(ENABLECUDAGDS)
02223 #include <cufile.h>
02224 #endif
02225 #endif
02226 
02227 /* get the time of day from the system clock, and store it (in seconds) */
02228 double time_of_day(void) {
02229 #if defined(_MSC_VER)
02230   double t;
02231 
02232   t = GetTickCount();
02233   t = t / 1000.0;
02234 
02235   return t;
02236 #else
02237   struct timeval tm;
02238   struct timezone tz;
02239 
02240   gettimeofday(&tm, &tz);
02241   return((double)(tm.tv_sec) + (double)(tm.tv_usec)/1000000.0);
02242 #endif
02243 }
02244 
02245 int main(int argc, char *argv[]) {
02246   molfile_timestep_t timestep;
02247   float *coords0=NULL, *aligncoords0=NULL;
02248   float *coords1=NULL, *aligncoords1=NULL;
02249 
02250   void *v;
02251   jshandle *js;
02252   int natoms, i;
02253   ptrdiff_t sz, blocksz;
02254   float sizeMB =0.0, totalMB = 0.0;
02255   double starttime, endtime, totaltime = 0.0;
02256   int do_io = 1;
02257   int verbose = 0;
02258   int overlapiogpu = 1;
02259 
02260   printf("Standalone tests for JS plugin:\n");
02261   
02262   if (getenv("VMDJSNOIO") != NULL)
02263     do_io = 0;
02264 
02265   if (getenv("VMDJSVERBOSE") != NULL || getenv("VMDJSTESTVERBOSE"))
02266     verbose = 1;    
02267 
02268   if (do_io)
02269     printf("  Timestep disk I/O enabled.\n");
02270   else
02271     printf("  Timestep disk I/O DISABLED.\n");
02272 
02273 #if defined(ENABLECUDATESTS)
02274   printf("  CUDA GPU support compiled in.\n");
02275 
02276   // If the code is compiled with CUDA support, we benchmark 
02277   // host I/O immediately followed by host-GPU copies of each timestep
02278   cudaError_t crc;
02279   cudaStream_t devstream;
02280   ptrdiff_t maxatomidx=-1;
02281   int devcount;
02282   float *devptr=NULL;
02283 
02284   crc = cudaGetDeviceCount(&devcount);
02285   printf("  GPU device count: %d\n", devcount);
02286   if (devcount==0)
02287     printf("  No GPU devices, continuing with host only...\n");
02288 
02289   // Only do the CUDA tests if asked to
02290   if (getenv("VMDJSCUDATESTS") == NULL) {
02291     devcount = 0;
02292     printf("  GPU tests disabled.\n");
02293     printf("  Enable GPU tests with VMDJSCUDATESTS env variable\n");
02294   } else {
02295     printf("  Disable GPU tests by unsetting VMDJSCUDATESTS env variable\n");
02296   }
02297 
02298 #if defined(ENABLEJSSHORTREADS)
02299   /* test code for an implementation that does short reads that */
02300   /* skip bulk solvent, useful for faster loading of very large */
02301   /* structures                                                 */
02302   if (getenv("VMDJSMAXATOMIDX") != NULL) {
02303     fio_size_t bszmask;
02304 
02305     maxatomidx = atoi(getenv("VMDJSMAXATOMIDX"));
02306     if (maxatomidx < 0)
02307       maxatomidx = 0;
02308     if (maxatomidx >= js->natoms)
02309       maxatomidx = js->natoms - 1;
02310 
02311     printf("jsplugin) Short-copies of GPU timesteps enabled: %ld / %ld atoms (%.2f%%)\n",
02312            maxatomidx, js->natoms, 100.0*(maxatomidx+1) / ((float) js->natoms));
02313   }
02314 #endif
02315 #endif
02316 
02317 
02318   while (--argc) {
02319     int syncframe;
02320     ++argv; 
02321     natoms = 0;
02322     v = open_js_read(*argv, "js", &natoms);
02323     if (!v) {
02324       printf("jsplugin) open_js_read failed for file %s\n", *argv);
02325       return 1;
02326     }
02327     js = (jshandle *)v;
02328     sizeMB = ((natoms * 3.0) * js->nframes * 4.0) / (1024.0 * 1024.0);
02329     totalMB += sizeMB; 
02330     printf("jsplugin) file: %s\n", *argv);
02331     printf("jsplugin)   %d atoms, %d frames, size: %6.1fMB\n", natoms, js->nframes, sizeMB);
02332 
02333     starttime = time_of_day();
02334 
02335     /* ensure we have a large enough allocation so we can align */
02336     /* the starting pointer to a blocksz page boundary          */
02337     blocksz = MOLFILE_DIRECTIO_MIN_BLOCK_SIZE;
02338     sz = 3L*sizeof(float)*natoms + blocksz;
02339 
02340     /* pad the allocation to an even multiple of the block size */
02341     size_t blockpadsz = (sz + (blocksz - 1)) & (~(blocksz - 1));
02342 
02343     /* allocate multi buffers so we can have concurrent work/transfers/reads */
02344     aligncoords0 = (float *) alloc_aligned_ptr(sz, blocksz, (void**) &coords0);
02345     aligncoords1 = (float *) alloc_aligned_ptr(sz, blocksz, (void**) &coords1);
02346 
02347 #if vmdplugin_ABIVERSION > 17
02348     int filepgalignsz = 1;
02349     read_js_timestep_pagealign_size(v, &filepgalignsz);
02350     if (filepgalignsz != blocksz) {
02351       printf("jsplugin) Plugin-returned page alignment size doesn't match!\n");
02352     } else {
02353       printf("jsplugin) Page alignment size: %d\n", filepgalignsz);
02354     }
02355 #endif
02356 
02357 #if defined(ENABLECUDATESTS)
02358 #if defined(ENABLECUDAGDS)
02359     int cufileinuse=0;
02360     CUfileHandle_t cfh;
02361     CUfileDescr_t cfhdesc;
02362     CUfileError_t cferr;
02363     memset(&cfh, 0, sizeof(cfh));
02364     memset(&cfhdesc, 0, sizeof(cfhdesc));
02365 #endif
02366 
02367     if (crc == cudaSuccess && devcount > 0) {
02368       cudaSetDevice(0);
02369       printf("jsplugin) allocating GPU memory buffer for CUDA tests...\n");
02370       crc = cudaMalloc((void**) &devptr, blockpadsz);
02371       if (crc != cudaSuccess) {
02372         printf("Failed to allocate GPU buffer!\n");
02373         return -1;
02374       }
02375 
02376       cudaStreamCreate(&devstream);
02377 
02378 #if defined(ENABLECUDAGDS)
02379       cuFileBufRegister(devptr, blockpadsz, 0); 
02380 
02381       cfhdesc.handle.fd = js->directio_fd; // typedef of Unix FD
02382       cfhdesc.type = CU_FILE_EXTERNAL_MEMORY_HANDLE_TYPE_OPAQUE_FD;
02383       cferr = cuFileImportExternalFile(&cfh, &cfhdesc);
02384       if (cferr.err != CU_FILE_SUCCESS) {
02385         printf("Failed to import file handle for use by cuFile APIs!\n");
02386         return -1;
02387       }
02388       cufileinuse=1;
02389 #endif
02390     }
02391 #endif
02392 
02393     /* loop over all timesteps ... */
02394     for (syncframe=0,i=0; i<js->nframes; i++) {
02395       if (do_io) {
02396         /* read even/odd frame into alternating buffers so     */
02397         /* that we can overlap a read with an ongoing GPU copy */
02398         if (i & 1) 
02399           timestep.coords = aligncoords1;
02400         else  
02401           timestep.coords = aligncoords0;
02402 
02403         /* disk I/O is synchronous */
02404         if (verbose) {
02405           printf("%sreading frame[%d]...", (i!=0) ? "\r" : "", i);
02406           fflush(stdout);
02407         }
02408         int rc=0;
02409 #if defined(ENABLECUDAGDS)
02410         if (cufileinuse) {
02411 #if 0
02412           /* read an even multiple of the block size */
02413           ptrdiff_t rsz = natoms * 3L * sizeof(float);
02414           rsz = (rsz + (blocksz - 1)) & (~(blocksz - 1));
02415           ptrdiff_t foffset = i * rsz;
02416 #endif
02417           ptrdiff_t startoffset, foffset, readlen;
02418           fio_fd directio_fd;
02419           read_js_timestep_index_offsets(v, natoms, i, 0, natoms,
02420                                          &directio_fd,
02421                                          &startoffset,
02422                                          &foffset,
02423                                          &readlen);
02424 
02425 printf("cuFileRead(): offset %ld  readlen: %ld\n", foffset, readlen);
02426           ptrdiff_t ret = 0;
02427           ret = cuFileRead(cfh, (char *) devptr, readlen, foffset);
02428           if (ret < 0) {
02429             const char *descp = "unknown error code";
02430 #if 0
02431             // XXX this requires linkage with libcuda.so, which is
02432             //     against convention, so we avoid it for now
02433             if (cuGetErrorName(ret, &descp) != CUDA_SUCCESS)
02434               descp = "unknown cuda error";
02435 #endif
02436 
02437             printf("Error: cuFileRead(): %ld, '%s'\n", ret, descp);
02438             return -1;
02439           }
02440         } else
02441 #else
02442           rc = read_js_timestep(v, natoms, &timestep);
02443 #endif
02444 
02445         if (rc) {
02446           printf("jsplugin) error in read_js_timestep on frame %d\n", i);
02447           /* return 1; */
02448         }
02449       }
02450 
02451 #if defined(ENABLECUDATESTS)
02452       if (crc == cudaSuccess && devcount > 0) {
02453 #if defined(ENABLECUDAGDS)
02454         if (!cufileinuse) {
02455 #endif
02456           /* allow overlap of async memcpy with next round of direct I/O */
02457           if (overlapiogpu) {
02458             if (verbose) {
02459               printf("sync frame[%d]...", syncframe);
02460               fflush(stdout);
02461             }
02462             cudaStreamSynchronize(devstream);
02463           }
02464 
02465           if (verbose) {
02466             printf("cudaMemcpyAsync() frame[%d]...", i);
02467             fflush(stdout);
02468           }
02469 
02470           size_t bsz = (maxatomidx >= 0) ? (maxatomidx+1) : natoms;
02471           bsz *= 3L*sizeof(float);
02472           crc = cudaMemcpyAsync(devptr, timestep.coords, bsz, cudaMemcpyHostToDevice, devstream);
02473           syncframe=i;
02474 
02475           if (!overlapiogpu) {
02476             if (verbose) {
02477               printf("sync frame[%d]...", syncframe);
02478               fflush(stdout);
02479             }
02480             cudaStreamSynchronize(devstream);
02481           }
02482   
02483           if (verbose) {
02484             printf("       ");
02485             fflush(stdout);
02486           }
02487 #if defined(ENABLECUDAGDS)
02488         }
02489 #endif
02490       }
02491 #endif
02492     }
02493 
02494 #if defined(ENABLECUDATESTS)
02495     /* wait for last GPU memcpy */
02496     cudaStreamSynchronize(devstream);
02497     printf("\n");
02498 
02499     /* wait for any pending GPU calls to complete */
02500     cudaDeviceSynchronize();
02501 
02502 #if defined(ENABLECUDAGDS)
02503     if (cufileinuse) {
02504       cuFileBufDeregister(devptr); 
02505       cuFileDestroyFile(&cfh);
02506     }
02507 #endif
02508 #endif
02509 
02510     endtime = time_of_day();
02511     close_js_read(v);
02512 
02513 #if defined(ENABLECUDATESTS)
02514     cudaStreamDestroy(devstream);
02515     if (crc == cudaSuccess && devcount > 0) {
02516       cudaFree(devptr);
02517     }
02518 #endif
02519     free(coords0);
02520     free(coords1);
02521 
02522     totaltime += endtime - starttime;
02523     printf("jsplugin)  Time: %5.1f seconds\n", endtime - starttime);
02524     printf("jsplugin)  Speed: %5.1f MB/sec, %5.1f timesteps/sec\n", sizeMB / (endtime - starttime), (js->nframes / (endtime - starttime)));
02525   }
02526   printf("jsplugin) Overall Size: %6.1f MB\n", totalMB);
02527   printf("jsplugin) Overall Time: %6.1f seconds\n", totaltime);
02528   printf("jsplugin) Overall Speed: %5.1f MB/sec\n", totalMB / totaltime);
02529   return 0;
02530 }
02531       
02532 #endif
02533 

Generated on Wed Nov 11 03:06:27 2020 for VMD Plugins (current) by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002