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

graspplugin.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: graspplugin.C,v $
00013  *      $Author: johns $       $:  $             $State: Exp $
00014  *      $Revision: 1.23 $       $Date: 2016/11/28 05:01:54 $
00015  *
00016  ***************************************************************************/
00017 
00018 /* 
00019  * Reader for GRASP binary surface files
00020  * The file format is briefly described at these two sites:
00021  *   http://honiglab.cpmc.columbia.edu/grasp/grasp_contents.html#A.1
00022  *   http://www.msg.ucsf.edu/local/programs/grasp/html/Appendix%20A.html
00023  */
00024 
00025 //Modified by Biol. Angel H. Jiménez Pardo and Luis Rosales Leòn
00026 //Visualisation Department
00027 //DGSCA
00028 //Universidad Nacional Autònoma de Mèxico UNAM
00029 //2006
00030 
00031 #include <stdio.h>
00032 #include <math.h>
00033 #include <string.h>
00034 #include "molfile_plugin.h"
00035 #include "endianswap.h"
00036 
00038 
00039 #define POTENTIALS    1
00040 #define CURVATURES    2
00041 #define DISTANCES     4
00042 #define PROPERTY_1    8
00043 #define PROPERTY_2   16
00044 #define VERTEXCOLOR  32
00045 #define GLOBALCOLOR  64
00046 
00048 typedef float COLOUR [3]; 
00049 typedef float COLOUR [3];
00050 typedef unsigned int SHORTWORD;
00051 typedef float VECTOR [3];
00052 
00054 typedef struct {
00055   SHORTWORD flag, nvert;
00056   COLOUR clo, cmd, chi, colors;
00057   VECTOR range;
00058 } GRASSP;
00059 
00061 typedef struct {
00062   FILE *fd;
00063   molfile_graphics_t *graphics;
00064 } grasp_t;
00065 
00066 
00068 void Set_Colour(float *c, float r, float g, float b) {
00069   *(c+0)=r;
00070   *(c+1)=g;
00071   *(c+2)=b;
00072 } //Set_Colour end
00073 
00074 
00076 void ClinComb2(float c[], float k1, float c1[], float k2, float c2[]) {
00077   (c)[0]=(k1)*(c1)[0] + (k2)*(c2)[0];
00078   (c)[1]=(k1)*(c1)[1] + (k2)*(c2)[1];
00079   (c)[2]=(k1)*(c1)[2] + (k2)*(c2)[2];
00080 } //ClinComb2 end
00081 
00082 
00084 void Get_Property_Values(GRASSP *grassp, float *properties, float *colores, int nvert) {
00085   long  i;
00086   const char *name[] = {"potential","curvature","distance","property1","property2"};
00087   SHORTWORD  index;
00088   int k=0, j=0;
00089   float  val, weight, min, mid, max, midmin, maxmid;
00090   val = weight = min = mid = max = midmin = maxmid = 0.0;
00091 
00092   //Range values
00093   grassp->range[0]=-1000.0;
00094   grassp->range[1]= 0.0;
00095   grassp->range[2]= 1000.0;
00096 
00098   index = (SHORTWORD)( (log((double) grassp->flag) / log(2.0)) + 0.5 );
00099 
00101   if ((grassp->flag)!=POTENTIALS) {
00102     if (index <= 4)
00103       printf("graspplugin) No data available for '%s' option\n", name[index]);
00104     else 
00105       printf("graspplugin) out of range property, flag: %d index: %d\n", grassp->flag, index);
00106     printf("graspplugin) Will use white color instead\n");
00107     grassp->flag = GLOBALCOLOR;
00108     Set_Colour(grassp->clo, 1, 1, 1);
00109   } else {
00110     printf("graspplugin) Getting %s values.\n", name[index]);
00111   }
00112 
00114   max=  0.01;  /* should be > 0 */
00115   min= -0.01;  /* should be < 0 */
00116 
00117   // Get values
00118   for(i=0; i < nvert; i++) {
00119     if (properties[i] < min)
00120       min= properties[i];
00121     else if(properties[i] > max)
00122       max= properties[i];
00123   }
00124 
00125   // cut properties that are out of specified range
00126   if (min < grassp->range[0] || max > grassp->range[2]) {
00127     for(i=0; i < nvert; i++) {
00128       if(properties[i] < grassp->range[0])
00129         properties[i]= grassp->range[0];
00130       else if(properties[i] > grassp->range[2])
00131         properties[i]= grassp->range[2];
00132     }
00133   } else {
00134     // or reset range
00135     grassp->range[0] = min;
00136     grassp->range[2] = max;
00137   }
00138 
00139   // check mid value 
00140   if (grassp->range[1] <= grassp->range[0] || 
00141       grassp->range[1] >= grassp->range[2])
00142    grassp->range[1] = (grassp->range[0] + grassp->range[2]) / 2;
00143 
00144   printf("graspplugin) Computing colors for range %g,%g,%g\n", 
00145          grassp->range[0], grassp->range[1], grassp->range[2]);
00146 
00147   // Prepare color interpolation parameters
00148   min = grassp->range[0];
00149   mid = grassp->range[1];
00150   max = grassp->range[2];
00151   midmin = mid-min;
00152   maxmid = max-mid;
00153 
00154   // Create color for each vertex and copies to a vector
00155   k=0;
00156   for (i=0; i < nvert; i++) {
00157     val = properties[i];
00158     if (val <= mid) {
00159       weight = (midmin) ? (val-min)/midmin : 0;  
00160       ClinComb2(grassp->colors, 1-weight, grassp->clo, weight, grassp->cmd);
00161       for(j=0; j<=2; j++) {
00162         *(colores+k)=grassp->colors[j];
00163         k++;
00164       }
00165     } else {
00166       weight = (maxmid) ? (val-mid)/maxmid : 0;
00167       ClinComb2(grassp->colors, 1-weight, grassp->cmd, weight, grassp->chi);
00168       for (j=0; j<=2; j++) {
00169         *(colores+k)=grassp->colors[j];
00170         k++;
00171       }       
00172     }
00173     // clean up 
00174   }
00175 } // Get_Property_Values end
00176 
00177 
00179 void line3 (FILE * infile, GRASSP * grassp) {
00180   char line3 [81];
00181   fread(line3, 1, 80, infile);
00182 
00183   grassp->flag=0;
00184   int i=0;
00185   if (line3 [0]==',') i++;
00186   while (i<80 && line3[i]!=' ') {
00187 #if 0
00188     // XXX the rest of the code doesn't properly process flags yet, so
00189     // there's no point in setting them currently.
00190     if(!strncmp(line3+i, "potentials", 10))   grassp->flag |= POTENTIALS;
00191     if(!strncmp(line3+i, "curvature",   9))   grassp->flag |= CURVATURES;
00192     if(!strncmp(line3+i, "distances",   9))   grassp->flag |= DISTANCES;
00193     if(!strncmp(line3+i, "gproperty",   9))   grassp->flag |= PROPERTY_1;
00194     if(!strncmp(line3+i, "g2property", 10))   grassp->flag |= PROPERTY_2;
00195     if(!strncmp(line3+i, "vertexcolor",11))   grassp->flag |= VERTEXCOLOR;
00196 #endif
00197     i++;
00198   }
00199 
00200   // Assign default property colors
00201   if (grassp->flag > 0 && grassp->flag < VERTEXCOLOR) {
00202     switch (grassp->flag ) {
00203       case POTENTIALS:
00204         Set_Colour(grassp->clo, 1.0, 0.0, 0.0 );
00205         Set_Colour(grassp->cmd, 1.0, 1.0, 1.0 );
00206         Set_Colour(grassp->chi, 0.0, 0.0, 1.0 );
00207         break;
00208 
00209       case CURVATURES:
00210         Set_Colour(grassp->clo, 0.5, 0.5, 0.5 );
00211         Set_Colour(grassp->cmd, 1.0, 1.0, 1.0 );
00212         Set_Colour(grassp->chi, 0.0, 1.0, 0.0 );
00213         break;
00214 
00215      case DISTANCES:
00216        Set_Colour(grassp->clo, 1.0, 1.0, 1.0 );
00217        Set_Colour(grassp->cmd, 0.0, 0.0, 1.0 );
00218        Set_Colour(grassp->chi, 1.0, 0.0, 0.0 );
00219        break;
00220 
00221      default: // Global color
00222        Set_Colour(grassp->clo, 1.0, 0.0, 0.0 );
00223        Set_Colour(grassp->cmd, 0.5, 0.0, 0.5 );
00224        Set_Colour(grassp->chi, 0.0, 0.0, 1.0 );
00225        break;
00226    }
00227  }
00228 
00229  if (!grassp->flag)
00230    grassp->flag = GLOBALCOLOR; 
00231 } //line3 end
00232 
00233 
00234 // check endianness
00235 static int is_little_endian(void) {
00236   int x=1;
00237   return *((char *)&x);
00238 }   
00239 
00240 
00241 static void *open_file_read(const char *filepath, const char *filetype,
00242     int *natoms) {
00243   FILE *fd;
00244   grasp_t *grasp;
00245   
00246   fd = fopen(filepath, "rb");
00247   if (!fd) 
00248     return NULL;
00249   grasp = new grasp_t;
00250   grasp->fd = fd;
00251   grasp->graphics = NULL;
00252   *natoms = 0;
00253   return grasp;
00254 }
00255 
00256 
00257 static int read_rawgraphics(void *v, int *nelem, 
00258     const molfile_graphics_t **data) {
00259   grasp_t *grasp = (grasp_t *)v;
00260   FILE *infile = grasp->fd;
00261 
00262   // Reverse engineering is your friend, and combined with FORTRAN code, voila!
00263   // od -c shows the header starts off:
00264   // \0  \0  \0   P   f   o   r   m   a   t   =   2
00265   // and according to ungrasp, this is a 1 for grasp versions 1.0
00266   // and 1.1, and 2 for grasp version 1.2
00267   // Also, the header lines are of length 80 characters + 4 header chars
00268   // + 4 trailer chars
00269   // The 4 bytes at the beginning/end are standard Fortran array trash
00270 
00272   GRASSP datax;
00273    
00274   char trash[4];
00275 #define TRASH fread(trash, 4, 1, infile)
00276   char line[81];
00277 
00278   // FIRST LINE OF HEADER; contains format type
00279   TRASH; 
00280   fread(line, 1, 80, infile); 
00281   // make sure it says 'format='
00282   if (strncmp(line, "format=", 7) != 0) {
00283     printf("graspplugin) First characters of file don't look like a GRASP file\n");
00284     return MOLFILE_ERROR;
00285   }
00286   TRASH;
00287 
00288   // next char should be a 0 or 1
00289   char gfiletype = line[7];
00290   if (gfiletype == '1') {
00291     gfiletype = 1;
00292   } else if (gfiletype == '2') {
00293     gfiletype = 2;
00294   } else {
00295     printf("graspplugin) GRASP file is in format %c, but only '1' or '2' is supported\n", gfiletype);
00296     return MOLFILE_ERROR;
00297   }
00298 
00299   // SECOND LINE: contains "vertices,accessibles,normals,triangles"
00300   TRASH; 
00301   fread(line, 1, 80, infile); 
00302   TRASH;
00303 
00304   // THIRD LINE: contains (0 or more of)?
00305   //  "potentials,curvature,distances,gproperty,g2property,vertexcolor
00306   TRASH; 
00307   line3(infile, &datax);
00308   TRASH;
00309 
00310   // FOURTH LINE stores vertices, triangles, gridsize, lattice spacing
00311   int nvert, ntriangles, gridsize;
00312   float lattice;
00313   TRASH; 
00314   fread(line, 1, 80, infile); 
00315   TRASH;
00316   sscanf(line, "%d%d%d%f", &nvert, &ntriangles, &gridsize, &lattice);
00317 
00319   float *colores = new float [3*nvert];
00320 
00321   // FIFTH LINE stores the center (x,y,z) position
00322   float center[3];
00323   TRASH; 
00324   fread(line, 1, 80, infile); 
00325   TRASH;
00326   sscanf(line, "%f%f%f", center, center+1, center+2);
00327 
00328   float *vertex = new float[3 * nvert];
00329   float *access = new float[3 * nvert];
00330   float *normal = new float[3 * nvert];
00331   int *triangle = new int[3 * ntriangles];
00332   float *properties = new float[3* nvert];
00333 
00334   if (!vertex || !access || !normal || !triangle || !properties) {
00335     delete [] vertex;
00336     delete [] access;
00337     delete [] normal;
00338     delete [] triangle;
00339     delete [] properties;
00340     printf("graspplugin) Failed vertex/access/normal/triangle allocations.\n");
00341     return MOLFILE_ERROR;
00342   }
00343 
00344   // ungrasp says:
00345   //    if (filetype.eq.1) then integer*2
00346   //    if (filetype.eq.2) then integer*4
00347 
00348   // And read them in.  Who needs error checking?
00349   TRASH; 
00350   fread(vertex, 3 * sizeof(float), nvert, infile); 
00351   TRASH;
00352   TRASH; 
00353   fread(access, 3 * sizeof(float), nvert, infile); 
00354   TRASH;
00355   TRASH; 
00356   fread(normal, 3 * sizeof(float), nvert, infile); 
00357   TRASH;
00358  
00359   if (is_little_endian()) {
00360     swap4_aligned(vertex, 3*nvert);
00361     swap4_aligned(access, 3*nvert);
00362     swap4_aligned(normal, 3*nvert);
00363   }
00364 
00365   if (gfiletype == 2) {
00366     TRASH; 
00367     fread(triangle, 3 * sizeof(int), ntriangles, infile); 
00368     TRASH;
00369     TRASH;
00370     fread(properties, 3 * sizeof(float), nvert, infile);
00371     if (is_little_endian()) {
00372       swap4_aligned(triangle, 3*ntriangles);
00373       swap4_aligned(properties, 3*nvert);
00374     }
00375   } else {
00376 #if 1
00377     int i;
00378     short *striangle = new short[3 * ntriangles];
00379     if (!striangle) {
00380       delete [] vertex;
00381       delete [] access;
00382       delete [] normal;
00383       delete [] triangle;
00384       delete [] properties;
00385       printf("graspplugin) Failed short triangle allocation.\n");
00386       return MOLFILE_ERROR;
00387     }
00388 
00389     TRASH;
00390     fread(striangle, sizeof(short), 3 * ntriangles, infile);
00391     TRASH;
00392     TRASH;
00393     fread(properties, sizeof(float), 3 * nvert, infile);
00394     
00395     if (is_little_endian()) {
00396     swap2_aligned(striangle, 3 * ntriangles);
00397     swap4_aligned(properties, 3*nvert);}
00398     
00399     for (i=0; i<3*ntriangles; i++) {
00400       triangle[i] = striangle[i];
00401     }
00402     delete [] striangle;  
00403     
00404 #else
00405     // do it the slow way (converting from short to int)
00406     int i;
00407     short tmp[3];
00408     TRASH;
00409     for (i=0; i<ntriangles; i++) {
00410       fread(tmp, sizeof(short), 3, infile);
00411       if (is_little_endian()) swap2_aligned(tmp, 3);
00412       triangle[3*i+0] = tmp[0];
00413       triangle[3*i+1] = tmp[1];
00414       triangle[3*i+2] = tmp[2];
00415     }
00416     TRASH;
00417     TRASH;
00418     fread(properties, sizeof(float), 3 * nvert, infile);
00419       if (is_little_endian())
00420       swap4_aligned(properties, 3*nvert);
00421 
00422 #endif
00423   }   
00424   
00426   Get_Property_Values(&datax, properties, colores, nvert);
00427   
00428   // And draw things
00429   grasp->graphics = new molfile_graphics_t[3*ntriangles];
00430   int vert1, vert2, vert3;
00431 
00432   for (int tri_count = 0; tri_count < ntriangles; tri_count++) {
00433     vert1 = triangle[3*tri_count+0] - 1;  // from 1-based FORTRAN
00434     vert2 = triangle[3*tri_count+1] - 1;  // to 0-based C++
00435     vert3 = triangle[3*tri_count+2] - 1;
00436 
00437     if (vert1 <      0 || vert2 <      0 || vert3 <      0 ||
00438         vert1 >= nvert || vert2 >= nvert || vert3 >= nvert) {
00439       printf("graspplugin) Error, out-of-range vertex index, aborting.\n"); 
00440       delete [] vertex;
00441       delete [] access;
00442       delete [] normal;
00443       delete [] triangle;
00444       delete [] properties;
00445       return MOLFILE_ERROR;
00446     }
00447 
00448     grasp->graphics[2*tri_count  ].type = MOLFILE_TRINORM;
00449     grasp->graphics[2*tri_count+1].type = MOLFILE_NORMS;
00450     grasp->graphics[2*tri_count+2].type = MOLFILE_COLOR;
00451     
00452     float *tridata =  grasp->graphics[2*tri_count  ].data;
00453     float *normdata = grasp->graphics[2*tri_count+1].data;
00454     float *colordata = grasp->graphics[2*tri_count+2].data;
00455         
00456     memcpy(tridata  , vertex+3*vert1, 3*sizeof(float));
00457     memcpy(tridata+3, vertex+3*vert2, 3*sizeof(float));
00458     memcpy(tridata+6, vertex+3*vert3, 3*sizeof(float));
00459     
00460     memcpy(normdata  , normal+3*vert1, 3*sizeof(float));
00461     memcpy(normdata+3, normal+3*vert2, 3*sizeof(float));
00462     memcpy(normdata+6, normal+3*vert3, 3*sizeof(float));
00463     
00464     memcpy(colordata  , properties+3*vert1, 3*sizeof(float));
00465     memcpy(colordata+3, properties+3*vert2, 3*sizeof(float));
00466     memcpy(colordata+6, properties+3*vert3, 3*sizeof(float));
00467   } 
00468 
00469   *nelem = 2*ntriangles;
00470   *data = grasp->graphics;
00471 
00472   delete [] triangle;
00473   delete [] normal;
00474   delete [] access;
00475   delete [] vertex;
00476   delete [] properties;
00477 
00478   return MOLFILE_SUCCESS;
00479 }
00480 
00481 
00482 static void close_file_read(void *v) {
00483   grasp_t *grasp = (grasp_t *)v;
00484   fclose(grasp->fd);
00485   delete [] grasp->graphics;
00486   delete grasp;
00487 }
00488 
00489 
00490 /*
00491  * Initialization stuff here
00492  */
00493 static molfile_plugin_t plugin;
00494 
00495 VMDPLUGIN_API int VMDPLUGIN_init(void) {
00496   memset(&plugin, 0, sizeof(molfile_plugin_t));
00497   plugin.abiversion = vmdplugin_ABIVERSION;
00498   plugin.type = MOLFILE_PLUGIN_TYPE;
00499   plugin.name = "grasp";
00500   plugin.prettyname = "GRASP";
00501   plugin.author = "Justin Gullingsrud, John Stone";
00502   plugin.majorv = 0;
00503   plugin.minorv = 8;
00504   plugin.is_reentrant = VMDPLUGIN_THREADSAFE;
00505   plugin.filename_extension = "srf,SRF,grasp";
00506   plugin.open_file_read = open_file_read;
00507   plugin.read_rawgraphics = read_rawgraphics;
00508   plugin.close_file_read = close_file_read;
00509   return VMDPLUGIN_SUCCESS;
00510 }
00511 
00512 VMDPLUGIN_API int VMDPLUGIN_register(void *v, vmdplugin_register_cb cb) {
00513   (*cb)(v, (vmdplugin_t *)&plugin);
00514   return VMDPLUGIN_SUCCESS;
00515 }
00516 
00517 VMDPLUGIN_API int VMDPLUGIN_fini(void) { return VMDPLUGIN_SUCCESS; }
00518 
00519 
00520 

Generated on Tue Jan 21 02:56:13 2020 for VMD Plugins (current) by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002