Hi PyMOL users,

MEAD's potential maps can be used in PyMOL, but it requires a little more
than a bit of tweaking, as the maps are in the AVS .fld format, which can
not be read by PyMOL. (These maps can be read by Dino, as Paulo just
mentioned). Also, PyMOL reads (as far as I can see) only big-endian
phi-maps. 

I have written a small program, which does the conversion from .fld to
big-endian-.phi, so the maps can be read by PyMOL. It just finished it
today, and it has only been tested on Linux (i386), so beware ;-)

To create a nice electrostatic surface, you need:


* Gromacs (not strictly necessary, but makes life easier) Download from
http://www.gromacs.org 
* MEAD. Download from http://www.scripps.edu/bashford/
* fld2phi, source code quoted below (sorry, but I don't have access to our
external web server, and it's only about 4kb). 
* PyMOL, of course (v0.88 or newer).


Here is an example how to create everything from scratch for 4PTI.pdb


Gromacs steps
-------------
1) use pdb2gmx to create .gro and .top files:
pdb2gmx -f 4PTI.pdb -o 4PTI.gro -p 4PTI.top 

2) use grompp to create .tpr file:
grompp -f 4PTI.mdp -p 4PTI.top -o 4PTI.tpr -c 4PTI.gro

The file 4PTI.mdp contains the parameters for the Gromacs simulation. But as
we are not going to do any simulation this time, an empty file is ok. It can
be created with touch 4PTI.mdp. 

3) use editconf to create a MEAD-readable pdb file: 
editconf -f 4PTI.tpr -mead -o 4PTI.pqr.pdb

The output pdb file must then be renamed to be recognized by MEAD: 
mv 4PTI.pqr.pdb 4PTI.pqr 


MEAD steps
----------

1) create a .ogm file which specifies the grid size. Notice that PyMOL can
only handle grids which are 65x65x65 points, so your only option is to
change the spacing between points. You can specify focussing options in the
.ogm file, but only the coarsest grid is written anyway, so you only need
one line in the 4PTI.ogm file:

ON_GEOM_CENT 65 1.0

See documentation for other centering options (first parameter). Next number
is the number of grid points on each side (must be 65 to be readable by
PyMOL). The last number is a real specifying the distance between grid
points. 

2) Run 'potential' to create the grid:
potential -epsin 2 -CoarseFieldOut 4PTI 4PTI

The epsin option is mandatory and specifies the internal (in the protein)
dielectric constant. The program will say something like:

WARINING from potential main program:
Could not open field point file, 4PTI.fpt, for reading. Exiting without
giving any potentials.

This can be ignored. The program still writes out a 4PTI.fld file with the
grid. Notice that this file in not overwritten, so you must delete it
manually if it already exists. 

Convert to PyMOL readable grid
------------------------------

1) Use the fld2phi to convert the .fld file to a big-endian .phi grid file:
fld2phi 4PTI.fld 4PTI.phi

This sould create a 4PTI.phi file which is readable by PyMOL. 


PyMOL steps
-----------

1) Load the structure including the hydrogens built by Gromacs:
load 4PTI.pqr, 4PTI 

2) Create a selection of the water:
select water, 4PTI and resn SOL 

3) Remove the water atoms:
remove water 

4) Show the surface of the 4PTI object:
show surface, 4PTI 

5) Load the electrostatic grid:
load 4PTI.phi, map

You can show the extent of the grid box by clicking on the object called
"map" in the object list to the right. 

6) Create a color ramp object:
ramp_new e_lvl, map, [-0.02,0.00,0.02] 

7) Color the surface according to the grid and map:
set surface_color, e_lvl, 4PTI 

Thats it. You can change the color scale on the fly by issuing another
ramp_new command with other numbers. The 3 numbers are red-point,
white-point and blue-point, respectively. The scale can also be changed by
ctrl+mid-click while you drag the color scale.

It is also possible to create one or more contour surfaces: 
isosurface contour1, map, -0.05

where contour1 is the object name of the surface (choose whatever you want),
"map" is the object name of the electrostatic potential map and the number
is the contour level. The commands isomesh and isodot have the same syntax
and do exactly what you think.



Best regards (and please forgive me for quoting source code here :)

Esben


**************************************************************
*  Source code for fld2phi                                   *
*  save as: fld2phi.c                                        *
*  compile with: gcc -o fld2phi fld2phi.c                    *
**************************************************************

/* This is fld2phi.c */
/* This program reads little-endian AVS fields (.fld-files) from MEAD and
*/
/* convert it to big-endian .phi maps readable by PyMOL
*/
/* Tested on i386-Linux only */
/* (c) Esben Friis, 2003 */

/* ------------ INCLUDES --------------- */

#include <stdio.h>
#include <stdlib.h>

/* ------------- DEFINES --------------- */

#define VERSION "1.0"
#define DATE    "2003-07-09"

/* ------------- FUNCTIONS ------------- */

int get_fld_header (char *filename, int *gridpts, float *xmin, float *ymin,
float *zmin, float *xmax, float *ymax, float *zmax )
{
  char buffer[256];
  int i;
  FILE *infile;
   
  if (!(infile = fopen(filename,"rb")))
    return 0;
  else
  {
    for (i=0; i<4; i++)
      fgets(buffer, 256, infile);

    fgets(buffer, 256, infile);
    *gridpts = atoi(buffer+5);   /* the number of grid points in each
direction */

    for (i=0; i<5; i++)
      fgets(buffer, 256, infile);

    fgets(buffer, 256, infile);  /* min extent line */
    sscanf(buffer+8, "%f %f %f", xmin, ymin, zmin);

    fgets(buffer, 256, infile);   /* max extent line */
    sscanf(buffer+8, "%f %f %f", xmax, ymax, zmax);
  
    fgets(buffer, 256, infile);

    fclose (infile);     

    return 1;
  }
} 
           

int get_fld_grid (char *filename, int gridpts, float *grid)
{
  char buffer[256];
  int i, x, y, z, nobj;
  FILE *infile;

  if (!(infile = fopen(filename,"rb")))
    return 0;
  else
  {  
    for (i=0; i<13; i++)
       fgets(buffer, 256, infile);

    /* --- end of header ---- */
  
    fgets(buffer, 3, infile);   /* junk bytes */
    nobj = fread(grid, 4, gridpts*gridpts*gridpts, infile);   /* data array
*/
    
    printf ("read %d data points from %s\n", nobj, filename);

    fclose (infile);     
  }
}


void write_big_endian_float(FILE *stream, float *f)
{
   char *p;
   char rev[4];
   
   p = (char*) f;
   rev[0] = p[3];
   rev[1] = p[2];
   rev[2] = p[1];
   rev[3] = p[0];
   fwrite(&rev, 4, 1, stream);
}


int write_phi(char *filename, int gridpts, float *grid, float xmin, float
ymin, float zmin, float xmax, float ymax, float zmax)
{
   float xmid, ymid, zmid, scale;
   FILE *outfile;
   int i, x, y,z;
      
   if (!(outfile = fopen(filename,"wb")))
     return 0;
   else
   {

      xmid = (xmax+xmin)/2.0;
      ymid = (ymax+ymin)/2.0;
      zmid = (zmax+zmin)/2.0;

      scale = ((float)gridpts-1.0)/(xmax-xmin); 

      printf ("midpoint %f, %f, %f, scale = %f\n", xmid, ymid, zmid, scale);

      fprintf(outfile,"%c%c%c%c",0,0,0,0x14);
      fprintf(outfile,"now starting phimap ");
      fprintf(outfile,"%c%c%c%c",0,0,0,0x14);
      fprintf(outfile,"%c%c%c%c",0,0,0,0x46);
      fprintf(outfile,"potential ");
      fprintf(outfile,"Map converted from MEAD AVS (fld) file
");
      fprintf(outfile,"%c%c%c%c",0,0,0,0x46);
      fprintf(outfile,"%c%c%c%c",0,0x10,0xc3,0x04);

      /* order is z y x */  
      for (z=0; z<gridpts; z++)
        for (y=0; y<gridpts; y++)
          for (x=0; x<gridpts; x++)
             write_big_endian_float(outfile,
grid+gridpts*gridpts*z+gridpts*y+x);

 
fprintf(outfile,"%c%c%c%c%c%c%c%c",0x13,0x80,0x00,0x10,0xc3,0x04,0x00,0x00);

      fprintf(outfile," end of phimap  ");

      fprintf(outfile,"%c%c%c%c",0x10,0,0,0);
      fprintf(outfile,"%c%c%c%c",0x10,0,0,0);

      write_big_endian_float(outfile, &scale); 
      write_big_endian_float(outfile, &xmid);
      write_big_endian_float(outfile, &ymid);
      write_big_endian_float(outfile, &zmid);
      fprintf(outfile,"%c%c%c%c",0x10,0,0,0);

      fclose(outfile);
  }
}

void usage(void)
{
   fprintf(stderr,"USAGE:\n");
   fprintf(stderr,"fld2phi <map.fld> <map.phi>\n\n");
   fprintf(stderr,"Converts an AVS .fld grid file <map.fld> written by
MEAD's 'potential'\n");
   fprintf(stderr,"program to a big-endian Delphi grid file <map.phi> which
is readable\n");
   fprintf(stderr,"by eg. PyMOL.\n");
   fprintf(stderr,"Notice that PyMOL (v0.90) can only handle grids of
65x65x65 points.\n");
   fprintf(stderr,"fld2phi is version %s from %s\n",VERSION,DATE);
  
   exit(0);
}

int main (int argc, char **argv)
{
   float xmin, ymin, zmin, xmax, ymax, zmax;
   float *grid;
   int gridpts;
   
   if (argc !=3)
     usage();

   get_fld_header(argv[1], &gridpts, &xmin, &ymin, &zmin, &xmax, &ymax,
&zmax);
   
   printf ("%d, %f, %f, %f, %f, %f, %f\n", gridpts, xmin, ymin, zmin, xmax,
ymax, zmax);
   
   printf ("Grid size is %d.\n", gridpts);
   if (gridpts != 65)
   {
     printf ("WARNING:\n");
     printf ("PyMOL version 0.90 only reads maps which are 65x65x65!\n");
   }
   
   grid = (float*) malloc(gridpts*gridpts*gridpts*4); 

   get_fld_grid(argv[1], gridpts, grid);
   
   write_phi(argv[2], gridpts, grid, xmin, ymin, zmin, xmax, ymax, zmax);

   return 0;
}  

Reply via email to